&EPA
United States
Environmental Protection
Agency
Research and Development
Environmental Research
Laboratory
Athens GA 30613
EPA/600 3-85 065
August 1985
Computer Program
Documentation for the
Enhanced Stream
Water Quality
Model QUAL2E
> I
-------
EPA/600/3-85/065
August 1985
COMPUTER PROGRAM DOCUMENTATION FOR THE ENHANCED
STREAM WATER QUALITY MODEL QUAL2E
by
Linfield C. Brown* and Thomas 0. Barnwell, Jr.**
*Department of Civil Engineering
Tufts University
MedfoYd, MA 02155
**Environmental Research Laboratory
U.S. Environmental Protection Agency
Athens, GA 30613
Cooperative Agreement No. 811883
t^-rtrt pro
ENVIRONMENTAL RESEARCH LABORATORY
OFFICE OF RESEARCH AND DEVELOPMENT
U.S. ENVIRONMENTAL PROTECTION AGENCY
ATHENS, GEORGIA 30613
-------
DISCLAIMER
The information in this document has been funded wholly or in part by
the United States Environmental Protection Agency under Cooperative Agreement
No. 811883 to Tufts University. It has been subject to the Agency's peer and
administrative review, and it has been approved for publication as an EPA
document. Mention of trade names or commercial products does not constitute
endorsement or recommendation for use by the U.S. Environmental Protection
Agency.
ii
-------
FOREWORD
As environmental controls become more costly to implement and the
penalties of judgment errors become more severe, environmental quality
management requires more efficient management tools based on greater
knowledge of the environmental phenomena to be managed. As part of this
Laboratory's research on the occurrence, movement, transformation, impact
and control of environmental contaminants, the Assessment Branch develops
management or engineering tools to help pollution control officials achieve
water quality goals.
The stream water quality model QUAL-II is widely used for waste load
allocations, discharge permit determinations, and other conventional
pollutant evaluations in the United States. Since its introduction in
1970, several different versions of the model have evolved. This manual
presents the most recent modifications in the form of an enhanced state-of-
the-art model called QUAL2E, which was developed by the National Council
for Air and Stream Improvement (NCASI). Distribution and maintenance of
the QUAL2E computer program, and training and assistance to model users,
will be provided by EPA's Center for Water Quality Modeling at this Labora-
tory.
Rosemarie C. Russo, Ph.D.
Director
Environmental Research Laboratory
Athens, Georgia
iii
-------
ABSTRACT
Presented in this manual are recent modifications and improvements to
the widely used stream water quality model QUAL-II. Called QUAL2E, the
enhanced model incorporates improvements in eight areas: (1) algal,
nitrogen, phosphorus, and dissolved oxygen interactions; (2) algal growth
rate; (3) temperature; (4) dissolved oxygen; (5) arbitrary non-conservative
constituents; (6) hydraulics; (7) downstream boundary concentrations; and
(8) input/output modifications.
QUAL2E, which can be operated either as a steady-state or as a dynamic
model, is intended for use as a water quality planning tool. The model,
for example, can be used to study the impact of waste loads on instream
water quality or to identify the magnitude and quality characteristics of
nonpoint waste loads as part of a field sampling program. The user also
can model effects of diurnal variations in meteorological data on water
quality (primarily dissolved oxygen and temperature) or examine diurnal
dissolved oxygen variations caused by algal growth and respiration.
This report was submitted in partial fulfillment of Cooperative Agree-
ment No. 811883 by Tufts University under the partial sponsorship of the
U.S. Environmental Protection Agency. This report covers a period from
August 1984 to June 1985, and work was completed as of June 1985.
IV
-------
CONTENTS
FOREWORD ill
ABSTRACT iv
1. INTRODUCTION 1
1.1 QUAL2E Development 2
1.1.1 History 2
1.1.2 Enhancements to QUAL-II 3
1.1.3 Information Sources 4
1.1.4 Organization of This Report 5
1.2 QUAL2E Computer Model 5
1.2.1 Prototype Representation 5
1.2.2 Model Limitations 6
1.2.3 Model Structure and Subroutines 6
1.2.4 Program Language and Operating Requirements 6
2. GENERAL MODEL FORMULATION 9
2.1 Introduction 9
2.2 Conceptual Representation 10
2.3 Functional Representation 10
2.4 Hydraulic Characteristics 14
2.4.1 Discharge Coefficients 14
2.4.2 Trapezoidal Cross Sections 15
2.4.3 Longitudinal Dispersion 15
2.5 Flow Augmentation 18
3. CONSTITUENT REACTIONS AND INTERRELATIONSHIPS 21
3.1 General Considerations 21
3.2 Chlorophyll a_ (Phytoplankton Algae) 21
3.2.1 Algal Respiration Rate 23
3.2.2 Algal Specific Growth Rate 23
3.2.3 Algal-Light Relationships 25
3.2.4 Algal-Nutrient Relationships 33
3.2.5 Temperature Dependence in Algae Simulation 34
3.3 Nitrogen Cycle 34
3.3.1 Organic Nitrogen 34
3.3.2 Ammonia Nitrogen 35
3.3.3 Nitrite Nitrogen 35
3.3.4 Nitrate Nitrogen 36
3.3.5 Inhibition of Nitrification at Low DO 36
-------
CONTENTS (cont'd)
Page
3.4 Phosphorus Cycle 37
3.4.1 Organic Phosphorus 37
3.4.2 Dissolved Phosphorus 38
3.5 Carbonaceous BOD 38
3.6 Dissolved Oxygen 39
3.6.1 Dissolved Oxygen Saturation Coefficient 40
3.6.2 Atmospheric Reaeration Coefficient Estimation .... 41
3.6.3 Ice Cover 47
3.6.4 K2 Default Values 47
3.6.5 Dam Reaeration 48
3.7 Coliforms 48
3.8 Arbitrary Nonconservative Constituent ... 49
3.9 Temperature 49
3.10 Temperature Dependence of Rate Coefficients 50
3.11 Reaction Rates and Physical Constants ... 52
4. FUNCTIONAL REPRESENTATION OF TEMPERATURE 55
4.1 Basic Temperature Equation 55
4.2 Definition of HN 56
4.3 Net Short-Wave Solar Radiation 58
4.3.1 Extraterrestrial Radiation 59
4.3.2 Radiation Scattering and Absorption 61
4.3.3 Cloudiness 63
4.3.4 Reflectivity 63
4.4 Long-Wave Atmospheric Radiation 64
4.5 Water Surface Back Radiation 64
4.6 Evaporation 65
4.7 Conduction 67
5. COMPUTATIONAL REPRESENTATION 68
5.1 Prototype Representation 68
5.2 Forcing Functions 69
5.3 Model Limitations 70
5.4 Numeric Solution Technique 71
5.4.1 Formulation of the Finite Difference Schmem 71
5.4.2 Method of Solution 74
5.4.3 Boundary Conditions 76
6. REFERENCES 78
APPENDIX A. QUAL2E USERS MANUAL 82
APPENDIX B. EXAMPLE INPUT/OUTPUT DATA SETS 118
-------
LIST OF FIGURES
No. Page
1-1 General Structure of QUAL2E 7
II-l Discretized Stream System 11
II-2 Stream Network of Computational Elements and Reaches 12
III-l Major Constituent Interactions in QUAL2E .- 22
III-2 QUAL2E Light Functions 28
IV-1 Heat Transfer Terms Associated with Interfacial Heat Transfer 57
V-l Classical Implicit Nodal Scheme 71
LIST OF TABLES
No. Page
II-l Values of Manning's "n" Roughness Coefficient 17
II-2 Experimental Measurements of Longitudinal Dispersion
in Open Channels 19
III-l Comparison of Dissolved .Oxygen Saturation Concentrations . . 42
111-2 Default Temperature Correction Values for QUAL2E 51
111-3 Typical Ranges for QUAL2E Reaction Coefficients 52
IV-1 Definition of Heat Transfer Terms Illustrated in Figure IV-1 58
IV-2 Empirical Coefficients for Determing Rs 64
vii
-------
ACKNOWLEDGMENT
Over the years, many investigators have contributed to the development
of what has become QUAL2E. The foundation upon which the model has been
built was laid by the Texas Water Development Board in the late 1960s in
the QUAL-I model. Many versions of the model emerged in the 1970s. The
lineage of QUAL2E can be traced to work done for the Southeast Michigan
Council of Governments (SEMCOG) by Water Resources Engineers, Inc. (now
Camp, Dresser, McKee Inc.). QUAL-II/SEMCOG was chosen for distribution
by the Center for Water Quality Modeling (CWQM) in the late 1970s and began
to receive wide use in water quality modeling and wasteload allocation
programs.
QUAL-II/SEMCOG was throughly reviewed, tested, and documented by the
National Council of the Paper Industry for Air and Stream Improvement, Inc.,
(NCASI) as discussed in NCASI Technical Bulletin No. 391. Changes arising
from this review were incorporated in a model called QUAL-II/NCASI, which
was adopted for distribution by the Center for Water Quality Modeling.
Because of a mutual interest in the program, CWQM partially sponsored an
NCASI review of other versions of the QUAL-II computer program and incor-
porated useful features of these versions in the program called QUAL2E.
Appendix A of this documentation report, the QUAL2E users manual, is
taken directly from NCASI Technical Bulletin No. 457, "Modifications to
the QUAL-2 Water Quality Model and User Manual for QUAL2E Version 2.2."
We express our appreciation to NCASI, for permission to reprint this
material in the QUAL2E documentation report.
The QUAL2E program also has been made available for IBM PC-compatible
microcomputer. The microcomputer installation of this program was performed
by Mr. Bruce Bartell of Computer Sciences Corporation, Inc. and was made
possible through the support of Mr. King Boynton of the US EPA's Office
of Water and through an agreement with the US-Spain Joint Committee for
Scientific and Technical Cooperation.
vm
-------
1. INTRODUCTION
QUAL2E is a comprehensive and versatile stream water quality model.
It can simulate up to 15 water quality constituents in any combination
desired by the user. Constituents which can be simulated are:
1. Dissolved Oxygen
2. Biochemical Oxygen Demand
3. Temperature
4. Algae as Chlorophyll a_
5. Organic Nitrogen as N
6. Ammonia as N
7. Nitrite as N
8. Nitrate as N
9. Organic Phosphorus as P
10. Dissolved Phosphorus as P
11. Coliforms
12. Arbitrary Nonconservative Constituent
13. Three Conservative Constituents
The model is applicable to dendritic streams that are well mixed. It
assumes that the major transport mechanisms, advection and dispersion, are
significant only along the main direction of flow (longitudinal axis of the
stream or canal). It allows for multiple waste discharges, withdrawals,
tributary flows, and incremental inflow and outflow. It also has the
capability to compute required dilution flows for flow augmentation to
meet any prespecified dissolved oxygen level.
Hydraulically, QUAL2E is limited to the simulation of time periods
during which both the stream flow in river basins and input waste loads are
essentially constant. QUAL2E can operate either as a steady-state or as a
-------
dynamic model, making it a very helpful water quality planning tool. When
operated as a steady-state model, it can be used to study the impact of
waste loads (magnitude, quality and location) on instream water quality
and also can be used in conjunction with a field sampling program to
identify the magnitude and quality chracteristics of nonpoint source
waste loads. By operating the model dynamically, the user can study the
effects of diurnal variations in meteorological data on water quality
(primarily dissolved oxygen and temperature) and also can study diurnal
dissolved oxygen variations due to algal growth and respiration. However,
the effects of dynamic forcing functions, such as headwater flows or point
loads, cannot be modeled in QUAL2E.
1.1 QUAL2E DEVELOPMENT
1.1.1 History
QUAL2E, a new release of QUAL-II, was developed jointly by the
National Council for Air and Stream Improvement, Inc. (NCASI) and the EPA
Center for Water Quality Modeling (CWQM), Environmental Research Labora-
tory, Athens, GA. It includes many modifications made in the model since
the revisions to the SEMCOG version of QUAL-II (NCASI, 1980). QUAL2E is
intended to supersede all prior releases of QUAL-II.
QUAL-II is an extension of the stream water quality model QUAL-I
developed by F. D. Masch and Associates and the Texas Water Development
Board (1971) and the Texas Water Development Board (1970). In 1972, Water
Resources Engineers, Inc. (WRE) under contract to the U.S. Environmental
Protection Agency, modified and extended QUAL-I to produce the first ver-
sion of QUAL-II. Over the next 3 years, several different versions of
the model evolved in response to specific user needs. In March 1976, the
Southeast Michigan Council of Governments (SEMCOG) contracted with WRE to
make further modifications and to combine the best features of the exist-
ing versions of QUAL-II into a single model. The significant modifica-
tions made in the SEMCOG version by WRE (Roesner et al_., 1981a and b)
were:
Option of English or metric units on input data
Option for English or metric outputchoice is independent of
input units
t Option to specify channel hydraulic properties in terms of
trapezoidal channels or stage-discharge and velocity discharge
curves
t Option to use Tsivoglou's computational method for stream
reaeration
Improvement in output display routines
t Improvement in steady-state temperature computation routines
-------
The SEMCOG version of QUAL-II was later reviewed, documented, and
revised (NCASI, 1982). The revised SEMCOG version has since been main-
tained and supported by the EPA Center for Water Quality Modeling (CWQM).
In 1983, EPA, through the CWQM, contracted with NCASI to continue the
process of modifying QUAL-II to reflect state-of-the-art water quality
modeling. Extensive use of QUAL-II/SEMCOG had uncovered difficulties that
required corrections in the algal-nutrient-light interactions. In addition,
a number of modifications to the program input and output had been suggested
by users. This report describes the most recent modifications made to
QUAL-II and fully documents the enhanced model, now named QUAL2E.
1.1.2 Enhancements to QUAL-II
QUAL2E, a modified version of QUAL-II/SEMCOG, incorporates improvements
in eight areas:
1. Algal, nitrogen, phosphorus, dissolved oxygen interactions
Organic nitrogen state variable
Organic phosphorus state variable
Nitrification inhibition at low DO
Algal preference factor for NH3
2. Algal growth rate
Growth rate dependent upon both NH3 and N03 concentrations
Algal self-shading
Three light functions for growth rate attenuation
Three growth rate attenuation options
Four diurnal averaging options for light
Temperature
Link to algal growth via solar radiation
Default temperature correction factors
Dissolved Oxygen (DO)
New Standard Methods DO saturation function
Traditional SOD units (g/m2-day or g/ft2-day)
Dam reaeration option
Arbitrary non-conservative constituent
First order decay
Removal (settling) term
Benthal source term
Hydraulics
t Input factor for longitudinal dispersion
-------
t Test for negative flow (i.e., withdrawal greater than flow)
t Capability for incremental outflow along reach
7. Downstream boundary
Option for specifying downstream boundary water quality
constituent concentrations
8. Input/output modifications
Detailed summary of hydraulic calculations
New coding forms
Local climatological data echo printed
Enhanced steady-state convergence
Five part final summary including components of DO
deficit and plot of DO and BOD
1.1.3 Information Sources
Major sources of information for this revised documentation were:
1. Roesner, L. A., Giguere, P. R. and Evenson, D« E. Computer
Program Documentation for Stream Quality Modeling (QUAL-II).
U.S. Environmental Protection Aghens, GA. EPA-600/9-81-014,
February 1981.
2. 'JRB Associates. Users Manual for Vermont QUAL-II Model.
Prepared for U.S. Environmental Protection Agency, Washington,
DC. June 1983.
3. 'National Council for Air and Stream Improvement. A Review of
the Mathematical Water Quality Model QUAL-II arid Guidance for
its Use. NCASI. Newg York, NY. Technical Bulletin No. 391,
Decemberr 1982.
This documentation of QUAL2E consolidates material from these and
other sources into a single volume. The basic theory and mechanics behind
the development of QUAL2E are described in this volume, which is intended
to supplement the users manual (NCASI, 1985). The QUAL2E user manual con-
tains a detailed description of data requirements, as well as the input
coding forms and an example input/output data file. Both reports and a
copy of the QUAL2E computer code are available from the Center for Water
Quality Modeling, U.S. Environmental Protection Agency, Athens, GA 30613.
-------
1.1.4 Organization of this Report
The general program structure, specifications, and limitations of
QUAL2E are discussed in the remainder of this chapter. Chapter 2 des-
cribes the conceptual and functional representation of QUAL2E as well as
the hydraulic characteristics of the model. The mathematical basis of
the water quality constituent formulations is presented in Chapter 3.
Chapter 4 presents the frame work for modeling temperature. It is ex-
tracted verbatim from Roesner jjt aj_., 1981. Chapter 5 describes the com-
putational representation of the model and the numerical solution algorithm,
A list of references used in this report is found in Chapter 6.
For the convenience of the majority of users, all of the units speci-
fications are given in the English system of measurement. QUAL2E, however,
will recognize either English or metric units.
1.2 QUAL2E COMPUTER MODEL
1.2.1 Prototype Representation
QUAL2E permits simulation of any branching, one-dimensional stream
system. The first step in modeling a system is to subdivide the stream
system into reaches, which are stretches of stream that have uniform hy-
draulic characteristics. Each reach is then divided into computational
elements of equal length. Thus, all reaches must consist of an integer
number of computational elements.
There are seven different types of computational elements:
1. Headwater element
2. Standard element
3. Element just upstream from a junction
4. Junction element
5. Last element in system
6. Input element
7. Withdrawal element
Headwater elements begin every tributary as well as the main river system,
and as such, they must always be the first element in a headwater reach.
A standard element is one that does not qualify as one of the remaining
six element types. Because incremental flow is permitted in all element
types, the only input permitted in a standard element is incremental
flow. A type 3 element is used to designate an element on the main stem
just upstream of a junction. A junction element (type 4) has a simulated
tributary entering it. Element type 5 identifies the last computational
-------
element in the river system; there should be only one type 5 element.
Element types 6 and 7 represent inputs (waste loads and unsimulated
tributaries) and water withdrawals, respectively. River reaches, which
are aggregates of computational elements, are the basis of most data
input. Hydraulic data, reaction rate coefficients, initial conditions,
and incremental flows data are constant for all computational elements
within a reach.
1.2.2 Model Limitations
QUAL2E has been designed to be a relatively general program; however,
certain dimensional limitations have been imposed during program develop-
ment. These limitations are:
t Reaches: a maximum of 50
Computational elements: no more than 20 per reach or a total
of 500
Headwater elements: a maximum of 10
Junction elements: a maximum of 9
t Input and withdrawal elements: a maximum of 50
QUAL2E incorporates features of ANSI FORTRAN 77 that allow these limita-
tions to be easily changed.
1.2.3 Model Structure and Subroutines
QUAL2E is structured as one main program supported by 49 different
subroutines. Figure 1-1 illustrates the functional relationships between
the main program and the subroutines. New state variables can be added
or modifications to existing relationships can be made with a minimum of
model restructuring through the simple addition of appropriate subroutines.
The structural framework of QUAL2E has been modified from prior ver-
sions. The large MAIN program and subroutine INDATA have been divided
into smaller groups of subroutines, each with a more narrowly defined
task. The new subroutines in QUAL2E include the algal light functions
(GROW/LIGHT), the steady state algal output summary (MRPT1), the organic
nitrogen and phosphorus state variables (NH2S, PORG), and the line printer
plot routine (PRPLOT). This reorganization of QUAL2E into smaller pro-
grammatic units is the first step in adapting the model to micro and
minicomputers, which have limited space for memory.
1.2.4 Program Language and Operating Requirements
QUAL2E is written in ANSI FORTRAN 77 and is compatible with mainframe
and personal computer systems that support this language. QUAL2E typically
-------
»w Augmentation O 1
u.
0
I
E
oc
2
a.
UAL2E
k
c Simulation
Program Return Loop for Dynam
Steady State ConOverge
1
H INDATA
HTDIUAT
M REAERC
k
Steady State
Convergence Loop
^ 1
t
ft
LGHT
ft| bos
1 ft| sbVMt
ft! WRPT1
M WRPTI M WRPT2 1
ftl FLOAUG
wiF^F
F|_
H
H
H
-H
M
ft|
L-H
^rj
^£|
^
H
ftfFFPLoT
" 'C1JDC
LMro
1- EATEX
TEMPSS
HEATER
CONSVT
BODS
COLIS
ANCS
/* )AIA/
\J )00
I )1A
> )02
« W
\ )0^
' )05
> )06
s )07
^ )l)9
i )10
s )1 1
* )12
» )13
o w
v R
A T
t 2
Figure 1-1. General Structure of QUAL2E
-------
requires 256K bytes of memory and uses a single system input device
(cards or disk file) and the system's line printer (or disk file) as the
output device.
If the system's normal FORTRAN input device unit is not unit 5 or
the output unit is not unit 6, then the variables "NI" and "NJ" in the
subroutine INDATA should be changed to reflect the system's I/O unit
identifiers.
-------
2. GENERAL MODEL FORMULATION
2.1 INTRODUCTION
The primary objective of any stream water quality model development
is to produce a tool that has the capability for simulating the behavior
of the hydrologic and water quality components of a stream system. The
development of this tool to simulate prototype behavior by applying
mathematical model on a digital computer proceeds through three general
phases (Water Resources Engineers, Inc., 1967):
1. Conceptual representation
2. Functional representation
3. Computational representation
Conceptual representation involves a graphic idealization of the
prototype by description of the geometric properties that are to be
modeled and by identification of boundary conditions and interrelation-
ships between various parts of the prototype. Usually, this process
entails dividing the prototype into discrete "elements" of a size com-
patible with the objectives that the model must serve, defining these
elements according to some simple geometric rules, and designating the
mode by which they are connected, either physically or functionally, as
integral parts of the whole. A part of this conceptual structuring is
the designation of those boundary conditions to be considered in the
simulation.
Functional representation entails formulation of the physical fea-
tures, processes, and boundary conditions into sets of algebraic equations.
It involves precise definition of each variable and its relationship to
all other parameters that characterize the model or its input-output
relationships.
Computational representation is the process whereby the functional
model is translated into the mathematical forms and computational proce-
dures required for solution of the problem over the desired time and
space continuum. It is concerned with development of a specific solution
technique that can be accommodated by the computer and with codification
of the technique in computer language.
In the remainder of this section the Conceptual Representation of
QUAL2E will be described together with its general Functional Representa-
-------
tion for mass transport, hydraulic characteristics, and longitudinal
dispersion. Chapter 3 will discuss specific constituent reactions and
interactions. Chapter 4 will develop the Functional Representation of
stream temperature as simulated in QUAL2E.
2.2 CONCEPTUAL REPRESENTATION
Figure II-l shows a stream reach (n) that Tias been divided into a
number of subreaches or computational elements, each of length AX. For
each of these computational elements, the hydrologic balance can be
written in terms of flows into the upstream face of the element (Q-j-i),
external sources or withdrawals (Qxi), and the outflow (Q-j) through the
downstream face of the element. Similarly, a materials balance for any
constituent C can be written for the element. In the materials balance,
we consider both transport (Q'C) and dispersion (A £L aC) as the movers
AX" 3x
of mass along the stream axis. Mass can be added to or removed from the
system via external sources and withdrawals (QxCx)-j and added or removed
via internal sources or sinks (Sj) such as benthic sources and biological
transformation. Each computational element is considered to be completely
mixed.
Thus, the stream can be conceptualized as a string of completely
mixed reactorscomputational elementsthat are linked sequentially to
one another via the mechanisms of transport and dispersion. Sequential
groups of these reactors can be defined as reaches in which the computa-
tional elements have the same hydrogeometric propertiesstream slope,
channel cross section, roughness, etc.--and biological rate constants
BOD decay rate, benthic source rates, algae settling rates, etc.--so that
the stream shown at the left of Figure 11-2 can be conceptually represented
by the grouping of reaches and computational elements shown on the right
of Figure 11-2.
2.3 FUNCTIONAL REPRESENTATION
2.3.1 Mass Transport Equation
The basic equation solved by QUAL2E is the one dimensional advection-
dispersion mass transport equation, which is numerically integrated over
space and time for each water quality constituent. This equation includes
the effects of advection, dispersion, dilution, constituent reactions and
interactions, and sources and sinks. For any constituent, C, this equation
can be written as:
3M 3(AyDi 3x) 3(AX u C) dC
_ = JJLi - 1 dx - - dx + (Axdx) ~ + s II-l
at 8x 3x dt
10
-------
Computational
Element i
(QC)S
FLOW
BALANCE
MASS
BALANCE
Figure II-l. Discretized Stream System
-------
Number
Coflipwlotlonol
Numbtr
Figure 11-2. Stream Network of Computational Elements and Reaches
12
-------
where
M = mass (M)
x = distance (L)
t = time (T)
C = concentration (M L~3)
Ax = cross-sectional area
D|_ = dispersion coefficient (L2 T'1)
"u = mean velocity (L T-1)
s = external source or sinks (M T~l)
Because M = VC, we can write
3M a(vc) ac 3V
_ = __ = v + C II-2a
3t 3t 3t 3t
where
V = AX dx = incremental volume
If we assume that the flow in the stream is steady, i.e., 8Q/3t = 0, then
the term 3V/3t = 0 and equation II-2a becomes
3M 3C
-- = V -- II-2b
at at
Combining equations II-l and II-2b and rearranging,
3T
8C a(AxDL lx) 8(AX "u C) dC s
at Ax ax AX ax dt v
The terms on the right-hand side of the equation represent, respec-
tively, dispersion, advection, constituent changes, external sources/sinks,
and dilution. The dC term refers only to constituent changes such as
dt 3C
growth and decay, and should not be confused with the term , the local
at
concentration gradient. The latter term includes the effect of constituent
changes as well as dispersion, advection, sources/sinks, and dilutions.
13
-------
Under steady-state conditions, the local derivative becomes equal to
zero; in other words:
8C
= 0 II-4
at
Changes that occur to individual constituents or particles independent of
advection, dispersion, and waste inputs are defined by the term
dC
= individual constituents changes 11-5
dt
These changes include the physical, chemical, and biological reactions and
interactions that occur in the stream. Examples of these changes are
reaeration, algal respiration and photosynthesis, and coliform die-off.
2.4 HYDRAULIC CHARACTERISTICS
QUAL2E assumes that the stream hydraulic regime is steady-state;
i.e., 3Q/^t = 0, therefore, the hydrologic balance for a computa-
tional element can be written simply as (see Figure II-l):
( - (QX). .1-6
where (Qx). is the sum of the external inflows and/or withdrawals to that
element, i
2.4.1 Discharge Coefficients
Once equation II-6 has been solved for Q, the other hydraulic
characteristics of the stream segments can be determined by equations of
the form:
u = aQb II-7
Ax = Q/U 1 1-8
and
d = aQB 1 1-9
where a, b, a and e are empirical constants, and d is the stream depth.
These constants usually can be determined from stage-discharge rating
curves.
14
-------
2.4.2 Trapezoidal Cross Sections
Alternatively, if the cross-sectional properties of the stream segment
are available as a function of the depth d, u can be obtained as a function
of discharge by the trial and error solution of Mannings equation:
1>486 J>n i/?
Q = _ Ax Rx2/3 Se !/2 11-10
n
where
AX = cross-section area of the channel or canal, ft2
Rx = mean effective hydraulic radius, ft
n = Manning roughness factor (usual range 0.010 to 0.10)
Se = slope of the energy grade line (unit!ess)
Q = discharge, ft^/sec
The value for "u is then determined from equation II-8.
2.4.3 Longitudinal Dispersion
Dispersion is basically a convective transport mechanism. The term
"dispersion" is generally used for transport associated with spatially
averaged velocity variation, as opposed to "diffusion," which is reserved
for transport that is associated primarily with time-averaged velocity
fluctuations.
Taylor (1956) was able to derive a predictive equation for the
longitudinal dispersion coefficient, DL, in long straight pipes, as
DL = 10 r0 u*, ft2/sec 11-11
where r0 is the pipe radius and u* is the average shear velocity defined
as
u* = /T0/p, ft/sec 11-12
where
TO = boundary shear stress, lb/ft2, and
p = mass fluid density, Ib-sec2/ft4
Some investigators have attempted to apply Taylor's expression to stream-
flow. Such applications are only approximate, however, because of the
difference between the geometry or velocity distributions in streamflow
and those in a pipe.
15
-------
Elder (1959).assumed that only the vertical velocity gradient was
important in streamflow and developed an expression analogous to Taylor's
expression:
DL = Kdu* 11-13-
where d is the mean depth in feet of the stream. Elder used a value of
5.93 for K in this equation.
Other investigators have derived similar expressions for D|_ and found
it to be extremely sensitive to lateral velocity profiles. Elder's
expression, however, seems adequate in one-dimensional situations where
the channel is not too wide. For very wide channels, Fisher (1964) has
shown that half-width rather than depth is the dominant scale and there-
fore is important to the definition of the longitudinal dispersion coeffi-
cient. Equations 11-11 and 11-13 can be written in terms of the Manning
Equation and other variables characteristic of stream channels.
As an example, for steady-state open-channel flow.
u* = C ^~RSe~ 11-14
where
C = Chezy's coefficient
R = the hydraulic radius
Se = the slope of the energy grade line
Chezy's coefficient is given by:
Rl/6
C = - 11-15
where n is the Manning roughness coefficient tabulated for different types
of channels in Table II-l.
Se, the slope of the energy gradient, is given by
u" n 9
S = ( O2 H-16
6 1.486 R2/3
where "u is the mean velocity. Substituting equations 11-14, 11-15 and II-
16 into equation 11-13 and letting R = d for a wide channel yields the
expression
DL = 3.82 K n u" d5/6 11-17
16
-------
TABLE II-l
VALUES OF MANNING'S "n" ROUGHNESS COEFFICIENT
After Henderson (1966)
Artificial Channels n
Glass, plastic, machined metal 0.010
Dressed timber, joints flush 0.011
Sawn timber, joints uneven 0.014
Cement plaster 0.011
Concrete, steel troweled 0.012
Concrete, timber forms, unfinished 0.014
Untreated gunite 0.015-0.017
Brickwork or dressed masonry 0.014
Rubble set in cement 0.017
Earth, smooth, no weeds 0.020
Earth, some stones, and weeds 0.025
Natural River Channels n
Clean and straight 0.025-0.030
Winding with pools and shoals 0.033-0.040
Very weedy, winding and overgrown 0.075-0.150
Clean straight alluvial channels 0.031 d!/6
(d = D-75 size in ft.
= diameter that 75
percent of parti-
cles are smaller
than)
17
-------
where
DL = longitudinal dispersion coefficient, ft2/sec
K = dispersion constant (dimensionless)
n = Manning's roughness coefficient (dimensionless)
u = mean velocity, ft/sec
d = mean depth, ft
Typical values for dispersion coefficients, D|_, and values of the
dispersion constant, K, cited by Fisher et al. (1979), are given in Table
11-2. Note that the dispersion constant, K, shown in this table is one
to three orders of magnitude greater than that used by Elder.
2.5 Flow Augmentation
When the DO concentration in a stream drops below some required target
level, such as the state water quality standard for DO, it may be desirable
to raise this DO concentration by augmenting the flow of the stream.
According to the originators of the flow augmentation routine in QUAL2E,
Frank D. Masch and Associates and the Texas Water Development Board
(1971), the amount of flow necessary to bring the DO concentrations up to
required standards cannot be calculated by an exact functional relationship.
A good approximation of the relationship is used in QUAL2E and has the
following quadratic form:
- DOT - D0min 11-18
and
DOR
Q = Q [ + 0.15 ( )2] 11-19
DOT DOT
where,
= dissolved xoygen concentration required to meet target
conditions, mg/L
DOj = required target level of DO, mg/L
D0min = minimum DO concentration (critical level) in the oxygen
sag curve, mg/L
QR = amount of flow augmentation required, ft3/sec
QQ = flow at the critical point in the oxygen sag curve,
ft3/sec
18
-------
TABLE II-2
EXPERIMENTAL MEASUREMENTS OF LONGITUDINAL DISPERSION IN OPEN CHANNELS
(After Table 5.3, Fisher et al., 1979)
Channel
Chicago Ship Channel
Sacramento River
River Derwent
South Platte River
Yuma Mesa A Canal
Trapezoidal Laboratory
Channel with roughened
sides
Green-Duwarmish River
Missouri River
Copper Creek (below gage)
Clinch River
Copper Creek (above gage)
Powell River
Cinch River
Coachella Canal
Bayon Anacoco
Nooksack River
Wind/Bighorn Rivers
John Day River
Comite River
Sabine River
Yadkin River
Depth
d
(ft)
26.5
13.1
0.82
1.5
11.3
0.115
0.154
0.115
0.115
0.069
0.069
3.61
8.86
1.61
2.79
1.61
2.79
6.89
6.89
1.31
2.79
1.90
5.12
3.08
2.98
2.49
3.61
7.09
1.90
8.10
1.41
6.69
15.6
7.71
12.6
Width
W
(ft)
160
--
--
--
1.31
1.41
1.31
1.12
1.08
0.62
66
66
52
59
52
154
197
174
62
112
118
79
85
121
210
194
226
82
112
52
341
417
230
236
Mean
Velocity
u
(ft/sec)
0.89
1.74
1.25
2.17
2.23
0.82
1.48
1.48
1.44
1.48
1.51
5.09
0.89
1.97
0.85
1.05
3.08
2.62
0.52
0.49
0.69
2.33-
1.12
1.31
2.20
2.89
5.09
3.31
2.69
1.21
1.90
2.10
1.41
2.49
Shear
Velocity
u*
(ft/sec)
0.063
0.17
0.46
0.23
1.13
0.066
0.118
0.115
0.114
0.108
0.127
0.16
0.24
0.26
0.33
0.26
022
034
0.35
0.38
0.18
0.16
0.14
0.22
0.22
0.89
0.39
0.56
0.46
0.59
0.16
0.16
0.26
0.33
0.43
Dispersion
Coefficient
(ft2/Sec)
32
161
50
174
8.2
1.3
2.7
4.5
0.8
4.3
2.4
70-92
16,000
215
226
102
151
581
506
97
102
87
103
355
420
377
452
1722
151
700
151
3390
7200
1200
2800
Dispersion
Constant
K
20
74
131
510
8.6
174
150
338
205
392
270
120-160
75000
500
250
245
235
245
210
220
200
280
140
524
640
170
318
436
172
146
650
3100
1800
470
520
19
-------
The model augments the stream flow by first comparing, after steady-
state conditions have been reached, the simulated DO concentration with
the prespecified target level of DO in each reach. If the calculated DO
is below the target level, the program finds those upstream sources that
the user has specified for dilution purposes, and adds water equally from
all these sources. The DO calculations are then repeated. This process
continues until the DO target level is satisfied. (NOTE: The flow
augmentation subroutine can be used for DO only.)
20
-------
3. CONSTITUENT REACTIONS AND INTERRELATIONSHIPS
3.1 GENERAL CONSIDERATIONS
One of the most important considerations in determining the waste-
assimilative capacity of a stream is its ability to maintain an adequate
dissolved oxygen concentration. Dissolved oxygen concentrations in streams
are controlled by atmospheric reaeration, photosynthesis, plant and animal
respiration, benthal demand, biochemical oxygen demand, nitrification,
salinity, and temperature, among other factors.
The most accurate oxygen balance would consider all significant fac-
tors. The QUAL2E model includes the major interactions of the nutrient
cycles, algae production, benthic oxygen demand, carbonaceous oxygen
uptake, atmospheric aeration and their effect on the behavior of dissolved
oxygen. Figure III-l illustrates the conceptualization of these interac-
tions. The arrows on the figure indicate the direction of normal system
progression in a moderately polluted environment; the directions may be
reversed in some circumstances for some constituents. For example, under
conditions of oxygen supersaturation, which might occur as a result of
algal photosynthesis, oxygen might be driven from solution, opposite to
the indicated direction of the flow path.
Coliforms and the arbitrary nonconservative constituent are modeled
as nonconservative decaying constituents and do not interact with other
constituents. The conservative constituents, of course, neither decay nor
interact in any way with other constituents.
The mathematical relationships that describe the individual reactions
and interactions are presented in the following paragraphs.
3.2 CHLOROPHYLL A (PHYTOPLANKTONIC ALGAE)
Chlorophyll ^ is considered to be directly proportional to the
concentration of phytoplanktonic algal biomass. For the purposes of this
model algal biomass is converted to chlorophyll ^ by the simple relation-
ship:
Chi 1 = a0 A III-l
21
-------
where
Chi ji = chlorophyll ^concentration, ug-Chl a/L
A = algal biomass concentration, mg-A/L
a0 = a conversion factor, ug Chi ^/mg A
The differential equation that governs the growth and production of algae
(chlorophyll aj is formulated according to the following relationship.
Atmospheric
Reaeration
Figure III-l. Major Constituent Interactions in QUAL2E
22
-------
dA oj
= yA - pA - A III-2
dt d
where
A = algal biomass concentration, mg-A/L
t = time, day
v = the local specific growth rate of algae as defined below,
which is temperature dependent, day*
p = the local respiration rate of algae, which is temperature
dependent, dayl
<*1 = the local settling rate for algae, which is temperature
dependent, ft/day
d = average depth, ft
3.2.1 Algal Respiration Rate
In QUAL2E, the single respiration rate parameter, p, is used to
approximate three processes: (a) the endogenous respiration of algae, (b)
the conversion of algal phosphorus to organic phosphorus, and (C) the con-
version of algal nitrogen to organic nitrogen. No attempt is made to use
separate rate coefficients for these three processes, as is done in the
State of Vermont, revised Meta Systems version of QUAL-II (JRB Associates,
1983; and Walker, 1981).
3.2.2 Algal Specific Growth Rate
The local specific growth rate of algae, y, is known to be coupled
to the availability of required nutrients (nitrogen and phosphorus) and
light. A variety of mathematical expressions for expressing multiple
nutrient-light limitations on algal growth, rate have been reported (De
Groot, 1983; Scavia and Park, 1976; and Swartzman and Bentley, 1979).
QUAL2E has the capability of modeling the interaction among these limiting
factors in three different ways.
Growth Rate Option 1. Multiplicative. The kinetic expressions used
to represent the effects of nitrogen, phosphorus, and light are multiplied
together to determine their net effect on the local algal growth rate.
This option has as its biological basis the multiplicative effects of
enzymatic processes involved in photosynthesis:
Umax (FD () (Fp) ni-3a
23
-------
where
Mmax = maxlnium specific algal growth rate, day"-*-
FL = algal growth limitation factor for light
FN = algal growth limitation factor for nitrogen
FP = algal growth limitation factor for phosphorus
This formulation is used in the SEMCOG version of QUAL-II.
Growth Rate Option 2. Limiting Nutrient. This option represents the
local algal growth rate as limited by light and either nitrogen or. phosphorus,
but not both. Thus, the nutrient/light effects are multiplicative, but the
nutrient/nutrient effects are alternate. This formulation mimics Liebig's
law of the minimum:
v = Umax (FL) Min (FN.FP) III-3b
Thus, the algal growth rate is controlled by the nutrient (N or P) with the
smaller growth limitation factor. This option is used in the State of
Vermont version of QUAL-II.
Growth Rate Option 3. Inverse Additive. This option, a compromise
between options 1 and 2, is a modification of an intuitive form suggested
by Scavia and Park (1976) and is mathematically analogous to the total
resistance of two resistors in series. In this option, an effective nutrient
limitation factor is computed as the average of the inverse reciprocals of
the individual nitrogen and phosphorus growth limitation factors, i.e.,
1/FN + 1/FP
Thus, the algal growth rate is controlled by a multiplicative relation
between light and nutrients, but the nutrient/nutrient interactions are
represented by an inverse average. This option has been used by Water
Resources Engineers in the application of a QUAL-II-like model, WREDUN, to
Lake Dunlap (Brandes and Stein, no date).
Walker (1983) has cautioned against using the inverse-additive option
in systems where one nutrient is in excess (say nitrogen, so that FN -» 1.0)
and the other is extremely limiting (say phosphorus, so that FP -» 0.0).
In this case the value of the nutrient attenuation factor approaches 2
FP, rather than FP, as expected.
24
-------
3.2.3 Algal -Light Relationships
3.2.3.1 Light Functions
A variety of mathematical relationships between photosynthesis and
light have been reported in the literature (Jassby and Platt, 1976;
Field and Effler, 1982). Although they differ in mathematical form, the
relationships exhibit similar characteristics. All show an increasing
rate of photosynthesis with increasing light intensity up to a maximum or
saturation value. At high light intensities, some of the expressions
exhibit photoinhibition, whereas others show photosynthetic activity
remaining at the maximum rate.
QUAL2E recognizes three options for computing the algal growth limi-
tation factor for light, FL in Equations III-3a,b,c. Light attenuation
effects on the algal growth rate may be simulated using a Monod half-
saturation method, Smith's function (Smith, 1936), or Steele's equation
(Steele, 1962).
Light Function Option 1. Half Saturation. In this option, the algal
growth limitation factor for light is defined by a Monod expression:
FL Z
Z ...... .. -
KL + iz
where
FLZ = algal growth attenuation factor for light at intensity I2
Iz = light intensity at a given depth (z), Btu/ft2-hr
KL = half saturation coefficient for light, Btu/ft2-hr
z = depth variable, ft
Light Function Option 2. -Smith's Function. In this option, the algal
growth limitation factor for light is formulated to include second order
effects of light intensity:
FL = _ - - III-4b
(KL2 + I,2)1/2
where
K|_ = light intensity corresponding to 71% of the maximum growth
rate, Btu/ft2-hr
with the other terms as defined in Equation III-4a.
25
-------
Light Function Option 3. .Steel's Equation. This option incorporates
an exponential function to model the effect of photoinhibition on the algal
growth rate:
Jz Iz
FLZ = ( ) exp (1 - «-) III-4c
KL KL
where
KL = saturation light intensity at which the algal growth rate is
a maximum, Btu/ft^-hr
with the other terms as defined in Equation III-4a.
Note: The parameter KL, which appears in all three light function equations
is defined differently in each.
All of the light functions in Equations III-4a,b,c express the value
of FL for an optically thin layer. In QUAL2E photosynthesis occurs throughout
the depth of the water column. Light intensity varies with depth according
to Beer's law:
Iz = I exp (-x z) III-5
where
Iz = light intensity at a given depth (z), Btu/ft2-hr
I = surface light intensity, Btu/ft2-hr
\ = light extinction coefficient, ft~l
z = depth variable, ft
When Equation III-5 is substituted into Equations III-4a,b,c and
integrated over the depth of flow, the depth-averaged light attenuation
factor is obtained. The resulting expressions for the three options are:
Option 1: Half Saturation
KL + I
FL = (I/Ad) In [ i] III-6a
KL + Ie'xd
KL = light intensity at which growth rate is 50%
of the maximum growth rate.
26
-------
Option 2: Smith's Function
+ (1 +
FL = (1/Xd) ln[ - - - . .] III-6b
I/KLe-xd + (1 +
K|_ = light intensity at which growth rate is
of the maximum growth rate.
Option 3: Steel's Equation
-I/K
- e
._
L]
KL = light intensity at which growth rate is
equal to the maximum growth rate.
where
FL = depth-averaged algal growth attenuation factor for light
KL = light saturation coefficient, Btu/ft2-hr
A = light extinction coefficient, ft~l
d = depth of flow, ft
I = surface light intensity, Btu/ft2-hr
The relative merits of these light functions are discussed by various
authors (Bannister, 1974; Platt et al_. , 1981; Swartzmann and Bentley, 1979;
and Field and Effler, 1982). The half saturation method is the form used
in the SEMCOG version of QUAL-II. Evidence shows that the use of Smith's
function is preferrable over the half saturation method if photoinhibition
effects are unimportant (Jassby and Platt, 1976). The mathematical forms
of Equations III-4a,b,c are compared graphically in Figure 1 1 1-2. All
three equations have a single parameter, Ki_; however, it is defined differ-
ently in each equation. In Figure II 1-2 trie values of KL are selected so
that each curve passes through a common point, namely FL = 0.5 at I = 5
intensity units (i.e., a half saturation rate equal to 5 light intensity
units).
3.2.3.2 Light Averaging Options
Steady state algal simulations require computation of an average value
of FL, the growth attenuation factor for light, over the diurnal cycle.
There are four options in QUAL2E for computing this average. The options
27
-------
arise from combinations of situations regarding two factors:
The source of the solar radiation data used in the computation,
i.e., whether it is supplied externally by the user or calcu-
lated internally in the temperature heat balance.
The nature of the averaging process, i.e., whether hourly values
of FL are averaged, or a single daylight average value of solar
radiation is used to estimate the mean value of FL.
The four daily light averaging options are defined below. In
each case, the half saturation light function is used as an example; in
practice any of the three light functions may be employed.
Option 1: FL is computed from one daylight average solar radiation
value calculated in the steady state temperature heat balance:
FL
= AFACT * f *
FL
Ad
Ia1
1.0-r
Saturation
Half Saturation
1 = Half Saturation ; KL = 5.0
2 = Smith's Function ; KL = 8.66
3 = Steele's Equation ; KL = 21.55
o
0.0
Light Intensity, I (arbitrary units)
Figure 111-2. QUAL2E Light Functions
28
-------
Ialg = TFACT* Iternp III-7C
where
FL = algae growth attenuation factor for light, adjusted for
daylight hours and averaging method
AFACT = a light averaging factor, used to provide similarity
between calculations using a single average value of
solar radiation and computations using the average of
hourly values of FL
f = fraction of daylight hours
Fl_i = growth attenuation factor_for light, based on daylight
average light intensity (Iaig
X = light extinction coefficient, ft'1
d = mean depth of stream, ft
KL = half saturation coefficient for light, Btu/ft2-hr
Taig = daylight average, photosynthetically active, light
intensity, Btu/ft2-hr
TFACT = fraction of solar radiation computed in the temperature
heat balance that is photosynthetically active
Itemp = daylight average light intensity as computed in the
temperature heat balance, Btu/ft2-hr
Option 2: FL is computed from one daylight average solar radiation
value supplied externally by the user. The calculations required to ob-
tain FL iji option 2 are the same as those for option 1, except that the
value of Ia]g is computed directly from user input of photosynthetically
active solar radiation:
Talg = kot/N n 1-8
where
Itot = total daily photosynthetically active solar radiation,
Btu/ft2
N = number of daylight hours per day, hr
Both Itot and N are supplied by the user as input information.
Equations III-8, III-7b, and III-7a are used to compute the value of FL.
Because the user input value of It0t is assumed to be the photosyntheti-
cally active radiation, the factor TFACT is not used in option 2.
29
-------
Option 3: FL is obtained by averaging the hourly daylight values of
FL that are computed from the hourly daylight values of solar radiation
calculated in the steady state temperature heat balance:
FL = f * FL2 III-9a
1 N 1 KL + lalg.i
FL2 = - I [ - 9' i] III-9b
Ni=l xd K + Ie-*d
, =TFACT*
where
FL2 = average of N hourly values of FL, based on
hourly values of light intensity (Iaig,i)
Ia]q -j = hourly value of photosynthetical ly active light
intensity, Btu/ft2-hr
Hemp i = hourly value of light intensity as computed in
the steady state temperature heat balance, Btu/
ft2-hr
with other terms are defined in Equations III-7a,b,c, and III-8.
Because the average FL computed in option 3 (and 4) is an average of
diurnal ly varying values of FL, the factor AFACT is not used in the
calculations.
Option 4: FL is obtained by averaging the hourly daylight values of
FL that are computed from the hourly daylight values of solar radiation
calculated from a single value of total daily, photosynthetical ly active,
solar radiation and an assumed cosine function. The calculations required
to obtain FL are the same as those for option 3, except that the values of
Jalg.i are computed from an internally specified cosine function:
COS 2 wi
(1 -- ) , 1 = 1,N IH-10
N + 1
As in the case of option 2, both I^0t and N are supplied by the user.
Equations 111-10, III-9b, and III-9a are then used to compute the value
of FL. Because the user specified value of Itot is assumed to be photo-
synthetical ly active, the factor TFACT is not used with option 4.
30
-------
Three empirical factorsdiurnal cosine function, AFACT, and TFACT--
are used in the formulations of the four light averaging options.
Two diurnal cosine functions were evaluated for use in QUAL2E : (1)
a modified form of the one in the SEMCOG version of QUAL-II, and (2) the
form used in QUAL-TX (Texas Water Development Board, 1984). The function
in SEMCOG was modified to produce non-zero solar radiation values for each
daylight hour, as given in Equation 111-10. The form used in QUAL-TX is:
Itot *(i-l) *1
- ' [COS(~ ) - COS(-)] . 1-l.N IH-H
2N N N
Equations 111-10 and III-ll were evaluated by comparing simulated
values of FL from modeling options 2 and 4 (i.e., in effect computing
values of AFACT). Simulations were performed over a range of values of
Ki_, A, d, I^ot, anc' N» as we^ as f°r eacn °f the three light functions.
The values of AFACT averaged 0.92 and 0.94 for the SEMCOG and Texas
equations, respectively. There was no compelling reason to include both
functions (the user specified the one to be used). The diurnal cosine
function used in QUAL2E, therefore, is the modified SEMCOG version given
in Equation 111-10.
AFACT is the adjustment factor accounting for the nonlinear averaging
inherent in computing a daily average value of FL. From the simulations
just described, a resonable value of AFACT is 0.92, with a range from 0.85
to 0.98. Zison et al_. (1978) report ah implied value of 1.0 (Eq. 3.33),
and Walker (1983j~suggests using a value of 0.85.
TFACT is the photosynthetically active fraction of total solar radia-
tion. When performing algae simulations, it is important that the value
of.light intensity and light saturation coefficient, Ki, be in units of
photosynthetically active radiation, PAR (Bannister, 1974; Field and
Effler, 1983; and Stefan et aU, 1983). Because the temperature heat
balance computes total radiation over a wide spectrum, this value must be
adjusted to PAR if it is to be used in the algae simulation. The ratio
of energy in the visible band (PAR) to energy in the complete (standard)
spectrum is approximately 0.43 to 0.45 (Bannister, 1974 and Stefan et
al., 1983). TFACT is a user input variable; thus a value to meet site
specific conditions may be used.
Summary of Daily Averaging Options: The selection of a light averag-
ing option depends largely on the extent to which the user wishes to
account for the diurnal variation in light intensity. Options 1 and 2
use a single calculation of FL based on an "average" daily solar radiation
value. Options 3 and 4 calculate hourly values of FL from hourly values
of solar radiation and then average the hourly FL values to obtain the
daily average value. Options 1 and 3 use the solar radiation from the
temperature heat balance routines. (Thus both algae and temperature
simulations draw on the same source for solar radiation.) Options 2 and
31
-------
4 use the solar radiation value provided by the user for algae simulation.
Thus, either option 2 or 4 must be selected when algae are simulated and
temperature is not. The light averaging factor (AFACT) is used to provide
similarity in FL calculations between options 1 and 2 versus options 3
and 4. The solar radiation factor (TFACT) specifies the fraction of the
solar radiation computed in the heat balance, which is photosynthetically
active. It is used only with options 1 or 3.
In dynamic algae simulations, photosynthetically active radiation is
computed hourly using Equation III-9c unless temperature is not simulated,
in which case photosynthetically active solar radiation data must be
supplied with the local climatology data.
3.2.3.3 Algal Self Shading
The light extinction coefficient, x, in Equations III-6a,b,c is
coupled to the algal density using the nonlinear equation
X = XQ + A1 agA + X2(aQA)2/3 111-12
where
XQ = non-algal portion of the light extinction coefficient, ft"1
\l = linear algal self shading coefficient, ft"1 (ug-Chlai/L)"1
\2 - nonlinear algal self shading coefficient, ft"1 (ug-Chla/L)~2/3
ofj = conversion factor, ug-Chl ji /mg A
A = algal biomasis concentration, mg-A/L
Appropriate selection of the values of X^ and Xg allows modeling of
a variety of algal self-shading, light-extinction relationships:
t No algal self shading (QUAL-II SEMCOG)
\l = X2 = 0
a Linear algal self shading (Meta Systems QUAL-II)
Xi i 0 , X2 = 0
Nonlinear algal self shading (Riley Eq., in Zison et^ al_., 1978)
\! = 0.00268, ft"1 (ug-ChU/L)"1
X2 = 0.0165, ft'1 (ug-ChU/L)"2/3
32
-------
or
\l = 0.0088, m'1 (ug-ChU/L)'1
\2 = 0.054, m'1 (ug-ChU/L)"2/3
3.2.4 Algal Nutrient Relationships
The algal growth limitation factors for nitrogen (FN) and for phos-
phorus (FP) are defined by the Monod expressions:
Ne
FN = . 111-13
N +K
and
P2
FP = 111-14
P2 + Kp
where
Ne = the effective local concentration of available inorganic
nitrogen, mg-N/L
KN = the Michaelis-Menton half-saturation constant for nitrogen,
mg-N/L
?2 = the local concentration of dissolved phosphorus, mg-P/L
Kp = the Michaelis-Menton half-saturation constant for
phosphorus, mg-P/L
Algae are assumed to use ammonia and/or nitrate as a source of in-
organic nitrogen. The effective concentration of available nitrogen is
given by:
Ne = N! + N3 111-15
where
HI = concentration of ammonia nitrogen, mg-N/L
N3 = concentration of nitrate nitrogen, mg-N/L
The empirical half-saturation constants for nitrogen, K^, and phos-
phorus, Kp, are used to adjust the algal growth rate to account for those
33
-------
factors that can potentially 1imit algal growth. Each constant is actually
the level at which that particular factor limits algal growth to half the
maximal or "saturated" rate (Zison et jil_., 1978). Table III-3 at the end
of this chapter lists typical values of the half-saturation constants for
nitrogen and phosphorus. If algal concentrations are simulated and
either nitrogen, phorphorus, or both are not simulated, the program
assumes that the parameter not simulated is not limiting.
3.2.5 Temperature Dependence in Algae Simulation
The algal growth rate and death rates are temperature dependent.
They are corrected within the model, as are all other temperature dependent
systems variables, according to the procedure explained in Section 3.10.
3.3 NITROGEN CYCLE
In natural aerobic waters, there is a stepwise transformation from
organic nitrogen to ammonia, to nitrite, and finally to nitrate. The
nitrogen cycle in QUAL2E contains all four of these components, as shown
in Figure III-l. The incorporation of organic nitrogen as a state variable,
an organic nitrogen settling term, and an algal nitrogen uptake preference
factor are the primary enhancements to the nitrogen cycle in QUAL2E com-
pared to the SEMCOG version of QUAL-II. The differential equations
governing transformations of nitrogen from one form to another are shown
below.
3.3*1 Organic Nitrogen
dN4
- = 01 P A - 33 N4 - o4 N4 111-16
dt
where
N4 = concentration of organic nitrogen, mg-N/L
33 = rate constant for hydrolysis of organic nitrogen to
ammonia nitrogen, temperature dependent, day*
a} = fraction of algal biomass that is nitrogen, mg-N/mg-A
p = Algal respiration rate, day-1
A = algal biomass concentration, mg-A/L
04 = rate coefficient for organic nitrogen settling, temperature
dependent, dayl
34
-------
3.3.2 Ammonia Nitrogen
» = foNA - BiNi + °3/d - FI onjiA 111-17
dt
where
Fi = PNNI/(PNNI + (i ~ PN)NS) ni-18
NI = the concentration of ammonia nitrogen, mg-N/L
N3 = the concentration of nitrate nitrogen, mg-N/L
N4 = the concentration of organic nitrogen, mg-N/L
B! = rate constant for the biological oxidation of ammonia nitrogen,
temperature dependent, day*
B3 = organic nitrogen hydrolysis rate, day'1
ai = fraction of algal biomass which is nitrogen, mg-N/mg-A
-------
where
HI = the concentration of ammonia nitrogen, mg-N/L
N2 = the concentration of nitrite nitrogen, mg-N/L
3} = rate constant for the oxidation of ammonia nitrogen,
temperature dependent, dayl
32 = fate constant for the oxidation of nitrite nitrogen,
temperature dependent, dayl
3.3.4 Nitrate Nitrogen
dN3
- (1 - F)cqyA 111-20
dt
where
F = fraction of algal nitrogen taken from ammonia pool, as
defined in Section 3.3.2
cq = fraction of algal biomass that is nitrogen, mg-N/mg-A
y = local specific growth rate of algae, dayl
3.3.5 Inhibition of Nitrification at Low Dissolved Oxygen
QUAL2E has the capability of inhibiting (retarding) the rate of
nitrification at low values of dissolved oxygen. This inhibition effect
has been reported by others (Department of Scientific and Industrial
Research, 1964; Texas Water Development Board, 1984).
Nitrification rates are modified in QUAL2E by computing an inhibition
correction factor (having a value between zero and one) and then applying
this factor to the values of the nitrification rate coefficients, P]_, and
32- The nitrification rate correction factor is computed according to
a first order equation:
CORDO = 1.0 - EXP(-KNITRF * DO) 111-21
where
CORDO = nitrification rate correction factor
EXP = exponential function
36
-------
KNITRF = first order nitrification inhibition coefficient, mg/L-1
DO = dissolved oxygen concentration, mg/L
The correction factor is applied to the ammonia and nitrite oxida-
tion rates by:
Ammonia: (3i)1nhib. = CORDO * Wlnput HI-22
Nitrite: (32)inhib. = C0RDO * (^input in~23
A value of 0.6 for KNITRF closely matches the inhibition formula-
tion in QUAL-TX, the Texas Water Development Board version of QUAL-II,
whereas, a value of 0.7 closely simulates the data for the Thames Estuary
3.4 PHOSPHORUS CYCLE
The phosphorus cycle operates like the nitrogen cycle in many respects.
Organic forms of phosphorus are generated by the death of algae, which
then converts to the dissolved inorganic state, where it is available to
algae for primary production. Phosphorus discharged from sewage treatment
plants is generally in the dissolved inorganic form and is readily taken up
by algae (Zison et aU , 1978). QUAL2E revises the SEMCOG version of QUAL-
II, which included only dissolved phosphorus, to simulate the interactions
between organic and dissolved phosphorus. Below are the differential
equations governing transformations of phosphorus from one form to another.
3.4.1 Organic Phosphorus
= a? p A - 34?! - acPi 111-24
dt
where
P! = the concentration of organic phosphorus, mg-P/L
a£ = phosphorus content of algae, mg P/mg-A
p = algal respiration rate, day-1
A = algal biomass concentration, mg-A/L
34 = organic phosphorus decay rate, temperature dependent,
day-1
05 = organic phorphorus settling rate, temperature dependent ,
day-1
37
-------
3.4.2 Dissolved Phosphorus
dP2
« = 34?! + 02/d - «2yA 111-25
dt
where
?2 ~ concentration of inorganic or dissolved phosphorus, mg-P/L
02 = benthos source rate for dissolved phosphorus, temperature
dependent, mg-P/ft^-day
d = mean stream depth, ft
p = algal growth rate, dayl
A = algal biomass concentration, mg-A/L
3.5 CARBONACEOUS BOD
The QUAL2E model assumes a first order reaction to describe deoxygen-
ation of ultimate carbonaceous BOD in the stream. The BOD function as
expressed in the model also takes into account additional BOD removal due
to sedimentation, scour and flocculation, which do not exert an oxygen
demand (Thomas, 1948):
dL
= - iqi - K3L 111-26
dt
where
L = the concentration of ultimate carbonaceous BOD, mg/L
KI = deoxygenation rate coefficient, temperature dependent, day"-*-
l<3 = the rate of loss of carbonaceous BOD due to settling,
temperature dependent, dayl
QUAL2E simulates ultimate BOD in the general case; however, the user
may choose to use 5-day BOD values for input and output. In this case,
the model will make the necessary coversions from 5-day to ultimate BOD.
The conversion equation is:
BOD5 = BODU (1.0 - EXP(5 * KBOD)) 111-27
38
-------
where
BOD5 = 5-day BOD, mg/L
BODU = ultimate BOD, mg/L
KBOD = BOD conversion rate coefficient, day1
The SEMCOG version of QUAL-II uses a value of 0.23 day1 for KBOD.
With QUAL2E, the user may specify the appropriate value for this conver-
sion. Note: when modeling 5-day BOD, the conversion coefficient is
applied to all input BODs forcing functions (headwaters, incremental
flows, point loads, and the downstream boundary condition).
3.6 DISSOLVED OXYGEN
The oxygen balance in a stream system depends on the capacity of the
stream to reaerate itself. This capacity is a function of the advection
and diffusion processes occurring within the system and the internal
sources and sinks of oxygen. The major sources of oxygen, in addition to
atmospheric reaeration, are the oxygen produced by photosynthesis and the
oxygen contained in the incoming flow. The sinks of dissolved oxygen
include biochemical oxidation of carbonaceous and nitrogenous organic
matter, benthic oxygen demand and the oxygen utilized by respiration
(Zison et aj_., 1978).
The differential equation used in QUAL2E to describe the rate of
change of oxygen is shown below. Each term represents a major source or
sink of oxygen.
dO
= K2(0*-0) + (03 y - 04?) A - K! L - K4/d - 05 Bj Nj. - 05 &2 N2 IH-28
dt
where
0 = the concentration of dissolved oxygen, mg/L
0* = the saturation concentration of dissolved oxygen at the
local temperature and pressure, mg/L
»
03 = the rate of oxygen production per unit of algal photo-
synthesis, mg-0/mg-A
04 = the rate of oxygen uptake per unit of algae respired,
mg-0/mg-A
05 = the rate of oxygen uptake per unit of ammonia nitrogen
oxidation, mg.->0/mg-N
39
-------
«5 = the rate of oxygen uptake per unit of nitrite nitrogen
oxidation, mg-0/mg-N
v = algal growth rate, temperature dependent, day"1
p = algal respiration rate, temperature dependent, day"1
A = algal biomass concentration, mg-A/L
L = concentration of ultimate carbonaceous BOD, mg/L
d = mean stream depth, ft
K! = carbonaceous BOD deoxygenation rate, temperature dependent,
day1
l<2 = the reaeration rate in accordance with the Fickian
diffusion analogy, temperature dependent, day-1
l<4 = sediment oxygen demand rate, temperature dependent,
g/ft2-day
B! = ammonia oxidation rate coefficient, temperature dependent,
day1
&2 = nitrite oxidation rate coefficient, temperature dependent,
day1
NI = ammonia nitrogen concentration, mg-N/L
N£ = nitrite nitrogen concentration, mg-N/L
3.6.1 Dissolved Oxygen Saturation Concentration
The solubility of dissolved oxygen in water decreases with increasing
temperature, increasing dissolved solids concentration, and decreasing
atmospheric pressure (Zison et al_., 1978). QUAL2E uses a predictive
equation for the saturation "(equilibrium) concentration of dissolved
oxygen (APHA, 1985).
InO* = -139.34410 + (1.575701 x 1Q5/T) - (6.642308 x 10?/T2)
+ (1.243800 x 1Q10/T3) - (8.621949 x lO11/!^) 111-29
where:
0* = equilibrium oxygen concentration at 1.000 atm, mg/L
T = temperature (°K) = (°C+273.150) and °C is within the
range 0.0 to 40.0°C
40
-------
For non-standard conditions of pressure, the equilibrium concentra-
tion of dissolved oxygen is corrected by the equation 1 1 1-30:
(1-Pwv/P) (HP)
Op = 0*P [.. ..] 111-30
(1-Pwv) (1-*)
where
Op = equilibrium oxygen concentration at non-standard pressure,
mg/L
0* = equilibrium oxygen concentration at 1.000 atm, mg/L
P = pressure (atm) and is within 0.000 to 2.000 atm
Pwv = partial pressure of water vapor (atm), which may be
computed from:
= 11.8571 - (3840. 70/T) - 216961/T2) 111-31
and
$ = 0.000975 - (1.426 x 10-5t) + (6.436 x 10-8t2) 111-32
where
t = temperature, °C
The equations in Standard Methods (1985) for computing dissolved oxy-
gen saturation concentrations also include corrections for salinity and
chloride. Because neither salinity nor chloride is explicitly modeled,
QUAL2E does not correct 0* for chloride or salinity. Furthermore, the
pressure correction to 0* (Equation 1 11-30) is made only when temperature
is modeled, because barometric pressure data are a primary requirement of
the heat balance equations.
The dissolved oxygen saturation concentrations computed from the
Texas and SEMCOG versions of QUAL-II are compared to those from the
Standard Methods formulations of QUAL2E in Table III-l.
3.6.2 Atmospheric Reaeration Coefficient Estimation
The reaeration coefficient (Ko) is most often expressed as a function
of stream depth and velocity. QUAL2E provides eight options for estimating
or reading in K£ values, which are discussed in the sections below. A
comparative study of reaeration prediction equation performance has been
reported by St. John et al_. (1984).
41
-------
TABLE III-l
COMPARISON OF DISSOLVED OXYGEN SATURATION CONCENTRATIONS
(Barometric Pressure = 1 atm, Chloride = 0.0 mg/L,
Equilibrium with Air Saturated with Water Vapor)
Temperature,
°C
0.
1.
2.
3.
4.
5.
6.
7.
8.
9.
10.
11.
12.
13.
14.
15.
16.
17.
18.
19.
20.
21.
22.
23.
24.
25.
26.
27.
28.
29.
30.
31.
32.
33.
34.
35.
36.
37.
38.
39.
40.
QUAL-II
SEMCOG
14.632
14.227
13.837
13.461
13.100
12.752
12.418
12.096
11.787
11.489
11.203
10.927
10.661
10.406
10.159
9.922
9.692
9.471
9.257
9.050
8.849
8.655
8.465
8.281
8.101
7.925
7.753
7.584
7.417
7.252
7.089
6.927
6.765
6.604
6.442
6.280
6.116
5.950
5.782
5.612
5.438
QUAL-TX
Texas
14.584
14.187
13.806
13.441
13.091
12.755
12.433
12.124
11.828
11.544
11.271
11.009
10.758
10.517
10.285
10.062
9.848
9.642
9.444
9.253
9.069
8.891
8.720
8.555
8.396
8.241
8.092
7.948
7.807
7.672
7.540
7.412
7.288
7.167
7.049
6.935
6.823
6.715
6.609
6.506
6.406
QUAL2E
Std. Meth.
14.621
14.217
13.830
13.461
13.108
12.771
12.448
12.139
11.843
11.560
11.288
11.027
10.777
10.537
10.306
10.084
9.870
9.665
9.467
9.276
9.093
8.915
8.744
8.578
8.418
8.264
8.114
7.969
7.828
7.691
7.559
7.430
7.305
7.183
7.065
6.949
6.837
6.727
6.620
6.515
6.413
42
-------
.K2 Option 1
Option 1 allows the user to read in K2 values that have been pre-
viously selected by the modeler. This option is useful in modeling
unusual situations such as ice cover (see Section 3.6.3).
K? Option 2
Using data collected in field measurements of stream reaeration,
Churchill, Elmore, and Buckingham (1962) developed the following expression
for K2 at 20°C.
K220 = 5.026 IT °-969 d-1-673 x 2.31 111-33
where
u = average velocity in the stream, ft/sec.
d = average depth of the stream, ft
K2 = reaeration coefficient, day"1
Kp Option 3
O'Connor and Dobbins (1958) proposed equations based on the turbulence
characteristics of a stream. For streams displaying low velocities and
isotropic conditions, Equation 111-34 was developed:
20 = (Dm ") ' | m_34
* dl.50
For streams with high velocities and nonisotropic conditions, the rela-
tionship is:
4800^*5 S^'25
Ko20 = I ..." ,. x 2.31 111-35
L dl-25
where
S0 = slope of the streambed, ft/ft
d = mean stream depth, ft
~\\ = mean velocity, ft/day
K2 = reaeration coefficient, day"1
43
-------
and Dm is the molecular diffusion coefficient (ft2/day), which is given
by:
Dm = 1.91 x 103 (1.037)1"-20 111-36
Equation 111-34 has been found to be generally applicable for most cases
and is the equation used in QUAL2E for Option 3. Equation II1-35 can
be used to calculate l<2 outside the model and input it directly under
Option 1.
K? Option 4
Based on the monitoring of six streams in England, Owens et al. (1964)
developed an equation for shallow, fast moving streams. Equation 111-37
can be used with streams that exhibit depths of 0.4 to 11.0 feet and velo-
cities of 0.1 to 5.0 ft/sec:
K220 = 9.4 Tj0.67/d1.85 x 2.31 111-37
where
u = mean velocity, ft/sec
d = mean depth, ft
K? Option 5 .
Thackston and Krenkel (1966) proposed the following equation based
on their investigation of several rivers in the Tennessee Valley Authority
system.
u*
K2,0 = 10.8 (1 + F0-5) «. x 2.31 111-38
d
where F is the Froude number, which is given by:
u*
F = . 111-39
and u* is the shear velocity, ft/sec.:
u n / g
u* = / d Seg = i 111-40
6 1.49 dl-167
44
-------
where
d = mean depth, ft
g = acceleration of gravity, ft/sec2
Se = slope of the energy gradient
u = mean velocity, ft/sec
n = Manning's coefficient
K Option 6
Langbien and Durum (1967) developed a formula for l<2 at 20°C:
K220 = 3.3 "u/d1-33 x 2.31 111-41
where
IF = mean velocity, ft/sec
d = mean depth, ft
K;? Option 7
This option computes the reaeration coefficient from a power function
of flow. This empirical relationship is similar to the velocity and depth
correlations with flow used in the hydraulics section of QUAL2E, i.e.,
K2 = aQb 1 1 1-42
where
a = coefficient of flow for K2
Q = flow, ft3/sec
b = exponent on fl ow for K2
Kg Option 8
The method of Tsivoglou and Wallace (1972) assumes that the reaera-
tion coefficient for a reach is proportional to the change in elevation
of the water surface in the reach and inversely proportional to the flow
time through the reach. The equation is:
45
-------
20
K220 = c 111-43
tf
where
c = escape coefficient, ft-1
Ah = change in water surface elevation in reach, ft
tf = flow time within reach, days
Assuming uniform flow, the change in water surface elevation is
Ah = Se AX 111-44
where
Se = slope of the energy gradient, ft/ft
AX = reach length, ft
and the time of passage through a reach is
AX
tf = «-- 111-45
u
where
IT = mean velocity in reach, ft/sec
Substituting the above in equation 111-43 gives
K220 = (3600 x 24) cSe u 111-46
Equation 111-46 is the form of Option 8 used in QUAL2E. The constants
3600 and 24 convert velocity to units of feet per day. The slope may be
input directly for computing K2 with this option, or it can be calculated
by:
_2
u n2
Se = ___ 111-47
(1.49)2 d4/3
46
-------
where
d = mean depth, ft
n = Manning's coefficient
The escape coefficient is usually treated as a variable and determined
empirically. TenEch (1978) recommends the following guideline in deter-
mining c values, analogous to that recommended for uncalibrated stream
segments by Tsivoglou and Neal (1976):
c = 0.054 ft-1 (at 20°C) for 15 _< Q <_ 3000 ft3/sec
c = 0.110 ft-1 (at 20°C) for 1 ^ Q ^ 15 ft3/sec
3.6.3 Ice Cover
Ice cover on streams during winter low flow conditions may signifi-
cantly affect reaeration. Reaeration rates are decreased because ice
cover reduces the surface area of the air-water interface through which
reaeration occurs (TenEch, 1978). Approaches recommended by TenEch
(1978) for estimating the extent of ice cover include:
o Statistical analyses of past records
o Steady state heat budget analysis (including the U.S. Army
Corps of Engineers differential equations)
o Extensive field observations
To adjust the reaeration rate for winter ice cover conditions in the
QUAL2E model, the calculated reaeration rate must be multiplied by an "ice
cover factor" and input under Option 1. TenEch recommends factors ranging
from 0.05 for complete ice cover to 1.0 for no ice cover. Depending on
the extent of cover, reaeration values can be greatly reduced.
3.6.4 K? Default Values
There are no default Kg values in QUAL2E. In some versions of QUAL-
II, a default value of K£ is computed, accounting for the influences of
wind-induced turbulence and diffusion under low-velocity conditions. In
those models, when the calculated values of Kg are less than two divided
by the depth of the reach (2/d), Kg is set equal to 2/d. This feature
has not always proved useful, particularly when simulating the very low
reaeration rates; thus it is not included in QUAL2E.
47
-------
3.6.5 Dam Reaeration
QUAL2E has the capability of modeling oxygen input to the system from
reaeration over dams. The following equation described by Zison et a1.
(1978) and attributable to Gameson is used to estimate oxygen input from
dam reaeration.
Da - Db = [1 - ] Da 111-48
1 + O.llab(l + 0.046T)H
where
Da = oxygen deficit above dam, mg/L
£>b = oxygen deficit below dam, mg/L
T = temperature, °C
H = height through which water falls, ft
a = empirical parameter
= 1.25 in clear to slightly polluted water
= 1.0 in polluted water
b = empirical parameter
. = 1.0 for weir with free fall
= 1.3 for step weir or cascades
The factors H, a and b are input for each dam. The model includes a
provision for bypassing some or all of the flow around the dams (e.g.,
through generators). The fraction of the total flow that spills over the
dam is supplied as an input variable.
3.7 COLIFORMS
Coliforms are used as an indicator of pathogen contamination in sur-
face waters. Expressions for estimating coliform concentrations are
usually first order decay functions, which only take into account coliform
die-off (Zison et aj_., 1978). The QUAL2E model uses such an expression:
dE
= - K5 E 111-49
dt b
48
-------
where
E = concentration of coliforms, colonies/100 ml
K5 = coliform die-off rate, temperature dependent, day"1
3.8 ARBITRARY NONCONSERVATIVE CONSTITUENT
QUAL2E has the provision for modeling an arbitrary nonconservative
constituent (ANC). In addition to a first order decay mechanism, there
are source and sink terms in the mass balance. The differential equation
describing the interactions for an arbitrary nonconservative constituent
is:
dR
= -K6 R - o6R + a7/d 111-50
dt
where
R = concentration of the nonconservative constituent, mg-ANC/L
Kg = decay rate for the constituent, temperature dependent, day"1
ag = rate coefficient for constituent settling, temperature
dependent, day-1
07 = benthal source for constituent, temperature dependent,
mg-ANC/ft2-day
d = mean stream depth, ft
3.9 TEMPERATURE
Temperature is modeled by performing a heat balance on each computa-
tional element in the system. The heat balance accounts for temperature
inputs and losses from the forcing functions as well as the heat exchanged
between the water surface and the atmosphere. The air-water heat balance
terms include long and short wave radiation, convection, and evaporation
using:
Hn = Hsn + Han - Hb - Hc - He 111-51
where
Hn = net heat flux passing the air water surface, Btu/ft^-day
49
-------
Hsn = net short wave solar radiation after losses from
absorption and scattering in the atmosphere and by
reflection at the interface, Btu/ft2-day
Han = net long wave atmosphere radiation after reflection,
Btu/ft2-day
Hb = outgoing long wave back radiation, Btu/ft2-day
HC = convective heat flux, Btu/ft2-day
He = heat loss by evaporation, excluding sensible heat loss,
Btu/ft2-day
In order for QUAL2E to perform the heat balance computations, the
user must supply a variety of data, including the longitude and latitude
of the basin, the time of year, evaporation coefficients, and a dust
attenuation coefficient. Local climatological information in the form of
time of day, wet and dry bulb air temperatures, atmospheric pressure,
cloud cover and wind velocity also must be provided. These data are
applied uniformly over the entire river basin.
In the dynamic mode, local climatological data must be supplied at
regular (typically 3 hour) intervals. In this manner the source/sink
term for the heat balance is updated in time to simulate the diurnal
response of the steady hydraulic system to changing temperature condi-
tions.
In the steady state mode, average local climatological data must
be supplied by the user. The program uses linear approximations for the
longwave back radiation and evaporation terms for solution of the steady
state heat balance. The reader is referred to Chapter 4 of this report
for a detailed treatment of the temperature simulation.
3.10 TEMPERATURE DEPENDENCE OF RATE COEFFICIENTS
The temperature values computed in QUAL2E are used to correct the rate
coefficients in the source/sink terms for the other water quality variables.
These coefficients are input at 20°C and are then corrected to temperature
using a Streeter-Phelps type formulation:
XT = X20 0 (T-20°) 111-52
where
Xj = the value of the coefficient at the local temperature (T)
X2Q = the value of the coefficient at the standard temperature
(20°C)
50
-------
TABLE 111-2
DEFAULT TEMPERATURE CORRECTION, 9, VALUES FOR QUAL2E
Rate Coefficient
BOD Decay
BOD Settling
Reaeration
SOD Uptake
Organic N Decay
Organic N Settling
Ammoni a Decay
Ammonia Source
Nitrite Decay
Organic P Decay
Organic P Settling
Dissol ved P Source
Algal Growth
Algal Respiration
Algal Settling
Coli form Decay
ANC
ANC
ANC
Symbol
K!
KS
K2
K4
S3
°4
31
°3
^2
34
°5
°2
y
P
°1
KS
<6
°6
°7
Default
SEMCOG
1.047
-
1.0159
-
-
-
1.047
-
1.047
-
-
-
1.047
1.047
-
1.047
1.047
-
-
Values
QUAL2E
1.047
1.024
1.024
1.060
1.047
1.024
1.083
1.074
1.047
1.047
1.024
1.074
1.047
1.047
1.024
1.047
1.000
1.024
1.000
Note: - = not temperature dependent in QUAL-II SEMCOG.
ANC = Arbitrary Nonconservative Constituent
51
-------
0 = an empirical constant for each reaction coefficient
The values of the temperature correction factors, 9, may be specified by
the user. In the absence of user specified values, the default values
shown in Table III-2 are employed. For comparison purposes, the 9 values
used in the SEMCOG version of QUAL-II are also listed in Table III-2.
If temperature is not simulated, the temperature value specified for
the initial condition is assumed to be the temperature for the simulation.
3.11 REACTION RATES AND PHYSICAL CONSTANTS
The chemical and biological reations that are simulated by QUAL2E
are represented by a complex set of equations that contain many system
parameters; some are constant, some are spatially variable, and some are
temperaturedependent. Table 111-3 lists these system parameters and
gives the usual range of values, units, and types of variation. Kramer
(1970), Chen and Orlob (1972), and Zison et al_. (1978) give detailed
discussions of the basic sources of data, ranges and reliabilities of
each of these parameters. Final selection of the values for many of
these system parameters or measurement of sensitive ones should be made
during model calibration and verification.
TABLE 111-3
TYPICAL RANGES FOR QUAL2E REACTION COEFFICIENTS
Range of Variable Temperature
Variable Description Units Values by Reach Dependent
ao
al
«2
a3
04
a5
Ratio of chlorophyl 1-a
to algal biomass
Fraction of algal biomass
that is Nitrogen
Fraction of algal biomass
that is Phosphorus
02 production per unit of
algal growth
02 uptake per unit of
algae respired
02 uptake per unit of
NH3 oxidation
ug-Chl a
mg A
mg-N
mg A
mg-P
mg A
mg-Q
mg A
mg-0
mg A
mg-0
mg N
10-100
0.07-0.09
0.01-0.02
1.4-1.8
1.6-2.3
3.0-4.0
No
No
No
No
No
No
No
No
No
No
No
No
52
-------
TABLE III-3 (cont'd)
TYPICAL RANGES FOR QUAL2E REACTION COEFFICIENTS
Variable Description
«6
^max
P
KL
KN
Kp
x°
M
*2
PN
°i
CT2
°3
04
°5
Q£ uptake per unit of
N02 oxidation
Maximum algal growth rate
Algal respiration rate
Michaelis-Menton half-
saturation constant
for light (Option 1)
Michaelis-Mention half-
saturation constant
for nitrogen
Michaelis-Menton half-
saturation constant
for phosphorus
Non-algal light extinc-
tion coefficient
Linear algal self-shading
coefficient
Nonlinear algal self-
shading coefficient
Algal preference factor
for ammonia
Algal settling rate
Benthos source rate for
dissolved phosphorus
Benthos source rate for
ammonia nitrogen
Organic nitrogen
settling rate
Organic phosphorus
settling rate
Units
mg-0
mg N
day'1
day1
Btu/ft2-
min
mg-N/L
mg-P/L
ft'1
I/ft
ug Chl-a/L
I/ft
Range of Variable Temperature
Values by Reach Dependent
1.0-1.14
1.0-3.0
0.05-0.5
0.02-0.10
0.01-0.20
0.01-0.05
Variable
0.002-0.02
0.0165
No
No
No
No
No
No
No
No
No
No
No
No
No
No
No
No
No
No
(ug Chl-a/L)Z/J (Riley)
-
ft/day
mg-P
ft^-day
mg-0
ft^-day
day'1
day'1
53
0.0-1.0
0.5-6.0
Variable
Variable
0.001-0.10
0.001-0.10
No
Yes
Yes
Yes
Yes
Yes
No
Yes
Yes
Yes
Yes
Yes
-------
TABLE III-3 (cont'd)
TYPICAL RANGES FOR QUAL2E REACTION COEFFICIENTS
Variabl
°6
°7
*1
K2
K3
K4
K5
<6
Pi
e Description
Arbitrary non-conserva-
tive settling rate
Benthal source rate for
arbitrary non-conserva-
tive settling rate
Carbonaceous deoxygenera-
tion rate constant
Reaeration rate constant
Rate of loss of BOD due
to settling
Benthic oxygen uptake
Col i form die-off rate
Arbitrary non-conserva-
tive decay coefficient
Rate constant for the
Units
day'1
mq-ANC
ft2-day
day"1
day'1
day'1
m-0
ft* -day
day -1
day'1
day -1
Range of Variable Temperature
Values by Reach Dependent
Variable
Variable
0.02-3.4
0.0-100
-0.36-0.36
Variable
0.05-4.0
Variable
0.10-1.00
Yes
Yes
Yes
Yes
Yes
Yes
Yes
Yes
Yes
Yes
Yes
Yes
Yes
Yes
Yes
Yes
Yes
Yes
B4
biological oxidation
of NH3 to N02
Rate constant for the
biological oxidation
of N02 to N03
Rate constant for the
hydrolysis of organic-
N to ammonia
Rate constant for the
decay of organic-P
to dissolved-P
day
'1
0.20-2.0
Yes
day-1 0.02-0.4 Yes
day1 0.01-0.7 Yes
Yes
Yes
Yes
54
-------
4. FUNCTIONAL REPRESENTATION OF TEMPERATURE
4.1 BASIC TEMPERATURE EQUATION
The basic mass transport equation for QUAL2E was given in Section
II as (see equation II-3):
_§C
8C 3(AXDL 3x) 3(AX u C) dC s
-- = , - , + - + - iv-l
3t Ax 3x Ax 3x dt V
In temperature modeling, C is taken as the concentration of heat (HL~3)
and can be equated to temperature through the relationship
C = p c (T - T0) IV-2
where
p = the density of water (M L~3)
c = the heat capacity of water (HM~1 D-1)
T = the water temperature
T0 = an arbitrary base temperature
M = mass
H = heat energy flux
D = degrees
The parameters p and c can be considered constant for practical purposes.
Also, the internal heat generation dC, which results from viscous dissi-
dt
pation of energy and boundary friction, is generally small enough to be
55
-------
considered negligible. Thus setting ^C = 0 in equation IV-1 and substituting
dt
equation IV-2 for C gives us (after some simplification):
*L
3T 3(AXDL 3x) 9(AX u T) 1 s
= _ - 'i + - IV-3
3t Ax 9x Ax 3x pC V
The source term s/V (with units of HL-3T-1) accounts for all heat
transferred across the system boundaries, i.e., heat transferred across
the air-water interface and heat conducted across mud-water interface.
Heat transfer across the mud-water interface is generally insignificant;
hence, s/V takes on the identify of the net rate of heat input per unit
volume of stream through the air-water interface.
It is most convenient to represent the interfacial heat transfer
rate as a flux (HN) having units of HL'M"1. For a stream element of
length dx and mean surface width W, H|\j is related to s/V as follows.
The total rate of heat input across the air-water interface is H^ dx
W. _This heat isjdistributed uniformly throughout the underlying volume
of Ax dx, where Ax is the mean cross-sectional area of the element. Thus
the rate of heat gain per unit volume of water, s/V, is computed as:
s s HN (Wdx) HN
V Ax dx Ax dx d
where d = AX/W is the hydraulic depth of the stream. Substituting equation
IV-4 into equation IV-3 gives the generalized form of the temperature
equation:
JL
9T 3(AXD[_ 3x) 8(AX u T) HN
_ = , _ | + «_. iv-5
3t Ax 3x Ax 3x pcd
4.2 DEFINITION OF HN
Heat is transferred across the air-water interface of a surface
water body by three difference processes: radiation exchange, evaporation,
and conduction. The individual heat terms associated with these processes
are shown in Figure IV-1 and are defined in Table IV-1 with the typical
ranges of their magnitudes in northern latitudes also listed.
56
-------
The expression that results from the sumnation of these various
energy fluxes is:
where
HN
H
sn
HN = Hsn + Han - (Hb ± Hc + He)
IV-6
net energy flux passing the air-water interface,
Btu/ft2-day
net short-wave solar radiation flux passing
through the interface after losses due to
absorption and scattering in the atmosphere
and by reflection at the interface, Btu/ft2-day
net long-wave atmospheric radiation flux passing
through the interface after reflection, Btu/ft^-day
outgoing long-wave back radiation flux, Btu/ft2-day
convective energy flux passing back and forth
between the interface and the atmosphere, Btu/ft2-day
energy loss by evaporation, Btu/ft2-day
These mechanisms by which heat is exchanged between the water surface and
the atmosphere are fairly well understood and are adequately documented in
the literature by Edinger and Geyer (1965). The functional representation
H.
Hb
4
4
nc
AIR-WATER
'INTERFACE
Figure IV-1.
Heat Transfer Terms Associated with
Interfacial Heat Transfer
57
-------
TABLE IV-1
DEFINITION OF HEAT TRANSFER TERMS
ILLUSTRATED IN FIGURE IV-1
Heat Term Units Magnitude
(BTU/ft2-dayl)
Hs = total incoming solar or HL^T'1 400-2800
short-wave radiation
H$r = reflected short-wave radiation HL"2T'1 40-200
Ha = total incoming atmospheric HL'2?"1 2400-3200
ratiation
Har = reflected atmospheric radiation HL"2^1 70-120
Hb = back radiation from the water HL"2^1 2400-3600
surface
He = heat loss by evaporation HL^T"1 150-3000
Hc = heat loss by conduction to HL^r1 -320 to +400
atmosphere
of these terms has been defined by Water Resources Engineers, Inc. (1967).
The formulations reported here were extracted from that more detailed work
by Frank D. Masch and Associates and the Texas Water Development Board
(1971).
4.3 NET SHORT-WAVE SOLAR RADIATION
The net incoming solar radiation is short-wave radiation which passes
directly from the sun to the earth's surface. Its magnitude depends on:
the altitude of the sun, which varies daily as well as seasonally for a
fixed location on the earth; -the dampening effect of scattering and
absorption in the atmosphere due to cloud cover, and the reflection from
the water surface.
58
-------
The net amount of solar radiation which reaches the surface of the
reach may be represented functionally on an hourly basis by:
Hsn = ^Jo jj^ (^- RSJ (1 - 0.6JSC) IV-7
(i) (ii) (111) (W)~~
where
Hsp = net short-wave solar radiation flux, Btu/ft2-hr
H0 = amount of radiation flux reaching the earth's
atmosphere, Btu/ft2-hr
at = atmospheric transmission term
Rs = Albedo or reflection coefficient
CL = cloudiness as a fraction of sky covered
It is appropriate for purposes of this discussion to identify and
treat separately the four components in equation IV-7 as (i) extra-
terrestrial solar radiation, (ii) radiation scattering and absorption,
(iii) reflectivity, and (iv) cloudiness.
4.3.1 Extraterrestrial Radiation
The short-wave solar radiation flux that strikes the earth's outer
atmosphere over a given period of time is given by Water Resources
Engineers, Inc. (1967) as:
Hsc ^
HQ = ,- { sin sin 6 (te -
r2 180
12 ir<|> Trte
M. cos * cos 6 [sin () - sin ()]} r IV-8
IT 180 12 12
where
Hsc = solar constant = 438.0 Btu/ft2-hr
r = normalized radius of the earth's orbit
$ = latitude of the site, degrees
59
-------
>te = hour angles corresponding to the beginning and end,
respectively, of any time interval between sunrise
and sunset
r = a correction factor for diurnal exposure to
radiation flux
Listed below are several parameters in equation IV-8 requiring
further definition as described by Water Resources Engineers, Inc. (1967).
a. Relative Earth-Sun Distance--
2*
r = 1.0 + 0.17 cos [ (186-Dy)] IV-9
365
where Dy is the number of the day of the year (beginning January 1)
b. Declination--
23.45 2ir
6 = 1 TT cos [ (173-Dy)] IV-10
180 365
c, Hour Angles--
tb = STb - Ats + ET - 12 IV-11
and
te = STe - Ats + ET - 12 IV-12
where ST|j, STe are the standard times at the beginning and end of the
time interval selected
ET = an expression for time from a solar ephemeris that
represents the difference in hours between "true
solar time" and that computed on the basis of a yearly
average. It is given for each day of the year, Dy, by
60
-------
2lT
ET = 0.000121 - 0.12319 sin [ (Dy-1) - 0.0714]
365
4*
0.16549 sin [--- (Dy-1) + 0.3088] IV-13
365
Ats = difference between standard and local civil time
in hours as determined from:
e
Ats = ~ (LSm - Mm) IV-14
J.
-------
in which a" is the mean atmospheric transmission coefficient after
scattering and absorption, given by:
a" = exp { - [0.465 + 0.0408 Pwc]
[0.179 + 0.421 exp (-0.721 9am)] 9am} IV-20
where 9am is the optical air mass given by the expression:
exp (-Z/2531)
9am = IV-21
sin a + 0.15 (180a + 3.885)-l-253
IT
in which
Z = elevation of the site in feet
a = sun's altitude in radians, given by:
1T(j) TT(J>
a = arc sin [sin * sin 6 + cos <^
180 180
Trt
cos <5 cos --] IV-22
12
in which t is the hour angle, described by an equation similar to equation
IV-11 and IV-12.
Pwc in equation IV-20 is the mean daily precipitable water content
in the atmosphere, given by the expression:
PWC = 0.00614 exp (0.0489Td) IV-23
where T^ is the dewpoint in °F, which can be obtained from the expression:
Td = In [(ea + 0.0837)/0.1001]/0.03 IV-24
where ea is the water vapor pressure of the air.
The mean atmospheric coefficient, a1, can also be represented by an
equation of the form of equation IV-20 as:
62
-------
a1 = exp { - [0.465 + 0.0408 Pwc]
[0.129 + 0.171 exp (-0.880 Oam)] 0am} IV-25
Dust attenuation of the solar radiation flux, which is represented
in equation IV-19 by the quantity d, varies with optical air mass, season
of the year, and geographic location. Water Resources Engineers, Inc.
(1967) gives a range of 0-0.13 for several locations.
4.3.3 Cloudiness
The dampening effect on the solar radiation flux is given by Water
Resources Engineers, Inc. (1967) as
Cs = 1.0 - 0.65 c£ IV-26
where C|_ is the decimal fraction of the sky covered. Water Resources
Engineers, Inc. (1967) reports that equation IV-26 gives satisfactory
results except for heavy overcast conditions, i.e., when C|_ approaches
1.0.
4.3.4 Reflectivity
The reflection coefficient, Rs, can be approximately computed as a
function of the solar altitude, a, by Anderson's (1954) empirical formula:
Rs = AaB IV-27
where a is in degrees, and A and B are functions of cloudiness, CL.
Values for A and B given by Anderson (1954) are shown in Table IV-2.
TABLE IV-2
EMPIRICAL COEFFICIENTS FOR DETERMINING R<
After Anderson (1954)
Cloudiness
CL
0
Clear
0.1 - 0.5
Scattered
0.6 - 0.9
Broken
1.0
Overcast
Coefficients A BAB AB AB
1.18 -0.77 2.20 -0.97 0.95 -0.75 0.35 -0.45
63
-------
4.4 LONG-WAVE ATMOSPHERIC RADIATION
The long-wave radiation emitted by the atmosphere varies directly
with the moisture content of the atmosphere. Although it is primarily
dependent on air temperature and humidity, it can also be affected by
ozone, carbon dioxide, and possibly other materials in the atmosphere.
Anderson (1954) indicated that the amount of atmospheric radiation is
also significantly affected by cloud height. The amount of long-wave
atmospheric radiation that is reflected is approximately a constant
fraction of the incoming radiation, found by Anderson (1954) to be
approximately 0.03.
The net atmospheric radiation flux can be expressed as:
Han = [2.89 x 10"6] o (Tfl + 460)6 (1.0 + 0.17Cf)(l-RL) IV-28
where
Han = net long-wave atmospheric radiation flux, Btu/ft^-hr
o = Stefan-Boltzman constant, 1.73 x 10-9 Btu/ft2/hr/
°Rankine4
Ta = air temperature at a level 6 feet above the water
surface, °F
RL = reflectivity of the water surface for atmospheric
radiation = 0.03
C[ = cloudiness, fraction of cloud cover
4.5 WATER SURFACE BACK RADIATION
The third source of radiation transfer through the air-water interface
is long-wave back radiation from the water surface, H^, which represents
a loss of heat from the water. It can be seen from Table IV-1 that back
radiation accounts for a substantial portion of the heat loss from a body
of water. This loss is expressed by the Stefan-Boltzman Fourth Power
Radiation Law for a blackbody as:
Hb = 0.97 a (Ts + 460)4 IV-29
where
Hb = water surface back radiation flux, Btu/ft2-hr
Ts = water surface temperature, °F
64
-------
Equation IV-29 can be linearized over a given temperature range as
Hb =
-------
W = wind speed, in mph, measured 6 feet above the water
surface
es = saturation vapor pressure of the air, in. of Hg,
at the temperature of the water surface, as
given by
es = 0.1001 exp (0.03 Ts) - 0.0837
and
ea = water vapor pressure, in. of Hg, at a height of
6 feet above the water surface, given as
ea = ewb - 0.000367 Pa (Ta - Twb)
IV-34
1571
where
saturation vapor pressure, in. of Hg, at the
wet bulb temperature from the expression
ewb = 0.1001 exp (0.03 T^) - 0.0837 IV-35
Pa = local barometric pressure, in. of Hg
Twb = wet bulb air temperature, °F
Ta = dry bulb air temperature, °F
The literature contains a wide range of values for the evaporation
constants a and b. Roesner (1969) reports that a good average value of a
would be 6.8 x 10-4 ft/hr-in. of Hg, while b would best be represented by
2.7 x l.O'4 ft/hr-in. of Hg.-mph.
To linearize the variation of evaporation rate with surface water
temperature Ts, equation IV-34 is approximated over 5°F intervals as:
es = 01 + BI Ts IV-36
Sets of 01, Bi are specified for twenty-one 5°F intervals between 35°F
and 135°F. The linearized evaporation expression is used in the steady-
state temperature solution.
66
-------
The sensible evaporative heat loss can be expressed simply as:
Hv = c y E (Ts - T0) IV-37
where
c = heat capacity of water = 1 Btu/lb-°F
T0 = reference temperature, °F
Sensible heat loss is very small compared to the other heat loss components
in the energy budget and thus is not included in the QUAL2E temperature
computation.
4.7 CONDUCTION
Heat that is transferred between the water and the atmosphere due to
a temperature difference between the two phases is normally called
conduction. .Using the fact that transfer by conduction is a function of
the same variables as evaporation, it is possible to arrive at a propor-
tionality between heat conduction and heat loss by evaporation. This
proportionality, known as Bowen's ratio, is expressed as:
B = - = CB [ $ «] ~ IV-38
He es - ea 29.92
where Cg is a coefficient = 0.01.
By using Bowen's ratio, the rate of heat loss to the atmosphere by
heat conduction, Hc, can be defined as:
Pa
Hc = y HL (a+bW) (0.01 ) (Ts - Ta) IV-39
.29.92
For practical purposes, the ratio (Pa/29.92) can be taken as unity.
67
-------
5. COMPUTATIONAL REPRESENTATION
5.1 PROTOTYPE REPRESENTATION
To expand upon the basic conceptual representation presented in Sec-
tions 1 and 2, QUAL2E permits any branching, one-dimensional stream system
to be simulated. The first step involved in approximating the prototype
is to subdivide the stream system into reaches, which are stretches of
stream that have uniform hydraulic characteristics. Each reach is then
divided into computational elements of equal length so that all computa-
tional elements in all reaches are the same length. Thus, all reaches
must consist of an integer number of computational elements.
There are seven different types of computational elements:
1. Headwater element
2. Standard element
3. Element just upstream from a junction
4. Junction element
5. Last element in system
6. Input element
7. Withdrawal element
Headwater elements begin every tributary as well as the main river system,
and as such, they must always be the first element in a headwater reach.
A standard element is one that does not qualify as one of the remaining
six element types. Because incremental flow is permitted in all element
types, the only input permitted in a standard element is incremental
flow. A type 3 element is used to designate an element on the mainstem
that is just upstream of a junction. A junction element (type 4), has a
simulated tributary entering it. Element type 5 identifies the last com-
putational element in the river system (downstream boundary); there
should be only one element type 5. Element types 6 and 7 represent
elements which have inputs (waste loads and unsimulated tributaries) and
water withdrawals, respectively.
68
-------
River reaches, which are aggregates of computational elements, are
the basis of most data input. Hydraulic data, reaction rate coefficients,
initial conditions, and incremental flow data are constant for all computa-
tional elements within a reach.
5.2 FORCING FUNCTIONS
Forcing functions are the user specified inputs that drive the
system being modeled. These inputs are specified in terms of flow, water
quality characteristics, and local climatology. QUAL2E accommodates four
types of hydraulic and mass load forcing functions in addition to local
climatological factorsheadwater inputs, point sources or withdrawals,
incremental inflow/outflow along a reach, and the (optional) downstream
boundary concentration.
1. Headwater Inputs - Headwater inputs are typically the upstream
boundary conditions at the beginning of the system. They are the condi-
tions required to generate the solution of the mass balance equations for
the first computational element in each headwater reach. Headwaters are
also the source of water for flow augmentation.
2. Point Sources and/or Withdrawals - These loads are used to repre-
sent point source discharges into the system (i.e., sewage and industrial
waste, or storm water runoff) and losses from the system resulting from
diversions. In QUAL2E point source discharges may represent either raw
or treated waste loads. If raw waste loads are used, the effect of
treatment can be simulated by applying a specific fractional removal for
carbonaceous BOD to each point source load.
3. Incremental Inflow - QUAL2E has the capability to handle flow
uniformly added or removed along a reach. The total flow increment along
a reach is apportioned equally to all computational elements in the
reach. This feature can be used to simulate the effects of non-point
source inputs to the system, or the effect of loss of stream flow to the
groundwater.
4. Downstream Boundary Concentration (optional) - QUAL2E has the
capability of incorporating known downstream boundary concentrations of
the water quality constituents into the solution algorithm. This feature
is useful in modeling systems with large dispersion in the lower reaches
(e.g., estuaries). -When downstream boundary concentrations are supplied,
the solution generated by QUAL2E will be constrained by this boundary
condition. If the concentrations are not provided, the constituent
concentrations in the most downstream element will be computed in the
normal fashion using the zero gradient assumption (see Section 5.4.3).
Local climatological data are required for the dynamic simulation of
algae and temperature. The temperature simulation uses a heat balance
across the air-water interface and thus requires values of wet and dry
bulb air temperatures, atmospheric pressure, wind velocity, and cloud
69
-------
cover. The algal simulation requires values of net solar radiation.
These climatological data must be input at regular time intervals over
the course of the simulation. For modeling steady-state temperature and
algae, average daily local climatological data are required. All clima-
data are applied uniformly over the entire river basin.
5.3 MODEL LIMITATIONS
QUAL2E has been developed to be a relatively general program; however,
certain dimensional limitations have been imposed upon it during program
development. These limitations are as follows:
Reaches: a maximum of 50
Computational elements: no more than 20 per reach or 500 in total
Headwater elements: a maximum of 10
Junction elements: a maximum of 9
Input and withdrawal elements: a maximum of 50 in total
(Note: These limitations may be modified, if necessary, by the user by
altering the PARAMETER statement specifications in the common blocks of
the program.)
QUAL2E can be used to simulate any combination of the following
parameters or groups of parameters.
1. Conservative minerals (up. to three at a time)
2. Temperature
3. BOD
4. Chlorophyll a^
5. Phosphorus cycle (organic and dissolved)
6. Nitrogen cycle (organic, ammonia, nitrite, and nitrate)
7. Dissolved oxygen
8. Coliforms
9. An arbitrary nonconservative constituent
All parameters can be simulated under either steady-state or dynamic
conditions. If either the phosphorus cycle or the nitrogen cycle are
not being simulated, the model presumes they will not limit algal growth.
70
-------
5.4 NUMERICAL SOLUTION TECHNIQUE
At each time step and for each constituent, Equation II-3 can be
written I times, once for each of the I computational elements in the
network, because it is not possible to obtain analytical solutions to
these equations under most prototype situations, a finite difference
method is usedmore specifically, the classical implicit backward differ-
ence method (Arden and Astill, 1970; Smith, 1966; and Stone and Brian,
1963).
The general basis of a finite difference scheme is to find the value
of a variable (e.g., constituent concentration) as a function of space at
a time step n+1 when its spatial distribution at the nth time step is
known. Time step zero corresponds to the initial condition. Backward
difference or implicit schemes are characterized by the fact that all
spatial derivatives (9/3x) are approximated in difference form at
time step n+1.
5.4.1 Formulation of the Finite Difference Scheme
The finite difference scheme is formulated by considering the consti-
tuent concentration, C, at four points in the mnemonic scheme as shown in
Figure V-l.
Three points are required at time n+1 to approximate the spatial
derivatives. The temporal derivative is approximated at distance step i.
Equation II-3 can be written in finite difference form in two steps.
First, the advection and diffusion terms are differentiated once with
respect to x, giving:
DOWNSTREAM4-
UPSTREAM
element i-H
element i
I
o-
N : t -
At
i
Figure V-l. Classical Implicit Nodal Scheme
71
-------
(ADL -) - (ADL ~)
aci a* i 3x i_! (A u
at
+ __ + _ v-i
dt Vi
where
Vi = Ai A Xi
Secondly, expressing the spatial derivative of the diffusion terms in finite
difference and thence the time derivative of C in finite difference, there
results:
T] c"+} - [(ADL)J] cf *
At V-j A Xj
s
. rn+1
) + ri Cn+1 + Pi +
In equation V-2, the term dC/dt is expressed as:
dC,-
^ = r. cn+l + D.
« r1 v,1 T PI
dt
where
^ = rate constant
Pi = internal constituent sources and sinks (e.g., nutrient
loss from algal growth, benthos sources, etc.)
Note that the dC/dt for every constituent modeled by QUAL2E can be expressed
in this form.
72
-------
If equation V-2 is rearranged in terms of the coefficients of
Cn+1, and C^}, we obtain the equation:
b1
where
At Qi_i At
V-3
= 1.0 + [(ADL)i +
At At
+ Qi - - r. At
Ci = -[(ADL)i
7. = r.n .
*i ui
At
At
At
The values of ai, b.j, Ci, and 1^ are all known at time n, and the C^
terms are the unknowns at time step n+1.
In the case of a junction element with a tributary upstream element,
the basic equation becomes:
a1 Cnt} + bi Cf1 + Ci Cn!j + dj Cn+1 = Zi
V-4
where
At
At
j = the element upstream of junction element i
Cn+1 = concentration of constituent in element j at time n+1
It can be seen that the dj term is analogous to the a^ term. Both
terms account for mass inputs from upstream due to dispersion and advec-
tion.
73
-------
oU-j
Under steady-state conditions, = 0 in equation V-l. Working
3t
through the finite difference approximations and rearranging terms as
before, the steady-state version of equation V-3 is derived:
n+l + h. rO+l + c- rn+l =
Di h + ci H-l
where
(ADL)l
Si
Note that equation V-5 is the same as equation V-3, with three
changes:
o At = 1.0
o the constant 1.0 in b-j = 0.0
o the initial concentration Cn in Z^ = 0.0
5.4.2 Method of Solution
Equations V-3 and V-5 each represent a set of simultaneous linear
equations whose solution provides the values of Cn+1 for all i's.
Expressed in matrix form, this set of equations appears as:
74
-------
bl cl
3o bo C2
a3 b3 C3
ai bi ci
al-l bl-l CI-1
al bl
X
V "
c||+1
cfi
T1
'PI
cf1
=
zl
Z2
Z3
Zi
ZM
zi
V-6
The left matrix is a tri-diagonal matrix. An efficient method that readily
lends itself to a computer solution of such a set of equations is:
Diyide through the first equation in V-6 by b^ to obtain:
Cn+1 + Wx
V-7
where
and G^ = I\l\>i.
Combine the expression for b^ (see V-3) and the second equation in
V-6 to eliminate ag and the result is:
+1
W2 c§+ = G2
V-8
where
W2 =
wl
and 62 =
- 82
Combine equation V-8 and the third equation in V-6 to eliminate 33
and the result is:
75
-------
Cf 1 + W3 Cf 1 = G3 V-9
where
c3 Z3 - a3 G2
W3 = i and 63 = .
b3 - a3 W2 b3 - a3 W2
Proceed through the equations, eliminating a-j and storing the values
of W-j and G-j given by:
ci
Wi = . i = 2, 3, .... I V-10
bi - ai Wi_i
and
zi - ai
Gi = - , i = 2, 3, .... I V-ll
bi - ^i wi-l
The last equation is solved for Cn+1 by
Cn+1 = Gj V-12
Solve for Cn^J, Cn+J, . . . , C}+1 by back substitution.
Cn+1 = Gi - W,. Cnlj, i = 1-1, 1-2, . . . , 1 V-13
5.4.3 Boundary Conditions
In most situations of interest, transport is unidirectional in nature,
i.e., there is no significant transport upstream. Therefore, the concen-
tration at some point just upstream from the beginning or end of the
stream reach of interest can be used as the boundary condition.
5.4.3.1 Upstream Boundary (Headwater Elements)
For headwater elements there is no upstream, i-1, element. Thus,
the headwater driving force is substituted in Equation V-3 for the upstream
concentration G-J_I. Because the headwater concentrations are fixed, they
are incorporated on the right hand side of Equation V-3 in the known term Z-j ,
for headwater elements as follows.
SiAt
Z, = Cin + + pi At - ai Cn V-14
76
-------
where C0 is the upstream boundary condition (headwater concentration).
5.4.3.1 Downstream Boundary (Last Element in the System)
QUAL2E has two options for modeling the downstream boundary. One
uses a zero gradient assumption; the other incorporates fixed downstream
constituent concentrations into the solution algorithm.
Zero Gradient Assumption (Arden and Astill, 1970)--For the last
computational element in the system, there is no downstream, i+1, element.
At this boundary, a zero gradient assumption is made that replaces C-j+i
with C-j_i. In this manner, the downstream boundary acts as a mirror to
produce a zero gradient for the concentration of the constituent variable.
The coefficient a-j, therefore, is modified to include the dispersion effect
normally found in the coefficient c-j for the last element in the system.
Thus, the equation for a-j in V-3 becomes:
At Qi_iAt
ai = -[((ADL)M + (ADL)I) + ] .v-is
VjAxj V!
and
G! = 0 V-16
where I = index of the downstream boundary element
Fixed Downstream Constituent ConcentrationsFor this boundary
option, the user supplies known downstream boundary concentrations C[_B
for each water quality constituent. Thus, the value of Cj+i in Equation
V-3 becomes
CI+1 = CLB V-17
Because the boundary concentrations are known in this option, they are
incorporated on the right hand side of Equation V-3 in the known term Z-j
for the downstream boundary element then results as
SjAt
zi = c! + + PiAt - ci CLB v-is
VI
77
-------
6. REFERENCES
American Public Health Association, Inc., Standard Methods for the
Examination of Water and Wastewater, American Public Health Associa-
tion, 1965 (12th edition), 1985 (16th edition).
Anderson, E.R., Energy Budget Studies in Water Loss InvestigationsLake
Hefner Studies, Technical Report, U.S. Geological Survey, Washington,
DC, Prof. Paper 269, 1954.
Arden, B.W. and K.N. Astill, Numerical Algorithms: Origins and Applica-
tions, Addison-Wesley, Reading, MA, 1970.
Bannister, T.T., "Production Equations in Terms of Chlorophyll Concentra-
tion, Quantun Yield, and Upper Limit to Production," Limnology and
Oceanography, Vol. 19, No. 1, pp 1-12, January 1974.
Brandes, R.J. and A.B. Stein, WREDUN Model Documentation Report, Water
Resources Engineers, Inc., prepared for Texas Department of Water
Resources, Construction Grants and Water Quality Planning Division,
no date.
Chen, C.W. and G.T. Orlob, Final Report, Ecologic Simulation of Aquatic
Environments, Water Resources Engineers, Inc., prepared for the
Office of Water Resources Research, U.S. Department of the Interior,
Washington, DC, October 1972.
Churchill, M.A., H.L. Elmore and R.A. Buckingham, "The Prediction of
Stream Reaeration Rates," International Journal of Air and Water
Pollution, Vol 6, pp 467-504, 1962.
DeGroot, W.T., "Modelling the Multiple Nutrient Limitation of Algal
Growth," Ecological Modelling, Vol. 18, pp 99-119, 1983.
Duke, James H., Jr., Provision of a Steady-State Version of the Stream
Model, QUAL, Water Resources Engineers, Inc., Austin, TX, prepared
for the U.S. Environmental Protection Agency, Washington, DC,
November 1973.
Edinger, J.E. and J.C. Geyer, Heat Exchange in the Environment, Johns
Hopkins Univ., Baltimore, MD, 1965.
78
-------
Elder, J.W., "The Dispersion of a Marked Fluid in Turbulent Shear Flow,"
Jour. Fluid Mech.. Vol. 5, Part 4, pp 544-560, May 1959.
Field, S.D. and S.W. Effler, "Photosynthesis-Light Mathematical Formula-
tions," Journal of Environmental Engineering Division, ASCE, Vol.
108, No. EE1, pp 199-203, February 1982.
Field, S.D. and S.W. Effler, "Light-Productivity Model for Onondaga, N.Y.,"
Journal Environmental Engineering Division, ASCE, Vol. 109, No. EE4,
pp 830-844, August 1983.
Fisher, H.B., Discussion to "Time of Travel of Soluble Contaminants in
Streams," by T.J. Buchanan, Proc. Sanitary Eng. Div., ASCE, v. 6,
1964.
Fisher, H.B., E.J. List, R.C.Y. Koh, J. Imberger, N.H. Brooks, Mixing in
Inland and Coastal Waters, Academic Press, New York, NY, 1979.
Frank D. Masch and Associates and the Texas Water Development Board,
Simulation of Water Quality in Streams and Canals, Theory and
Description of the QUAL-I Mathematical Modeling System, Report 128,
the Texas Water Development Board, Austin, TX, May 1971.
Henderson, F.M;, Open Channel Flow, Macmillan Co., New York, NY, 1966.
JRB Associates, "Users Manual for Vermont QUAL-II Model," prepared for
U.S. Environmental Protection Agency, Washington, DC, June 1983.
Jassby, A.D., and T. Platt, "Mathematical Formulation of the Relationship
Between Photosynthesis and Light for Phytoplankton," Limnology and
Oceanography. Vol. 21, No. 4, July 1976, pp 540-547.
Kramer, R.H., A Search of the Literature for Data Pertaining to
Bioenergetics and Population Dynamics of Freshwater Fishes, Desert
Biome Aquatic Program, Utah State University, Logan, UT, August
1970.
Langbien, W.B. and W.H. Durum, The Aeration Capacity of Streams, U.S.
Geological Survey, Washington, DC, Circ. 542, 1967.
National Council for Air and Stream Improvement, Inc., A Review of
the Mathematical Water Quality Model QUAL-II and Guidance for its
Use. NCASI, New York, NY, Technical Bulletin No. 391, December 1982.
National Council for Air and Stream Improvement, Inc., QUAL2E User Manual,
NCASI, New York, NY, Technical Bulletin No. 457, April 1985.
O'Connor, D.J. and W.E. Dobbins, "Mechanism of Reaeration in Natural Streams,"
Trans. ASCE, Vol. 123, pp 641-684, 1958.
79
-------
Owens, M., R.W. Edwards and J.W. Gibbs, "Some Reaeration Studies in Streams,"
International Journal of Air and Water Pollution, Vol. 8, No. 8/9,
pp 469-486, September 1964.
Platt, T., K.H. Mann, and R.E. Ulanowicz (eds.), Mathematical Models in
Biological Oceanography, Unesco Press, Paris, 1981.
Roesner, L.A., Temperature Modeling in Streams, Lecture notes, water
quality workshop, Tennessee Valley Authority, Knoxville, TN, 1969.
Roesner, L.A., P.R. Giguere, and D.E. Evenson, Computer Program Documentation
for Stream Quality Modeling (QUAL-II). U.S. Environmental Protection
Agency, Athens, GA, EPA-600/9-81-014, February 1981b.
Roesner, L.A., P.R. Giguere, and D.E. Evenson, User's Manual for Stream
Quality Model (QUAL-II). U.S. Environmental Protection Agency,
Athens, GA, EPA-600/9-81-015, February 1981b.
Scavia, D. and R.A. Park, "Documentation of Selected Constructs and
Parameter Values in the Aquatic Model CLEANER," Ecological Modeling,
Vol. 2, pp 33-58, 1976.
Smith, E.L., "Photosynthesis in Relation to Light and Carbon Dioxide,"
Proceedings, National Academy of Sciences, Vol. 22, pp 504-510, 1936.
Smith, J.D., Solutions to Partial Differential Equations, Macmillan Co.,
New York, NY, 1966.
St. John, J.R., T.W. Gallagher, and P.R. Paquin, "The Sensitivity of the
Dissolved Oxygen Balance to Predictive Reaeration Equations," in Gas
Transfer at Water Surfaces, W. Brutsaert and G. Jirka, eds., D. Reidl
Publishing Company, Dordrecht, Holland, 1984.
Steele, J.H., "Environmental Control of Photosynthesis in the Sea,"
Limnology and Oceanography, Vol. 7, pp 137-150, 1962.
Stefan, H.G., J.J. Cardoni, F.R. Schiebe, and C.M. Cooper, "Model of Light
Penetration in a Turbid Lake," Water Resources Research, Vol. 19,
No. 1, pp 109-120, February 19BT.
Stone, H.L. and P.O.T. Brian, "Numerical Solution of Convective Transport
Problems," Journal American Institute of Chemical Engineers, Vol. 9,
No. 5, pp 681-688, 1963.
Streeter, H.W. and E.B. Phelps, Study of the Pollution and Natural Purifi-
cation of the Ohio River. U.S. Public Health Service, Washington, DC,
Bulletin No. 146 (reprinted 1958), 1925.
Swartzman, G.L. and R. Bentley, "A Review and Comparison of Phytoplankton
Simulation Models," Journal of the International Society for Ecological
Modelling. Vol. 1, Nos. 1-2, pp 30-81, 1979.
80
-------
Taylor, G.I., "The Dispersion of Matter in Turbulent Flow Through a Pipe,"
Proceedings. Royal Society of London, Vol. 234A, No. 1199, pp 456-475,
March 6, 1954.
Texas Water Development Board, Simulation of Water Quality in Streams
and Canals, Program Documentation and User's Manual, Austin, TX,
September 1970.
Texas Water Development Board, QUAL-TX Users Manual, Version 2.5, Water
Quality Management Section, Austin, TX, November 1984.
TenEch Environmental Consultants, Inc. Waste Load Allocation
Verification Study: Final Report. Prepared for Iowa Department
of Environmental Quality, July 1978.
Thackston, E.L. and P.A. Krenkel, Reaeration Prediction in Natural Streams,
Journal of the Sanitary Engineering Division, ASCE, Vol. 95, No. SA1,
pp 65-94, February 1969.
Thomas, H.A., Jr., "Pollution Load Capacity o.f Streams." Water and Sewage
Works. Vol. 95, No. 11, pp 409-413, November 1948.
Tsivoglou, E.C. and J.R. Wallace, Characterization of Stream Reaeration
Capacity, Prepared for U.S. Environmental Protection Agency,
Washington, DC, 1972.
Tsivoglou, E.C. and L.A. Neal, "Tracer Measurement of Reaeration: III.
Predicting the Reaeration Capacity of Inland Streams," Jour.
WPCF, Vol. 48, No. 12, pp 2669-2689, December 1976.
Walker, W.W., QUAL2 Enhancements and Calibration to the Lower Winooski.
Prepared for the Vermont Agency of Environmental Conservation,
Montpelier, VT, December 1981.
Walker, W.W. Personal Communication, 1983.
Water Resources Engineers, Inc., Prediction of Thermal Energy Distribution
in Streams and Reservoirs, Prepared for the California Dept. of Fish
and Game, 1967.
Water Resources Engineers, Inc. Progress Report on Contract No. 68-01-0713,
Upper Mississippi River Basin Model Project, Sponsored by the
Environmental Protection Agency, submitted to Environmental Protection
Agency, September 21, 1972.
Wunderlich, W.O., The Fully-Mixed Stream Temperature Regime. ASCE Specialty
Conf., Utah State Univ., Logan, Utah, 1969.
Zison, S.W., D.B. Mills, D. Deimer, and C.W. Chen, Rates, Constants, and
Kinetics Formulations in Surface Water Quality Modeling. Prepared
for USEPA Environmental Research Lab., Athens, Georgia, September 1978,
revised 1985, EPA-600/3-85-040.
81
-------
APPENDIX A
QUAL2E User Manual*
The following sections illustrate the coding of input data
forms.
A. Title Data
All 16 cards are required in the order shown. The first
two are title cards, and columns 22 through 80 may be used to
describe the base, date of simulation, etc. Title cards 3
through 15 require either a "yes" or "no" in columns 10 through
12 and are right justified. Note that each of the nitrogen and
phosphorus series must be simulated as a group.
For each conservative substance (up to three) and the
arbitrary non-conservative, the constituent name must be
entered in columns 49 through 52. Corresponding input data
units are entered in columns 57 through 60 (e.g.. mg/L).
QUAL2E simulates ultimate BOD in the general case. If the
user wishes to use 5-day BOD for input and output, the program
will internally make the conversions to ultimate BOD. This
conversion is based upon first order kinetics and a decay rate
that can be specified by the user (Type 1 Data, line 8). If no
value is specified, the program uses a default value of 0.23
per day, base e. NCASI recommends that users work only with
ultimate BOD unless they have detailed knowledge of the river
water and point source BOD kinetics. To use the 5-day BOD
input/output option, write "5-Day Biochemical Demand" in mg/L
on the title 7 card beginning in column 22.
Card 16 must read ENDTITLE beginning in column 1.
*From: Modifications to the QUAL-2 Water Quality Model and User
Manual for QUAL-2E Version 2.2. National Council of the Paper
Industry for Air and Stream Improvement, Inc., New York, NY.
NCASI Tech. Bulletin No. 457. April 1985. Used by permission.
82
-------
B. Data Type 1 - Program Control
Type 1 Data define the program control options and the
characteristics of the stream system configuration, as well as
some of the geographical/meteorological conditions for modeling
temperature. There are a maximum of 17 Data 1 cards. The
first 13 are required, while the last four are necessary only
if temperature is being simulated.
The QUAL2E program recognizes Type 1 Data by comparing the
first four characters (columns 1-4) of each data card with a set
of internally fixed codes. If a match between the code and
characters occurs, then the data are accepted as supplied on
the card by the user. If a match does not occur, then the
program control options will revert to default values and the
system variables for the unmatched codes will be assigned the
value zero (0.0).
The first seven cards control program options. If any
characteristics other than those shown below are inserted in
the columns 1 through 4, the actions described will not occur.
LIST - Card T. list the input data.
WRIT - Card 2. write the intermediate output report,
WRPT2 (see SUBROUTINE WRPT2 in the documentation
manual, NCASI Technical Bulletin No. 391).
FLOW - Card 3, use flow augmentation.
STEA - Card 4 shows this is a steady-state simulation. If it
is not to be a steady-state, write DYNAMIC SIMULATION,
or NO STEADY STATE, and it is automatically a dynamic
simulation.
TRAP - Card 5, cross-sectional data will be specified for each
reach. If discharge coefficients are to be used for
velocity and depth computations, write DISCHARGE
COEFFICIENTS, or NO TRAPEZOIDAL CHANNELS, beginning in
column 1.
PRIN - Card 6. local climatological data specified globally for
the basin simulation will appear in the final output
listing.
PLOT - Card 7. dissolved oxygen and BOD will be plotted in final
output listing.
The next two cards provide further program flags and coeffi-
cients. This information is supplied in two data fields per
card; columns 26-35. and 71-80. Note that the character codes
83
-------
in columns 1-4 must occur as shown in order for the data to be
accepted by the program.
FIXE - Card 8. specifies: (a) whether the downstream boundary
water quality constituent concentrations are fixed (user
specified), and (b) the value of the rate coefficient
for converting input 5-day BOD to ultimate BOD. A value
of 1.0 (or larger) in columns 26-35 specifies that the
downstream boundary water quality constituent
concentrations will be supplied in Data Types 13 and
13A. A value less than 1.0 (usually 0.0 or blank) in
these columns means that the downstream boundary
concentrations are not user specified. In this case,
the concentrations in the most downstream element (Type
5) will be computed in the normal fashion using the zero
gradient assumption (1, 11). The second value on this
card, columns 71-80, is the rate coefficient for
converting 5-day to ultimate BOD. It is used only when
5-day BOD is being modeled (Title Card 7). If the
columns are left blank, the model uses a default value
of 0.23 per day, base e. Note that this conversion
factor is applied to all input BODg forcing functions
(headwaters, incremental flows, point loads, and the
downstream boundary condition).
INPU - Card 9, specifies whether the user will input and/or
output in metric or English units. The value of 1.0 (or
larger) in card columns 26-35 specifies metric input.
The value of 1.0 (or larger) in card column 71-80
specifies metric units for output. Any value less than
1.0 (usually 0.0 or blank) will specify English units.
The next four cards describe the stream system. There are
two data fields per card, columns 26-35 and 71-80. The program
restrictions on the maximum number of headwaters, junctions,
point loads, and reaches are defined by PARAMETER statements in
the Fortran code. These statements may be modified by the user
to accommodate a particular computer system or QUAL2E
simulation application. The values of the maximum constraints
are as follows:
Maximum number of headwaters 10
Maximum number of junctions 9
Maximum number of point loads 50
Maximum number of reaches 50
Maximum number of computational elements 500
NUMB - Card 10, defines the number of reaches into which the
stream is segmented and the number of stream junctions
(confluences) within the system.
84
-------
NUM_ - Card 11. shows the number of headwater sources and the
number of inputs or withdrawals within the system. The
inputs can be small streams, wasteloads, etc.
Withdrawals can be municipal water supplies, canals.
etc. NOTE: Withdrawals must have a minus sign ahead of
the flow in Data Type 11 and must be specified as
withdrawals in Data Type 4 by setting IFLAG = 7 for that
element. Note, the code for Card 11 is 'NUM_' (read:
NUM space) to distinguish it from the code for Card 10.
NUMB.
TIME - Card 12. contains the time step interval in hours and the
length of the computational element in miles
(kilometers). The time step interval is used only for a >
dynamic simulation, thus it may be omitted if the
simulation is steady-state.
MAXI - Card 13. provides information with different meanings de-
pending on whether a dynamic or steady-state simulation
is being performed. For a dynamic simulation, the
maximum route time is specified in columns 26-35. This
value represents the approximate time in hours required
for a particle of water to travel from the most upstream
point in the system to the most downstream point. The
time increment in hours for intermediate summary reports
of concentration profiles is specified in columns
71-80. For a steady-state simulati-on. the maximum
number of iterations allowed for solution convergence is
entered in columns 26-35. The value in columns 71-80
may be left blank because it is not required in the
steady-state solution.
The next four cards provide geographical and meteorological
information and are required only if temperature is being
simulated. There are two data fields per card, columns 26-35
and 71-80. Note: the character codes in columns 1-4 must
occur as shown in order for the data to be accepted by the
program.
LATI - Card 14, contains the basin latitude and longitude and
represent mean values in degrees for the basin.
STAN - Card 15. shows the standard meridian in degrees, and the
day of the year the simulation is to begin.
EVAP - Card 16. specifies the evaporation coefficients. Typical
values are AE = 6.8 x 10~4 ft/hr-in Hg and BE = 2.7 x
10~4ft/hr-in Hg-mph of wind for English units input.
or AE = 6.2 xlO~6 m/hr-mbar and BE = 5.5 x 10~6
m/hr-mbar-m/sec of wind for metric units input.
85
-------
ELEV - Card 17. contains the mean basin elevation in feet
(meters) above mean sea level, and the dust attenuation
coefficient (unitless) for solar radiation. The dust
attenuation coefficient generally ranges between zero
and 0.13. Users may want to consult with local
meteorologists for more appropriate valiies.
Data Type 1 must end with an ENDATA1 card.
C. Data Type 1A - Global Algal. Nitrogen, Phosphorus, and Light
Parameters
These parameters and constants apply to the entire
simulation and represent the kinetics of the algal, nutrient.
light interactions. It is important to note that proper use of
all options in QUAL2E requires detailed knowledge of the algal
growth kinetics appropriate for the water body being simulated.
These data cards are required only if algae, the nitrogen
series (organic, ammonia, nitrite, and nitrate), or the
phosphorus series (organic and dissolved) are to be simulated.
Otherwise they may be omitted, except for the ENDATA1A card).
Information is supplied in two data fields per card, columns
33-39 and 74-80. As with Type 1 Data, the QUAL2E program
recognizes Type 1A Data by comparing the first four characters
(columns 1-4) of each card with a set of internally fixed
codes. If a match between the codes and the characters occurs.
then the data are accepted as supplied on the card by the
user. If a match does not occur, then the system variables for
the unmatched codes will be assigned the value zero (0.0).
Note: the spaces (under bars) are an integral (necessary) part
of the four character code.
O_UP - Card 1, specifies the oxygen uptake per \rnit of ammonia
oxidation, and oxygen uptake per unit of nitrite
oxidation.
O_PR - Card 2. contains data on oxygen production per unit of
algae growth, usually 1.6 mg O/mg A. with a range of 1.4
to 1.8. It also contains data on oxygen uptake per unit
of algae, usually 2.0 mg O/mg A respired, with a range
of 1.6 to 2.3.
N_CO - Card 3. concerns the nitrogen content and phosphorus
content of algae in mg per mg of algae. The fraction of
algae biomass which is nitrogen is about 0.08 to 0.09.
and the fraction of algae biomass which is phosphorus is
about 0.012 to 0.015.
86
-------
ALG_ - Card 4. specifies the growth and respiration rates of
algae. The maximum specific growth rate has a range of
1.0 to 3.0 per day. The respiration value of 0.05 is
for clean streams, while 0.2 is used where the NE and
P2 concentrations are greater than twice the half
saturation constants.
N_HA - Card 5. contains the nitrogen and phosphorus half
saturation coefficients. The range of values for
nitrogen is from 0.2 to 0.4 mg/L and for phosphorus the
value normally used is 0.04 mg/L.
LIN_ - Card 6, contains the linear and nonlinear algal self
shading light extinction coefficients. The coefficients
\! and X-2 are defined below.
\1 = linear algae self-shading coefficient
(l/ft)/(ug chla/L), or (l/m)/(ug chla/L)
\2 = nonlinear algae self-shading coefficient
(l/ft)/ug chla/L)2/3. or (l/m)/(ug chla/L)2/3
These two self-shading coefficients are used with
\o. the non-algal light extinction coefficient (Type
6B Data) in the general light extinction eguation shown
below:
\ = \0 + X^A +
where X is the total light extinction coefficient and
A is the algae biomass concentration in mg A/L and
-------
LIGH - Card 7. contains the solar light function option for
computing the effects of light attenuation on the algal
growth rate, and the light saturation coefficient.
QUAL2E recognizes three different solar light function
options. The light saturation coefficient is coupled to
the selection of a light function, thus care must be
exercised in specifying a consistent pair of values.
The depth integrated form of the three light
functions and the corresponding definitions of the light
saturation coefficient are given in Figure A-l and
outlined in the following table.
Light Function Option
(Columns 33 - 39)
1 (Half Saturation)
2 (Smith's Function)
3 (Steele's Function)
Light Saturation Coefficient*
(Columns 74 - 80)
Half Saturation Coefficient
Light intensity corresponding
to 71% of maximum growth rate
Saturation Light Intensity
* Units of the Light Saturation Coefficient are as
follows:
English: BTU/ft2-min and Metric: Langleys/min
Light Function Option 1 utilizes a Michaelis-Menton
half saturation formulation for modeling the algal
growth limiting effects of light (FL). It is the method
used in the SEMCOG version of QUAL-2. Option 2 is
similar to Michaelis-Menton, but utilizes a second order
rather than first order light effect. Both options 1
and 2 are monotonically increasing functions of light
intensity. Option 3 includes a photo-inhibition effect
at high light intensities and has been reported in
(9.13).
DAIL - Card 8. contains the light averaging option (columns
33-39) and the light averaging factor (columns 74-80).
These values are used only in a steady-state simulation.
The light averaging option allows the user to specify the
mannerin which the light attenuation factor (FL in Figure
A-l) is computed, from the available values of solar
radiation. A summary of these options is given below.
88
-------
Option Description
1 FL is computed from one daily average solar
radiation value calculated in the steady-
state temperature subroutine (HEATER).
2 FL is computed from one daily average solar
radiation read from Data Type 1A.
3 FL is obtained by averaging the 24 hourly
values of FL. that are computed from the
24 hourly values of solar radiation calcu-
lated in the steady-state temperature sub-
routine (HEATER).
4 FL is obtained by averaging the 24 hourly
values of FL. that are computed from the
24 hourly values of solar radiation
computed from the total daily solar
radiation (Data Type 1A) and an assumed
cosine function.
Note: that if options 1 or 3 are selected, temperature
must be simulated.
The light averaging factor (columns 74-80) is used to
make a single calculation using daylight average solar
radiation (Option 1 or 2) agree with average of
calculations using hourly solar radiation values (Option
3 or 4). The factor has been reported to vary from 0.85
to 1.00.
The selection of a daily (diurnal) light averaging option
depends largely on the detail to which the user wishes to
account for the diurnal variation in light intensity. Options
1 and 2 utilize a single calculation of FL based on an average
daylight solar radiation value. Options 3 and 4 calculate
hourly values of FL from hourly values of solar radiation and
then average the hourly FL values to obtain the average
daylight value. Options 1 and 3 use the solar radiation from
the temperature heat balance routines (thus both algae and
temperature simulations draw on the same source for solar
radiation). Options 2 and 4 use the solar radiation value in
Data Type 1A for the algae simulation. Thus either option 2 or
4 must be selected when algae are simulated and temperature is
not. The light averaging factor is used to provide similarity
in FL calculations between options 1 and 2 versus options 3 and
4. The solar radiation factor (Data Type 1A. card 11)
specifies the fraction of the solar radiation computed in the
heat balance which is photosynthetically active. It is used
only with options 1 or 3.
89
-------
Light Functions (LFNOPT)
Option 1: Half saturation (SEMCOG)
FL = (-i) In [ KL + \ ,]
Ad K L + I e
KL = light intensity at which growth rate is 50% of the
maximum growth rate.
Option 2: Smith's Function (14)
1 T /If T -L /I J. t T /VT \ \ ^f ^"
FL = (-i) In [ ,H H , iV,]
Ad I/KLe-Xd + (1 + (I/KLe-Xd)2)1/2
KL = light intensity at which growth rate is 71% of the
maximum growth rate.
Option 3: Steele's Equation (9)
FL = 2.718 [e-(e"Xd
Ad
KL = light intensity at which growth rate is equal to the
maximum growth rate.
Notation: FL = light attention factor
X = extinction coefficient *
d = depth
I = surface light intensity
KL = light saturation coefficient
* Algal Self-shading
General Equation X . Xo f Xx oo A + X2 (a0A)2/3
where: X » light extinction coefficient
X9 » non-algal extinction
\l - linear algal self shading coefficient
\2 » non-linear algal self shading coefficient
A « algal biomass concentration (mg/L)
°0 - chlorophyll a to algae biomass ratio (ug chla/ragA)
Special Cases
No Self-shading (SEMCOG)
Linear Self-shading (META)
Non-linear Self-shading (TetraTech)
xl f X2 f °
e.g. X » XQ + 0.0088^ + 0.054(aQA)2/3 (Riley Eq.,
metric units)
FIGURE A-l. ALGAL GROWTH RATE
90
-------
In dynamic algae simulations option 3 is used (default)
unless temperature is not simulated in which case solar
radiation data are read in with the local climatology data.
NUMB - Card 9. contains the number of daylight hours (columns
33-39). and the total daily radiation (BTU/ft2, or
Langleys) (columns 74-80). This information is used if
light averageing options 2 or 4 are specified for the
simulation.
ALGY - Card 10, contains the light-nutrient option for computing
the algae growth rate (columns 33-39). and the algal
preference factor for ammonia nitrogen (columns 74-80).
The light-nutrient interactions for computing algae
growth rate are as follows (see also Figure A-l):
Option
Description
Multiplicative: (FL) * (FN) * (FP)
Limiting Nutrient: FL * [minimum (FN. FP)]
Inverse Additive: FL * 2
1/FN + 1/FP
Option 1 is the form used in QUAL-2 SEMCOG, while option 2
is used in the META Systems Version of QUAL-2. Option 3 is
described by Scavia and Park (13).
The algal preference factor for ammonia (columns 74-80)
defines the relative preference of algae for ammonia and
nitrate nitrogen (see also Figure A-2). The user defines
this preference by specifying a decimal value between 0
and 1.0, for example:
Algal Preference
factor
for Ammonia
Interpretation
0.0
0.5
1.0
Algae will use only nitrate for growth.
Algae will have equal preference for
ammonia and nitrate.
Algae will use only ammonia for growth,
91
-------
Nutrient Attenuation Factors
Nitrogen: NE = NH3 + NC>3 (mg-N/L)
FN = NE/(KN + NE)
where: FN = Nitrogen attenuation factor
KN = Nitrogen half saturation coefficient (mg-N/L)
Phosphorus: FP = ?2/(KP + P£)
?2 = Dissolved Phosphorus (mg-P/L)
where: FP = Phosphorus attenuation factor
KP = Phosphorus half saturation coefficient (mg-P/L)
Algal Preference foe Ammonia
F =
(PN)
(PN)(NH3) + (1-PN)N03
where: PN = Algal preference for ammonia nitrogen (0-1.0)
PN= 0.0; algae will use only nitrate for growth
PN= 0.5; algae will have equal preference for
ammonia and nitrate for growth
PN= 1.0; algae will use only ammonia for growth
F = Fraction of algal nitrogen uptake from ammonia pool.
FIGURE A-2. ALGAL GROWTH RATE
ALG/ - Card 11. contains the factor for converting the solar
radiation value from the heat balance to the solar
radiation value appropriate for the algae simulation
(columns 33-39) and the value of the first order
nitrification inhibition coefficient (columns 74-80).
The solar radiation factor specifies the fraction of
the solar radiation computed in the heat balance
(subroutine HEATER) that is photosynthetically active
(i.e.. used by algal cells for growth). It is required
only in steady-state simulations when light averaging
options 1 or 3 (Data Type 1A, card 8) are selected. A
decimal value between 0 and 1.0 specifies the value of
this fraction. Typically the value of this fraction is
about 0.45 (14).
The first order nitrification inhibition coefficient
is the value of KNITRF in the following equation (see
Figure A-3).
CORDO = 1.0 - EXP (-KNITRF *DO)
92
-------
where:
DO = dissolved oxygen concentration (mg/L), and
CORDO = correction factor applied to ammonia and nitrite
oxidation rate coefficients.
Nitrification Rate Correction Factor (CORDO)
CORDO = 1.0 - EXP(-KNITRF * DO)
The value of KNITRF is supplied by the user,
the default value in QUAL-2E is 10.0
Applied to Ammonia and Nitrite Oxidation Rates
Ammonia: ^l^inhib = ^l^ino t * (CORD°)
Nitrite: (Sn),_u,u = (B0).__..1. * (CORDO)
'Vinhib.
Magnitude of Correction Factor
'Vinput
The following table contains values of CORDO as a function
of DO (row) and KNITRF (column).
DO
(mg/L)
0.1
0.2
0.3
0.4
0.5
0.7
1.0
1.5
2.0
3.0
4.0
5.0
7.0
10.0
0.5
.05
.10
.14
.18
.22
.30
.39
.53
.63
.78
.86
.92
.97
.99
0.7
.07
.13
.19
.24
.30
.39
.50
.65
.75
.88
.94
.97
.99
1.00
KNITRF
1.0
.10
.18
.26
.33
.39
.50
.63
.78
.86
.95 1
.98 1
.99 1
1.00 1
1.00 1
2.0
.18
.33
.45
.55
.63
.75
.86
.95
.98
.00
.00
.00
.00
.00
5.0
.39
.63
.78
.86
.92
.97
.99
1.00
1.00
1.00
1.00
1.00
1.00
1.00
10.0
.63
.86
.95
.98
.99
1.00
1.00
1.00
1.00
1.00
1.00
1.00
1.00
1.00
KNITRF = 0.6, closely matches the nitrification inhibition
function in QUAL-2 Texas (4).
KNITRF
FIGURE A-3.
0.7, closely matches the nitrification inhibition
data of the Thames Estuary (12).
NITRIFICATION INHIBITION AT LOW DO
93
-------
A value of 0.6 for KNITRF closely matches the
inhibition formulation in QUAL-Z Texas, while a value of
0.7 closely matches the data for the Thames Estuary (12).
The default value of KNITRF is 10.0, i.e., no inhibition
of nitrification at low dissolved oxygen. A table of
CORDO values as a function of KNITRF and DO is given in
Figure A-3.
ENDA - The last card in Data Type 1A must be an ENDATA1A card,
regardless of whether algae, nitrogen, or phosphorus are
simulated.
D. Data Type IB - Temperature Correction Factors
Several of the processes represented in QUAL2E are affected
by temperature. The user may elect to input specific
temperature correction factors. In the absence of such
information, default values are used as noted in Figure A-4.
The user need supply only those values that are to be changed.
Data Type IB information is supplied as follows:
Alphanumeric code for each temperature
coefficient as noted in Figure A-4; Columns 10-17
User specific temperature coefficient Columns 19-26
The last card in Data Type IB must be an ENDATA1B card,
regardless of whether any of the default values are modified.
E. Data Type 2 - Reach Identification and River Mile/Kilometer
Data
The cards of this group identify the stream reach system by
name and river mile/kilometer by listing the stream reaches
from the most upstream point in the system to the most
downstream point. When a junction is reached, the order is
continued from the upstream point of the tributary. There is
one card per reach. The following information is on each card:
Reach Order or Number Columns 16-20
Reach Identification or Name Columns 26-40
River Mile/Kilometer at Head of Reach Columns 51-60
River Mile/Kilometer at End of Reach Columns 71-80
94
-------
DEFAULT VALUES
INDEX RATE COEFFICIENT
1 BOD Decay
2 BOD Settling
3 Reaeration
4 SOD Uptake
5 Organic N Decay
6 Organic N Settling
7 Ammonia Decay
8 Ammonia Source
9 Nitrite Decay
10 Organic P Decay
11 Organic P Settling
12 Dissolved P Source
13 Algae Growth
14 Algae Respiration
15 Algae Settling
16 Coliform Decay
17 Non-cons Decay
18 Non-cons Settling
19 Non-cons Source
FIGURE A-4. DEFAULT THETA VALUES FOR QUAL2E
SEMCOG
1.047
-
1.0159
-
-
-
1.047
-
1.047
-
-
-
1.047
1.047
-
1.047
1.047
-
-
QUAL-2E
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
.047
.024
.024
.060
.047
.024
.083
.074
.047
.047
.024
.074
.047
.047
.024
.047
.000
.024
.000
CODE
BOD DECA
BOD SETT
OXY TRAN
SOD RATE
ORGN DEC
ORGN SET
NH3 DECA
NH3 SRCE
N02 DECA
PORG DEC
PORG SET
DISP SRC
ALG GROW
ALG RESP
ALG SETT
COLI DEC
ANC DECA
ANC SETT
ANC SRCE
A very useful feature of QUAL2E pertaining to modifications
of reach identification once the system has been coded is that
existing reaches may be subdivided (or new reaches added) with-
out renumbering the reaches for the whole system. If. for
example, it is desired to divide the river reach originally
designated as REACH 3 into two reaches, the division is made by
calling the upstream portion REACH 3 and the "new reach"
downstream REACH 3.1. Up to nine such divisions can be made
per reach (3.1-3.9); thus REACH 3 (or any other reach) can be
divided into as many as 10 reaches numbered 3, 3.1-3.9. This
option of dividing a reach is useful particularly when new
field data indicate a previously unknown or a change in
geomorphology, or when the addition of a new or proposed load
alters the biochemistry in the downstream portion of the
reach. If this option is invoked, the number of reaches
specified in Data Type 1 must be changed to the new total
number of reaches.
Note: It is important to realize that this option cannot
be used to subdivide a reach into more (and thus smaller)
95
-------
computational elements, in an attempt to provide greater detail
to the simulation. All computational elements must have the
same length (as specified in Type 1 Data).
This option also will allow the user to add a new reach to
the system; for example, taking a tributary that was initially
modeled as a point source and changing it to a modeled reach
(or reaches) in the basin. This type of modification adds a
junction to the system and thus the junction information in
Data Types i. 4. and 9 must be modified accordingly.
This group of cards must end with ENDATA2.
F. Data Type 3 - Flow Augmentation Data
These cards, except ENDATA 3. are required only if flow
augmentation is to be used. The cards in this group contain
data associated with determining flow augmentation requirements
and available sources of flow augmentation. There must be as
many cards in this group as in the reach identification group.
The following information is on each card.
Reach Order or Number Columns 26-30
Augmentation Sources (the number of Columns 36-40
headwater sources which are avail-
able for flow augmentation)
Target Level (minimum allowable Columns 41-50
dissolved oxygen concentration (mg/L)
in this reach)
Order of Sources (order of avail- Columns 51-80
able headwaters, starting at most
upstream points
This card group must end with ENDATA3, even if no flow
augmentation is desired.
G. Data Type 4 - Computational Elements Flag Field Data
This group of cards identifies each type of computational
element in each reach. These data allow the proper form of the
routing equations to be used by the program. There are seven
element types allowed, they are listed below.
96
-------
I FLAG Type
1 Headwater source element.
2 Standard element, incremental inflow/
outflow only.
3 Element on mainstream immediately upstream
of a junction.
4 Junction element.
5 Most downstream element.
6 Input (point source) element.
7 Withdrawal element.
Each card.in this group (one for each reach), contains- the
following information:
Reach Order <5r Number Columns 16-20
Number of Elements in the Reach Columns 26-30
Element Type (these are numbers Columns 41-80
of a set, identifying each
element by type).
Remember that once a system has been coded, reaches can be
divided or new ones added without necessitating the renumbering
of the entire system (see Data Type 2 - Reach Identification
and River Mile/Kilometer Data for application and constraints).
When this option is invoked, the element types and number of
elements per reach for the affected reaches must be adjusted in
Data Type 4 to reflect the changes.
This card group must end with ENDATA4.
H. Data Type 5 - Hydraulics Data
Two options are available to describe the hydraulic
characteristics of the system. The first option utilizes a
functional representation while the second option utilizes a
geometric representation. The option desired is specified in
Data Type 1. card 5. The code "TRAPEZOIDAL" specifically
denotes the geometric representation. Any other code, such as
"NO TRAPEZOIDAL", or "DISCHARGE COEFFICIENTS", specifies the
functional representation.
97
-------
Note: With either option, the effect is global (for the
entire system). This option is not reach variable.
If the first option is selected, velocity is calculated as
V = aQb and depth is found by D =aQ,3 . Each card repre-
sents one reach and contains the values of a. b, a, and 6.
as described below.
Reach Order or Number Columns 16-20
Dispersion Constant Columns 23-30
a. coefficient for velocity Columns 31-40
b. exponent for velocity Columns 41-50
a. coefficient for depth Columns 51-60
|3. exponent for depth Columns 61-70
Mannings "n" for reach (if Columns 71-80
not specified, the program
default value is 0.02)
The dispersion constant is the value of K in the general
expression relating the longitudinal dispersion coefficient to
the depth of flow and shear velocity (10).
DL = Kdu*
where:
DL = longitudinal dispersion coefficient.
(ft2/day. m2/day)
K = dispersion constant, dimensionless
d = mean depth of flow, (ft.m)
u* = shear velocity, (ft/sec, m/sec) = (gdS)-1/2
g = gravitational constant (ft/sec2, m/sec2)
S = slope of the energy grade line (ft/ft, m/m)
Substitution of the Manning equation for S, leads to the
following expression for the longitudinal dispersion
coefficient. DL.
98
-------
DL = 3.82 KnVd5/6
where:
n = Mannings roughness coefficient, and
V = Mean stream velocity (ft/sec, m/sec).
Typical values of K range from 6 to 6000. A value of 5.93 leads
to the Elder equation for longitudinal dispersion, which is the
one used in the SEMCOG version of QUAL-2.
The coefficients a, b. a. and p should be expressed to
relate velocity, depth and discharge units as follows:
System
0
m3/sec
ft3/sec
V
m/sec
ft/sec
D
m
ft
Metric
English
If the second option is selected, each reach is represented as
a trapezoidal channel. This data form is also used to specify
the trapezoidal cross-section (bottom width and side slope), the
channel slope and the Manning's "n" corresponding to the reach.
The program computes the velocity and depth from this data using
Manning's Equation and the Newton-Raphson (iteration) method.
One card must be prepared for each reach as follows:
Reach Order or Number Columns 16-20
Dispersion Constant, K Columns 23-30
Side Slope 1 (run/rise; ft/ft, m/m) Columns 31-40
Side Slope 2 (run/rise; ft/ft, m/m) Columns 41-50
Bottom Width of Channel. Columns 51-60
(feet, meters)
Channel Slope (ft/ft, m/m) Columns 61-70
Mannings "n" (Default = 0.020) Columns 71-80
This group of data must end with an ENDATA5 card.
99
-------
I. Data Type 6 - BOD and DO Reaction Rate Constants Data
This group of cards includes reach information on the BOD
decay rate coefficient and settling rate, sediment oxygen
demand, as well as the method of computing the reaeration
coefficient. Eight options for reaeration coefficient
calculation are available. These are listed below.
K2 OPT Method (9)
1 Read in values of K2.
2 Churchill.
3 O'Connor and Dobbins.
4 Owens and Gibbs.
5 Thackston and Krenkel.
6 Langbien and Durum.
7 Use equation K2 = aQb
8 Tsivoglou-Wallace.
One card is necessary for each reach, and contains the
following information:
Reach Order or Number Columns 16-20
BOD Decay Rate Coefficient (I/day) Columns 21-28
BOD Removal Rate by Settling (I/day) Columns 29-36
Sediment Oxygen Demand Columns 37-44
(gm/ft2-day, gm/m2-day)
Option for K2 (1-8, as above) Columns 45-48
K2 (Option 1 only) Reaeration Columns 49-56
Coefficient, per day, base e. 20C
a. Coefficient for K2 (Option 7) Columns 57-64
or Coefficient for Tsivoglou
(Option 8)
b. Exponent for K2 (Option 7) or Columns 65-72
Slope of the Energy Gradient, Se,
(Option 8)
100
-------
The units of a and b vary depending on whether option 7 or
8 is used and on whether the input data are in English or
Metric units, as follows:
Units of a: English Metric
Option 7 (Coefficient) Consistent with Consistent with
flow in cfs flow in cms
Option 8 (Coefficient) I/ft 1/m
Units of b: English Metric
Option 7 (Exponent) Consistent with Consistent with
flow in cfs flow in cms
Option 8 (Se) Dimensionless Dimensionless
For option 8 (Tsivoglou's option), the energy gradient.
Se need not be specified if a Manning "n" value was assigned
under Hydraulic Data Type 5. Se will be calculated from
Manning's Equation using the wide channel approximation for
hydraulic radius.
This group of cards must end with ENDATA6.
J. Data Type 6A - N and P Coefficients
This group of cards is required if algae, the nitrogen
series (organic nitrogen, ammonia, nitrite, and nitrate), or
the phosphorus series (organic and dissolved) are to be
simulated. Otherwise, they may be omitted. Each card of this
group, one for each reach, contains the following information:
Reach Order or Number Columns 20-24
Rate Coefficient for Organic-N Columns 25-31
Hydrolysis (I/day)
Rate Coefficient for Organic-N Columns 32-38
Settling (I/day")
Rate Coefficient for Ammonia Columns 39-45
Oxidation (I/day)
Benthos Source Rate for Ammonia Columns 46-52
(mg/ft2-day, mg/m2-day)
101
-------
Rate Coefficient for Nitrite Columns 53-59
Oxidation (I/day)
Rate Coefficient for Organic Columns 60-66
Phosphorus Decay (l/day)
Rate Coefficient for Organic Columns 67-73
Phosphorus Settling (I/day)
Benthos Source Rate for Dissolved Columns 74-80
Phosphorus (mg/ft2-day. mg/m2-day)
Note that the benthos source rates are expressed per unit
of bottom area. Other versions of QUAL-2 use values per length
of stream. To convert to the areal rate, divide the length
value by the appropriate stream width.
This group of cards must end with ENDATA6A, even if algae.
nitrogen, or phosphorus are not simulated.
K. Data Type 6B - Algae/Other Coefficients
This group of cards is required if algae, the nitrogen
series, the phosphorus series, coliform, or the arbitrary
non-conservative is to be simulated. Otherwise, they may be
deleted. Each card of the group, one per reach, contains the
following information:
Reach Order or Number Columns 20-24
Chlorophyll a. to Algae Ratio Columns 25-31
(ug chla/mg algae)
Algal Settling Rate (ft/day, m/day) Columns 32-38
Non-Algal Light Extinction Columns 39-45
Coefficient (I/ft, 1/m)
Coliform Decay Coefficient (I/day) Columns 46-52
Arbitrary Non-Conservative Decay Columns 53-59
Coefficient (I/day)
Arbitrary Non-Conservative Settling Columns 60-66
Coefficient (I/day)
Benthos Source Rate for Arbitrary Columns 67-73
Non-Conservative (mg/ft2-day. mg/m2-day)
102
-------
This group of cards must end with ENDATA6B. even if algae.
nitrogen, phosphorus, coliform. or the arbitrary
non-conservative are not simulated.
L. Data Type 7 - Initial Conditions - 1
This card group, one card per reach, establishes the initial
conditions of the system, with respect to temperature, dissolved
oxygen concentration, BOD concentrations, and conservative
minerals. Initial conditions for temperature must always be
specified whether it is simulated or not. The reasons for this
requirement are: (a) when temperature is not simulated, the
initial condition values are used to set the value of the
temperature dependent rate constants; (b) for dynamic
simulations the initial condition for temperature, and every
other quality constituent to be simulated, defines the state of
the system at time zero; and (c) for steady state simulations
of temperature, an initial estimate of the temperature between
freezing and boiling is required to properly initiate the heat
balance computations. Specifying 68°F or 20°C for all
reaches is a sufficient initial condition for the steady state
temperature simulation case. The information contained is as
follows:
Reach Order or Number Columns 20-24
Temperature (F or C) Columns 25-31
Dissolved Oxygen (mg/L) Columns 32-38
BOD (mg/L) Columns 39-45
Conservative Mineral I* Columns 46-52
Conservative Mineral II* Columns 53-59
Conservative Mineral III* Columns 60-66
Arbitrary Non-Conservative* Columns 67-73
Coliform (No./lOO ml) Columns 74-80
* - Units are those specified on the Title Card.
This group of cards must end with ENDATA7.
103
-------
M. Data Type 7A - Initial Conditions - 2
This group of cards is required if algae, the nitrogen
series, or the phosphorus series are to be simulated. The
information is coded as follows:
Reach Order or Number Columns 20-24
Chlorophyll a. (ug/L) Columns 25-31
Organic Nitrogen as N (mg/L) Columns 32-38
Ammonia as N (mg/L) Columns 39-45
Nitrite as N (mg/L) Columns 46-52
Nitrate as N (mg/L) Columns 53-59
Organic Phosphorus as P (mg/L) Columns 60-66
Dissolved Phosphorus as P (mg/L) Columns 67-73
This group of cards must end with ENDATA7A. even if algae,
nitrogen, or phosphorus are not simulated.
N. Pita Type 8 - Incremental Inflow - 1
This group of cards, one per reach, accounts for the
additional flows into the system not represented by point
source inflows or headwaters. These inflows, which are assumed
to be uniformly distributed over the reach, are basically
groundwater inflows and/or distributed surface runoff that can
be assumed to be approximately constant through time.
An important new feature to QUAL2E is that incremental
outflow along a reach may be modeled. This option is useful
when field data show a decreasing flow rate in the downstream
direction indicating a surface flow contribution to groundwater.
Each card, one for each reach, contains the following
information:
Reach Order or Number Columns 20-24
Incremental Inflow (cfs. Columns 25-31
m3/sec) outflows are indicated
with a minus "-" sign.
Temperature (F. C) Columns 32-38
104
-------
Dissolved Oxygen (mg/L) Columns 39-44
BOD (mg/L) Columns 45-50
Conservative Mineral I Columns 51-56
Conservative Mineral II Columns 57-62
Conservative Mineral III Columns 63-68
Arbitrary Non-Conservative Columns 69-74
Coliform (No./lOO ml) Columns 75-80
This group of cards must end with ENDATA8.
O. Data Type 8A - Incremental Inflow - 2
This group of cards is a continuation of Data Type 8. and
is required only if algae, the nitrogen series or the
phosphorus series are to be simulated. Each card, one per
reach, contains the following information:
Reach Order or Number Columns 20-24
Chlorophyll
-------
indicated in Figure A-5; that is. the junctions must be ordered
so that the element numbers lust downstream of the junction are
specified in ascending order. In Figure A-5. the downstream
element numbers for Junction 1. 2 and 3 are 29, 56. and 64,
respectively. There is one card per junction, and the
following information is on each card:
Junction Order or Number Columns 21-25
Junction Names or Identification Columns 35-50
Order Number of the Last Element Columns 56-60
in the reach immediately
upstream of the junction
(see Figure A-5). In the
example, for Junction 1.
the order number of the
last element immediately
upstream of the junction
is number 17. For
Junction 2, it is number
49. For Junction 3. it
is number 43.
Order Number of the First Element Columns 66-70
in the reach immediately down-
stream from the junction.
It is these numbers that must
be arranged in ascending
order. Thus, for Figure A-5
these order numbers are as
follows:
Downstream
Junction Element No.
1 29
2 56
3 64
Order Number of the Last Element Columns 76-80
in the last reach of the tribu-
tary entering the junction.
For Figure A-5 these order
numbers for Junctions 1, 2,
and 3 are 28, 55, and 63,
respectively.
This group of cards must end with ENDATA9. even if there
are no junctions in the system.
106
-------
Most Upstream
Point
Reach
Number
Computational
Element Number
FIGURE A-5. STREAM NETWORK EXAMPLE TO ILLUSTRATE DATA INPUT
107
-------
Q. Data Type 10 - Headwater Sources Data - 1
This group of cards, one per headwater, defines the flow.
temperature, dissolved oxygen. BOD. conservative mineral.
algae, nutrient, coliform, and arbitrary nonconservative
concentrations of the headwater. The following information is
on each card:
Headwater Order or number
Starting at Most Upstream Point
Headwater Name or Identification
Flow (cfs, m3/sec)
Temperature (F. C)
Dissolved Oxygen Concentration (mg/L)
BOD Concentration (mg/L)
Conservative Mineral I
Conservative Mineral II
Conservative Mineral III
This group of cards must end with ENDATA10.
R. Data Type 10A - Headwater Sources Data - 2
This group of cards supplements the information in Data
Type 10. and is required if algae, the nitrogen series, the
phosphorus series, coliform. or arbitrary non-conservative are
to be simulated. Each card, one per headwater, contains the
following data:
Columns 15-19
Columns 20-35
Columns 36-44
Columns 45-50
Columns 51-56
Columns 57-62
Columns 63-68
Columns 69-74
Columns 75-80
Headwater Order or Number
Arbitrary Non-Conservative
Coliform. (No./lOO ml)
Chlorophyll a. (ug/L)
Organic Nitrogen as N (mg/L)
Ammonia as N (mg/L)
Nitrite as N (mg/L)
Columns 16-20
Columns 21-26
Columns 27-32
Columns 33-38
Columns 39-44
Columns 45-50
Columns 51-56
108
-------
Nitrate as N (mg/L)
Organic Phosphorus as P (mg/L)
Dissolved Phosphorus as P (mg/L)
Columns 57-62
Columns 63-68
Columns 69-74
This group of cards must end with ENDATA10A, even if algae.
nitrogen, phosphorus, coliform. or arbitrary non-conservative
are not simulated.
S. Data Type 11 - Point Load - 1
This group of cards is used to define point source inputs
and point withdrawals from the stream system. Point sources
include both wasteloads and unsimulated tributary inflows. One
card is required per inflow or withdrawal. Each card describes
the percent of treatment (for wastewater treatment), inflow or
withdrawal, temperature, and dissolved oxygen, BOD, and
conservative mineral concentrations. They must be ordered
starting at the most upstream point. The following information
is on each card:
Point Load Order or Number
Point Load Identification or Name
Percent Treatment (use only if in-
fluent BOD values are used)
Point Load Inflow or Withdrawal
(cfs. m3/sec) (a withdrawal must
have a minus {"-") sign
Temperature (F. C)
Dissolved Oxygen Concentration (mg/L)
BOD Concentration (mg/L)
Conservative Mineral I
Conservative Mineral II
Conservative Mineral III
Columns 15-19
Columns 20-31
Columns 32-36
Columns 37-44
Columns 45-50
Columns 51-56
Columns 57-62
Columns 63-68
Columns 69-74
Columns 75-80
This group of cards must end with ENDATA11,
109
-------
T. Data Type 11A - Point Load - 2
This group of cards supplements Data Type 11 and contains
the algal, nutrient, coliform and arbitrary non-conservative
concentrations of the point source loads. This information is
necessary only if algae, the nitrogen series, the phosphorus
series, coliform, or the arbitrary non-conservative are to be
simulated. Each card, one per waste load (withdrawal) contains
the following information:
Point Load Order or Number Columns 16-20
Arbitrary Non-Conservative Columns 21-26
Coliform (No./lOO ml) Columns 27-32
Chlorophyll a. (ug/L) Columns 33-38
Organic Nitrate as N (mg/L) Columns 39-44
Ammonia as N (mg/L) Columns 45-50
Nitrite as N (mg/L) Columns 51-56
Nitrate as N (mg/L) Columns 57-62
Organic Phosphorus as P (mg/L) Columns 63-68
Dissolved Phosphorus as P (mg/L) Columns 69-74
This group of cards must end with ENDATA11A, even if algae,
nitrogen, phosphorus, coliform, or arbitrary non-conservative
are not simulated.
U. Data Type 12 - Dam Reaeration
This group of cards is required if oxygen input from re-
aeration over dams is to be modeled as a component of the dis-
solved oxygen simulation. Dam reaeration effects are estimated
from the empirical equation attributed to Gameson (9). The
following inputs are required:
Dam Order or Number Columns 20-24
Reach Number of Dam Columns 25-30
Element Number Below Dam Columns 31-36
110
-------
ADAM Coefficient: Columns 37-42
ADAM = 1.25 for clear to
slightly polluted waters
ADAM = 1.00 for polluted water
BDAM Coefficient: Columns 43-48
BDAM = 1.00 for weir with
free fall
BDAM = 1.30 for step weirs or
cascades
Percent of Flow Over Dam Columns 49-54
(as a fraction 0.0-1.0)
Height of Dam (ft, m) Columns 55-60
This group of cards must end with ENDATA12. even if oxygen
input from dam reaeration is not to be modeled.
V. Data Type 13 - Downstream Boundary - 1
This data card supplies the constituent concentrations at
the downstream boundary of the system. It is required only if
specified in Data Type 1. card 8. This feature of QUAL2E is
useful in modeling systems with large dispersion in the lower
reaches (e.g.. estuaries). When downstream boundary
concentrations are supplied the solution generated by QUAL2E
will be constrained by this boundary condition. If the
concentrations are not provided, the constituent concentrations
in the most downstream element will be computed in the normal
fashion using the zero gradient assumption (I. 11).
Downstream boundary values for temperature, dissolved
oxygen. BOD, conservative mineral, coliform, and arbitrary
non-conservative are required as follows:
Temperature (F, C) Columns 25-31
Dissolved Oxygen (mg/L) Columns 32-38
BOD Concentration (mg/L) Columns 39-45
Conservative Mineral I Columns 46-52
Conservative Mineral II Columns 53-59
Conservative Mineral III Columns 60-66
111
-------
Arbitrary Non-Conservative Columns 67-73
Coliform (No./lOO ml) Columns 74-80
This data group must end with an ENDATA13 card, even if the
fixed downstream boundary concentration option is not used in
the simulation.
W. Data Type 13A - Downstream Boundary - 2
This group of data (one card) is a continuation of Data
Type 13. It is required only if the fixed downstream boundary
condition is used and if algae, the nitrogen series, the
phosphorus series are to be simulated. This card contains the
downstream boundary concentrations for algae, nitrogen, and
phosphorus as follows:
Chlorophyll a. (ug/L) Columns 25-31
Organic Nitrogen as N (mg/L) Columns 32-38
Ammonia as N (mg/L) Columns 39-45
Nitrite as N (mg/L) Columns 46-52
Nitrate as N (mg/L) Columns 53-59
Organic Phosphorus as P (mg/L) Columns 60-66
Dissolved Phosphorus as P (mg/L) Columns 67-73
This data group must end with an ENDATA13A card, even if
the fixed downstream boundary condition is not used, and if
algae, nitrogen, or phosphorus are not simulated.
X. Climatological Data
Climatological data are required for the following cases:
1. Temperature simulations, both steady-state and dynamic,
2. Dynamic simulations where algae is being simulated, and
temperature is not.
If neither temperature nor dynamic algae are being simulated,
these cards may be omitted.
112
-------
For steady-state temperature simulation, only one card is
required which gives average values of the climatological
data. For dynamic simulation, each card represents readings at
three hour intervals, chronologically ordered. There must be a
sufficient number of cards to cover the time period specified
for the simulation (Data Type 1, card 13, maximum route time).
The following information is on each card.
Month Columns 18-19
Day Columns 21-22
Year (last two digits) Columns 24-25
Hour of Day Columns 26-30
Net Solar Radiation* Columns 31-40
(BTU/ft2-hr. Langleys/hour)
Cloudiness**, fraction in Columns 41-48
tenths of cloud cover
Dry Bulb Temperature** (F, C) Columns 49-56
Wet Bulb Temperature** (F, C) Columns 57-64
Barometric pressure Columns 65-72
(inches Hg. millibars)
Wind speed** (ft/sec, m/sec) Columns 73-80
Required only if dynamic algae is simulated and
temperature is not.
Required if temperature is simulated.
There is no end card for the climatological data.
Y. Plot Reach Data
This data type is required if the plotting option for
DO/BOD is selected (Data Type 1, card 7. PLOT DO/BOD). The
following information is required for QUAL2E to produce a line
printer plot.
113
-------
1. Card I - BEGIN RCH
Reach number at which plot Columns 11-15
is to begin
2. Card 2 - PLOT RCH
a. Reach numbers in their Columns 11-15
input order (1, 2. 3..NREACH) Columns 16-20
21-26
b. If a reach is not to be etc.
plotted, (i.e.. a tributary) 76-80
replace the reach number
with a zero.
c. Use additional PLOT RCH cards
if there are more than 14
reaches in the system.
3. Additional plots can be obtained by repeating the
seguence of BEGIN RCH and PLOT RCH cards.
As an example of the plotting option, suppose that for the
river system shown in Figure A-5, one.wishes to obtain two
DO/BOD plots: one for the main stream (Reaches 1, 2, 5. 6. 10,
and 11) and one for the second tributary (Reaches 7 and 9).
The plot data would appear in the following order.
BEGIN RCH 1
PLOT RCH 120056000 10 11
BEGIN RCH 7
PLOT RCH 00000070900
No ENDATA card is required for the PLOT information.
Z. Summary
Constructing a consistent and correct input data set for a
QUAL2E simulation must be done with care. This user's guide is
designed to assist the user in this process. It has been
NCASI's and EPA's experience that two of the most frequently
made errors in constructing a QUAL2E input data set are:
(a) Using a numerical value that is inconsistent with the
units option selected, and
(b) Not adhering to the 4 character input codes for Data
Types 1 and 1A.
114
-------
As an aid to the units problem. Table A-l is included in
this report. It provides a complete summary of all the input
variables whose dimensions are dependent on whether English or
metric units are selected. Finally, the user is encouraged to
check and recheck the input codes in Data Types 1 and 1A for
accuracy, especially the codes for cards 10 and 11 of Data Type
1 (i.e.. "NUMB" and "NUM_").
REFERENCES - APPENDIX A
1. "A Review of the Mathematical Water Quality Model QUAL-2
and Guidance for its Use", NCASI Technical Bulletin No.
338. October 1980 (Revised No. 391. December. 1982).
2. "Computer Program Documentation for the Stream Quality
Model QUAL-2". prepared for Southeast Michigan Council
of Governments (SEMCOG). Detroit, Michigan (July, 1977).
3. "Users Manual for the Stream Quality Model QUAL-2".
prepared for Southeast Michigan Council of Governments
(SEMCOG). Detroit. Michigan (July. 1977).
4. "QUAL-TX User's Manual (Draft)". Texas Water Development
Board. Austin. Texas (June. 1981).
5. "Calibration and Application of QUAL-2 to the Lower Win-
ooski River: Preliminary Studies", prepared for state
of Vermont by Meta Systems. Inc. (July. 1979).
6. Norton. W.R.. et al. "Computer Program Documentation for
the Stream Quality Model-QUAL-2". for EPA Contract 68-
01-1869 (August. 1974).
7. Patterson. D.. e_t a_l, "Water Pollution Investigation
Lower Green Bay and Lower Fox River". Wisconsin DNR.
(EPA 68-01-1572). (June. 1975).
8. "A Study of the Selection, Calibration, and Verification
of Mathematical Water Quality Models". NCASI Technical
Bulletin No. 367. (March 1982).
9. "Rates, Constants, and Kinetic Formulations in Surface
Water Quality Modeling". USEPA. Athens. Georgia. NTIS
PB-290938 (March. 1978).
10. Fisher, H.B.. E.J. List. R.C.Y. Koh, J. Imberger. N.H.
Brooks, Mixing in Inland and Coastal Waters. Academic
Press, 1979.
115
-------
TABLE A-l.
LIST OF QUAL2E INPUT VARIABLES THAT ARE
ENGLISH/METRIC UNIT DEPENDENT
CARD
DATA OR
TYPB klM
1 8
8
1 11
1 IS
15
1 16
1A 6
6
1A 7
1A 9
2 all
11
5 all
(Discharge
Coefficient)
S all
(Trapezoidal)
6 all
6 all
6 all
6A all
CB all
7 all
B all
10 all
11 all
12 all
13 1
LCD all
VARIABLE DESCRIPTION
Input Units specification
Output Units Specification
Length of Computational Element
Evaporation coeftlcient
Evaporation Coefficient
Ban In Elevation
Linear Algal Extinction Coeff.
Non-linear Algal Extinction
Coefficient
Light Saturation Coefficient
Total Daily Solar Radiation
River Mile/km to Head of Reach
River Mile/da to End of Reach
Coefficient Flow for Velocity
Exponent Flow tot Velocity
Coeftlcient on Flow for Depth
Exponent on Flow for Depth
Bottom Width of Channel
SOD Rate
Option 7 for kj
Coefficient on flow for Kj
Exponent on flow for k2
Option 8 for Kj
Coefficient for Tsivoglou Eq.
Slope of Energy Gradient
Benthal Source Rate for
An>monla-N
Benthal Source Rate for
Phosphorus
Algal settling Rate
Non-algal Extinction Coefficient
Arbitrary Nonconeervative
Benthal Source Rate
Initial Condition-Temperature
Incremental Inflow
Flow Rate
Temperature
Headwater Conditions
Flow Rate
Temperature
Point Source/Withdrawal
Flow Rate
Temperature
Height of Dam
Downstream Boundary-Temperature
Solar Radiation
Dry Bulb Temperature
Het Bulb Temperature
Barometric pressure
Wind Speed
FORTRAN
CODE NAH8
METRIC
METOUT
DELX
AE
BE
ELEV
EXALG1
EXALG2
CKL
SONET
RMTHOR
RMTEOR
COEFQV
Exrogv
COEFOH
EXPOQH
WIDTH
CM
COEQK2
EXPQK2
COEQK2
EXPQK2
GNH3
SPHOS
ALOSET
EXCOEF
8RCANC
TINIT
01
TI
HWKLOH
HHTEHP
HSFLOH
HFTEHP
HDAH
LBTEHP
6OI.HR
DRYBI.B
METBLB
ATMPR
HIND
UNITS
ENGLISH
O
0
mile
ft/hr-ln Hq
ft/hr-in Hg-mph
ft
1/ft-ug chla/L
l/ft-(ug chla/L2/3
BTU/tt2
BTU/ft2
mile
mile
Consistent with
flow, velocity
and depth in
cts. fps. ft
respectively
ft
gm/ft2-day
Consistent with
flow in cfs
I/ft
ft/ft
Mg/ft2-day
mg/ft2-day
ft/day
I/day
»g/ft2-day
F
cfs
F
cfs
F
cfs
F
ft
F
BTO/ft2-hr
F
F
in Hg
ft/sec
METRIC
1
1
kilometer
m/hr-mbar
m/hr-mbar-m/sec
meters
1/m-ug chla/L
l/m-(ug chla/L2/3)
langley/mln
langleys
kilometer
kilometer
Consistent with
flow, velocity, and
depth in cms, mps.
respectively
meters
gm/m2-day
Consistent with
flow in cms
I/meter
meter/mater
mg/m2-day
mg/mz-day
m/day
I/meter
g/m2-day
C
cms
C
cms
C
cms
C
meters
C
langleys/hr
C
C
mbar
/sec
116
-------
11. Arden. B.W. and K.N. Astill. Numerical Algorithms: Ori-
gins and Applications, Addison-Wesley. 1970.
12. Department of Scientific and Industrial Research. Ef-
fects of Polluting Discharge on the Thames Estuary. Wa-
ter Pollution Research Technical Paper No. 11. Her Maj-
esty's Stationery Office. London. 1964.
13. Scavia. D. and R.A. Park, "Documentation at Selected
Constructs and Parameter Values in the Aquatic Model
CLEANER". Ecological Modelling. Vol. 2. pp. 33-58, 1976.
14. Bannister. T.T., "Production Equations in Terms of Chlo-
rophyll Concentration. Quantum Yield, and Upper Limit to
Production". Limnology and Oceanography. Vol. 19, No. 1.
pp. 1-12 (Jan., 1974).
117
-------
EH
H
cn
EH
<
Q
EH
D
ft
EH
J3
O
\
EH
D
Oi
2
M
X
W
CQ
X
H
Q
CO
s
or
o
Q
O
O
O
H
^
Q.
UJ
OJ
_J
O
I
_l
UJ
o
o
o
O
o:
UJ
a:
co
co
<
o
Q.
UJ
^
UJ
_J
1-
1
(*)
«.
z
1?
co
2
\
\
UJ >
1- Ul
UJ ^
s <
Q;
s
/
UJ
2
<
^
O
UJ
13
Q.
_J
V
\^
< 2
o to
SUJ
r UJ
Q.
or
o
tt
£
00
h-
fc
*
K
£
K
f!
£
F:
£
(-
or
Ul
CO
z
o
o
Kl
O
111
_l
-
_l
<
a:
Ul
z
2
Ul
>
_
l-
<
>
cc
Ul
CO
z
o
0
*
o
UJ
-1
_l
<
cc.
Ul
z
2
UJ
>
_
(-
<
>
a:
UJ
CO
z
0
u
in
o
UJ
-J
Ul
a:
3
(_
<
a:
Ul
a.
2
Ul
h-
CD
O
UJ
_l
O
z
<
2
Ul
O
z
UJ
cs
>.
X
o
_l
<
u
_
2
UJ
X
o
o
CD
f-
o
Ul
_J
_l
\
CD
3
Z
<
1
_l
X
O
CO
<
Ul
<
CD
_J
<
CD
O
Ul
_l
_l
v^
CD
2
Z
_
a.
CO
<
UJ
_i
o
>-
o
co
3
tr
0
X
a.
CO
0
X
a.
O)
O
m
_i
a.
i
a
UJ
>
_i
o
CO
CO
a
..
a.
i
CJ
z
«t
U)
tr
o
O
111
_J
_l
X.
CD
2
Z
Z
CO
<
UJ
_l
o
x
o
z
UJ
CD
O
tr
i-
z
UJ
_)
z
1
UJ
t-
tr
i
z
..
z
1
Ul
h-
cc
1-
z
..
z
1
<
z
o
2
2
<
-
z
1
o
z
X
o
a
UJ
>
-J
o
CO
CO
o
f)
Ul
_J
_
2
o
o
\
o
z
z
CO
2
tr
o
u.
_l
0
o
_l
<
o
UJ
u.
«
Ul
_>
L
>
1-
<
>
cc
Ul
CO
z
0
0
1
z
o
z
>-
tr
A
s
fe
3!
S
s
R
Si
in
g
ffi
(0
$
"
n
J
10
nj
;
^
S
D
£
u>
(1
s
o
s
n
S
M
S
s
R
f
V
S
N
M
A
5
t
(0
£
S
2?
u
=
£
9
0
-
<0
A
»
to
a-
o
O"
0
>
0
3
.0
f
0
+_
3
E
M
A
O
>
*
U
u
3
O
f
Q,
O
f
a.
j=
c
w
o
2
^
*
u
~ £
* o»
< o
CVJ C
"^
5-
* £
0 0
u. z
118
-------
C/)
^
or
o
o
o
o
UJ
CM
Z)
o
LU
Q
O
O
Z
O
or
UJ
or
CO
<
o
2
UJ
o
Q:
O
O
or
CO
o
10
a
£
e
u t
< *"
Q: «>
LJ 8
(- £
LU ^
2U>
< «
or s
2 =
O *°
g s
|_ K
"Z.
CJ
^E
f-r
K
CD
or
Q.
UJ S
3 -
rf "
> s
K ~ v
l»J s a:
fc s H
g 9 <0<
to . ^
O - =-"z
Q: E Q. < UJ
t- ~ ZZ S UJ
§ r-OOH
8 £ -=><
» « o o>-
0 w_.o
£ *HHU.<
n (ft Itl
a troi-
- _ii z
i M u H u
~
x en o
U. 0 KUU
uj~ to i o:
0 0 -
- - SO-OTJO
> I- H uj o:
023
H I tcio: IH u
-J'3 UJ UJ C
ZJIQ. OD CD XIUJIU
i KlSlZKlZZ >
Q Z) 3 3Z' O <
im o z:z _j H -ilc
n i > i
E
-u -xujar < ZHK < HOS
CO UJ Z
oxa. S2Sx»-<
ee _i- z z> = - « H
a. Q.U. ZZI-SJOT
e
J £
E u. S
UJ £
-oS
o S
ID
c z ;
tui o S
-IOD - s
Oi**. t
^ -
E = S
lu. Z £
Ju UJ S
-0 K!
io H
<
Q. h-' ?
;s
31UJIO
i S
II H
> S
UJ S
_i in
-UJ t
UJ- S
< 5
z. s
.CO f]
u.< =
U. CD 2
UJ a
OU. »
0 O - -
< (0
. H *>
Q. >< !
< uo » :
>-.z~ 1
UJ UJ UJ - i
I I
= O
S E J
E
£
E
E 3
I
10
- o «
_!
(9
< j, uj
< a
to a.
. E
X »
o ».
o
2 I
E £ £ £
UJ UJ
< CO
ISH
^ 11:
° ? s §
pill
u. *
****
119
-------
CO
s
o:
o
u.
o
Q
O
O
I
0
I-
Z)
a.
LJ
(VI
o
UJ
Q
O
O
O
cr
*
en
cr
UJ
i-
UJ
IT
S
Q_
Q
m
0
UJ
a.
UJ
cr
H-
cn
eo
o
a.
UJ
I
VALUE
PARAMETER
UJ
Z>
_l
%
PARAMETER
II
z
(9
Z
\
o
o
z
o
X
o
CM
Q
Z
l>-
DD
UJ
n
<
O
s
\
o
o
s
«-J
UJ
<
o
_l
-
ffl
ii
~
<
0
5
\
0.
(9
s
1
UJ
<
II
.
>
<
o
^
~-
UJ
H
<
tr
z
190
_l_
4lh-
<
i_
0
cc.
a.
1-105
ujlziuj
? *l*juj|a:
o
t
10
«
«r
2
~
r
o
01
a
t^
to
m
t
10
-
£
\
O
o
s
o
X
o
1C
X
z
>-
CO
UJ
*
4
t-
a.
3
O
D
O
n
_i
\
II
1C
\
UJICM
Z
_
1-
tn
z
*
*
_i
v
<
X
00
(J
z
o
e>
D
^
I
II
~
z
£
V.
I
Z
_
_
u.
UJ
o
o
z
xio
"KK»
<
ec
z>
(-UJ
-~
<-l
ok
viol
-1
UJ
t-
<
Z
t-
(/)
itz
I
h-
s
o
K
(9
O
UJ
CL
to
X
<
2
t9
_j
<
O
o
z
0
_
1-
<
CE
Z)
1-
<
to
u.
_i
<
X
z
_l
_
1-
<
a:
3
II
_
t-
tn
i-
X
C9
_l
~
t-
N.Q.
<
X
0
o
D
I
X
\
0
u
UJ
o
<
X
(0
(9
_J
<
Z
_l
4
K
X
0
_J
cr
n
z
u.
UJ
a:
a.
^
n
U.
Ul
o
o
z
o
1-
z|-
.|a>
IOI
X X
<|zz
_l
o
in
>-
_i
_
4
O
_1
4
1-
O
I-
H
Q-IX
O
0|>
Z
u.
_i
^
z
o
K
Q.
0
z
o
1
(J
z
3
U.
1-
X
o
-J
4
_l
Z
O
t-
Q.
O
e
z
_
(9
<
K
UJ
>
<
>-
_!
_
4
0
_l
a
^.
a:
3
O
X
1-
X
0
_l
>-
4
a
u.
0
or
UJ
m
s
Z)
z
!
a:
oiz
u.lo
u.|i-
LJI4
cclo
CL
U.
_J.
40:
C9K
_1
4IZ
-)
Kr
0-1 (-)
0
cc
o
_l
~-
z
o
1-
Q-
o
u
4
U.
1-
cc
o
K
o
4
u.
a
4
jce
4
O
X
t-
$
o
cc
w
>-
_J
4
cc
_i
o
S
IS
s
tf>
fl-
lfi
3
x>
Jl
n
n
In
m
K
iTi
A
O
T
OD
^
T
(C
V
if}
T
V
V
T
t
S1
£
1"
4
4
K
4
a
z
UJ
3
£
u?
8
^
10
5o
s
8
ffi
(NJ
M
W
in
«
Kl
s
w
y
2
(D
«
V
K)
N
-
0
£
a
1
a.
o
c"
o
'c
o
o>
0
tn
t-
z
o
K
UJ
^
0
£.
O
\
E
G
o
JC
u
3
4
c
«n
i
C
o
_J
»
a>
c
o
E
<
O o
t 5
O .2 £
o r
'Ice
^52
^ S 5
*; o:
-31 2
o « £^
u.
4
00
4
S
or
0
u.
o
£
4
t-
o
o
o
*
a.
ii
X
c
o
Z
'c
3
CC
UJ
tu
4
tc.
2
«>
"
t/>
1
4
5
'J
£
tn
tn
o
o>
4
o
1
z
£
3
§
|
in
3
0>
J
o
o
IT
o
0
tn
I
O
h-
120
-------
tf\
^
2
oc
o
LL
^
O
o
o
<
<
Q
H
z>
O.
z.
UJ in
C\J UJ
-J 3
< <
=> §
0 *r
. 5
UJ
-I I
UJ H-
Q
O
2 m
& LU
Q_
ri >
t 1-
_j
/-» «r
§ £
a <
i- 0
>
r
_l
<
13
^
1 1 1
fy
frt
CO
^
^r f9\
^ \fi/
CL
UJ "S
5
cc
£
LJ
:D
H
a
P
g
a
ffi
£
t
£
£
v
Kl
£
=
0
ID
h-
ID
if)
V
rt
-
^^ pX1*1
X
I
4
jJ
r
s
tv
S
0
a
t
<£
tf-
t
£
<\
=
s
01
Q 00
s
< u
t- o
-------
CO
2
cc
o
Q
O
O
O
H
D
CL
UJ
CM
a t
UJ
Q
O
O
z
I-
o
cr
O
UJ
cc
O)
CO
-------
en
S
o:
O
Q
O
O
<
Q
ID
CL
LU
CM
ID
O
LU
Q
O
O
z
o
cr
ID
O
LU
tr
CO
en
<
o
Q.
LU
1
<
o
n
f
£
_>
J
%
o
LJ_
to
UJ
n
>-
i-
|^
O
"Z.
<
1-
LJ
S
CP
D
LJ
m
<
_i
<
<
u.
o
or
UJ
Q
o
2
d^
^>
y?
. r-
m£
£-
^[7
Hc^
C
1
. "S.
LJ i.
Qn
or w
LJ
a.
h-
2
^
o
§
^ 00 1 !
^ 1
S i 1
s
»-
TJ "
^~ 2 '
3 j ' '
c S '
CVJ
tf) \
\
^ *
> )
JO
S ;
;
CO ? ,
LJ _J «,
oc§ B
?=! "
en
=
S 1-
0 S
oo e>
" ^
1C <
0
» 5
> O
~ _i
- u.
: ;
i
i
-J- T- ~r
i
: ' !
1
j
's ' , -
|
-* '. fc
S
01
as
r-
Is
K
,
]
~
\ f3
£
]
;
to
ID
i
;
j
j
i
i i
S
S
In
i
J
S
I
H
K
£
*
?
^
u
£
$
5
v
;
f
0
a
t
u
m
10
5
S
S
s
S
s
cc
w
s
s
I S
u "
tr a
N
8
O)
to s
UJ t
o s
CC E
O »
O E
CO S
I- S
2 c.
(5 CD
3 fO r-
< < »
t f
S < »
0 0 «
_1 2 IM
U- UJ -
-------
CO
or
o
Q
O
O
I-
Q
D
Q.
UJ
C\J
13
O
UJ
Q
O
Z
H
O
or
UJ
QL
h-
co
CO
o
Q.
UJ
Q
Q
LJ
LL.
0
LJ
UJ
_l
LJ
Z
O
s
o
o
LJ
a.
s
©
I
0
OJ
CT>
00
(^
CO
in
< °°
£ 1-
Q.
5 ,0
o ^
LJ_ .
0|5
2 UU
ZO UJ
UJ,, O
QO<
o t
LJ
Q_
V
(-
£
<
00
M
X
o>
£
f»
£
£
S
n
M
^
9
91
S.
«
m
<
_i
u.
V
4
t-
*t
a
z
UJ
g
M
R
y
m
£
f~
tt
A
J
IO
~
=
£
O)
t-
10
n
»
K)
CM
~
)
f
c
«>
g
3
C
"3
i
"5
>
o
V
A <
_ h-
0 <
n ^
' ^
1 s
c £
8 >
0 (C
li
E "*
1 z
n O
» 5
P~ g
~IP
b. 5
O «»
~Z o
x IS
°'° 5
- » i»
/^ « ^
O > uj
»' ° K
in » "
u. c t
x"?
£ s =
o" ° 2
to «> .,
"-- s
x ° "
*> § 6
>-
« '
-------
CO
S
cr
o
o
z
o
o
o
o
I-
o
Q.
z
UJ
C\J
o
I
_l
Ul
Q
O
O
cc
UJ
en
i-
co
CO
<
o
2
UJ
-
T
ID
Ll1
n
>
h-
I
X CC
°b
z o
O -1
U.UJ
Ld >
Of\
U.
z>-
o.sy
iJio
QO
(V -p
LUi. O
Qo<
CO
idl-
es 2
tf r\
i-
o H-
^m
°o
o
L
a
h
m
V
to
«j
-
,
'
X
O
a:
tn
o
_j
u
<
a:
o
>
X
!
|
'
[
1
1
i
1
1
1
!
»v^
j
,
[
'
1
1
^
I
i
I
-.
i
T
i
I
1
1
,
j
S
!
i
|
|
1
| |
I i
1
;
I
, J
1
i
I
1
!
j
*^. _
i '
i
!
'
>
t
1
!
s11
\
\
~^.
i
i
s
\
\
1
1
1
^
'
1
1
*«.
! 1
|
i
j
|
1
| [
I
L |
!
i
i
-
!
I
i
'
i
i
i
i
t
i
1
i
i
n
X
o
rr
o
_
_i
r>
<
IT
O
>l
X
\
j
in
<
l-
<
o
z
UJ
8
1
"
!S
~
Sci
S
S
t*-
VI
£
»
E
S
-
Q
Or
a>
t*-
_
§
»
«
V
0 TJ
(O
in w
- *r
o ?
CO *
t «
S 1
0 «
in £
u. -G
x £
10 8
OJ o>
< ?
* °
T u
< »>
PJ T>
^ «>
in r
Q
125
-------
CO
IV
u_
\3
Z
Q
O
o
<
\-
^
Q
\-
*~>
Q.
Z CO
~~ 1
Z
CM O
u.
^ H-
2 ^
O 0
O
o
-J Q
UJ "
9 °
O 0
2 co
O
Z «>
l_ UJ
ID ^
f"\
£ K
a <
,. H
£ <
r- o
_J
^
5
1 1 1
fy
fjO
^
)
Q.
UJ "5
5
cr
£
in
3
O
* °°
M°- o
oO
fecOQ
UJ U_
*
OLJ co
C\J ^ I
^^ Q.
OKU.
*
o ~~
< m*1
CCO°
UJ ~~~
o:
C\j j~
* o
Q JjJ
C/; ^
o:
(_ ^_
i j_
C/5UJ**
f-vO
o
DO
X
QUJ**
QOQ
QO^
m
UJ^O
o ^
UJ
Q.
1
^
1-
a
£
5
t
!£
£
S
m
~
-
g
0,
o>
|~
»
«i
*
10
-
,1
x.
o
a:
u.
UJ
o
o
H
U
<
UJ
cc
>
1X^
^
,
N
-^
S*
_ '
*V^
1
__
-
ii
X
o
tc
u.
LJ
O
o
(.-
u
<
111
c:
(0
<
t-
<
Q
Z
Ul
g
s
CJ
(T)
£
t
l£
12
V
2
~
=
S
0.
eo
-
u>
m
»
«
OJ
-
U
w
J"
;>
E
K -S
1 *
Q:
K >. *
5 ? I
(X) *-
E o \ p
V
^
in "
I- c
§ §
~- i *"
O V) £
» -1 >. 't
U. ID O
ro ?, ? c
W CM «
v '» j J
u. g, 8 c
o
u.
o
u, _
x" ^
« ^_
., O 00
evi * .
< UJ CM a
*" ^
a: o ui > o
U- {/) O ^~ W
126
-------
CO
s
oc
£
o
2
o
o
o
<
Q
Z)
QL
Z
UJ
LU O
C\J u.
_J LL.
< UJ
Is
a-
O z
o
z
I-
o
tr
o
LjJ
CC
CO
CO
<
o
CO
<
Q
ILJ
Uj
o cc
«
s
UJ
CL
>
s
s J
Q
»
' I
0
IT
i
U.
UJ
O
o
o.
o
z
<
z
^.x
1 !
I 1
1
J
1
{
' *" "^ -* ""*
1
1
1
'
"^* "^- "*"" "^* ^>« _- **^ "*" "*^ X
t
S
"
S
K
i
1
]
<
£
!
u
I
S
*
5
c
S
1
K
S
£
j
E
K
!
§
?
?
$
5
in
v
*
t
?
?
^
5
tr
a
fe
S
s
V
n
Jf
R
g
S
a
w
*
n
»
N
a
(U
S
ii 2
I 2
O t
IE - S
V U-
2 cT
UJ -11.
O 2 x
0 « »
< M
0- T
I CM ~
r f
, Q.
?<2
Z Q
o o
o
in
(D OD
127
-------
cn
2
o:
£
O
O
i-
Z)
a.
UJ
CVJ
o
I
UJ
a
o
O
a:
ID
O
UJ
cc
h-
co
CO
<
o
Q_
UJ
c/)
I-
UJ
o
u.
u.
UJ
o
o
tr
ui
X
t-
o
^
UJ
00
CD
UJ
a_
1- I-
<
a
.._ -o
0 «>
3
_)
< LU
fee"
z=>i
UJQ4
*
1-
1 ^^
n^ "*^
LL) u. >-
1 QQ
2 ~_
<^
,
O < tr
ZC>S
<£°
<-> o ^
C^
_1|_
rn X
_) UJ U.
"^ \ o
5^O
* '
(
I
to ^
>- <
oS
^
*
UJ
0
^
W
n
ii
X
o
a:
u.
UJ
o
o
IT
UJ
I
o
\
0
_J
4
ffl
<0
<
I-
4
0
Z
UJ
s
SI
z
ti
i!
s
t
IP
£
J
5
~
=
y
0>
CD
>-
1C
fl
»
m
~
-
O
f-
Q
IT)
U.
X
in
CO
*
^>
Ir>
»
E
CC
u.
C/)
t
z
e o
o i
-1 S
i o <"
s
4 Z ffi
128
-------
en
2
o:
o
CD
z
Q
O
O
ID
Q.
UJ
CVJ
O
UJ
Q
O
ID
O
o:
LU
CC
h-
co
en
<
o
2
UJ
I
0)
z
o
o
o
UJ
Q.
<
I-
Q
ir^
o~^
o
o
z
ro
CO
z
o
o
CVJ
CO
z
o
o
CO
o
o
cc
UJ
£L
CO
o"
X
10
zz
<
t
o
129
-------
5
t
£
£
S
2
«
-
10
«
»
«
rj
-
O
N
n
X
cJ
-------
CO
5
cc
o
u.
Q
o
o
^
1
^
O
1-
ID
Q- .
Z >
^E
_l
UJ U.
CM
1
/rj 00
UJ
1- £
ID I-
0
cc <
1
>- Q
1- °
'
^
«r
^
^
LJ
DC
H-
CO
1
CO
^
o
z
Q- _
UJ °
VCX
s
ir
o
LL.
- i o
or "^
f "^
; 2 O K
O' '^ *
(j .n
(J £
ro 5
CO
O '
o 5
W z
CO o
O "
1~
Ift
1
s
CO
o
?
Q x! ?
mC3 S
S if
^ s
Q _ I I
->LjJ ' v
O> O 5
C/)Q^[|
Q 1
Q_ *" _
5 * .
1- S
^ > °
rr O ~
O.-J S
Z £
* j»
ft' T w
UJQ.0 S _
0 ^ S
»
S D
*^ <
« O
E
UJ ~ -
Q. E -
>- »;
*~ =3
< e c
Q « U
K 2
<0 -
0
» 0
IO {_
(M 2
1
j
1
C
3
i
>
J
3
j
r^-
""»-.-
,
1
1
!
1
-*^ *
i
\ '
i
S.-J" ~S_
*' *
1
v.---
'
1
1
I
C
£C
i
$
p> c
_
U.
z
o:
o
z
-
g
?
;
i
10 "*
ID K ^" ^
ri °
.< S , . c
H-"~_| 03
< « H- ffi C 0
°- = I g s
2~ S < S 6
u.-£ > Sf
131
-------
CO
^
o:
o
Q
O
O
0
\-
^>
a.
LU
OJ
O
i
_l
UJ
Q
O
ID
O
o:
UJ
cr
I-
co
(O
<
o
a
LU
CM
i
O
_1
<
UJ
2
UJ
o:
CO
UJ
Q.
<£
Q
Ol
o
c
Q. _1
i \
-------
fjft
LL
.-
o
o
o
<
H
^
O
t-
^
Q.
Z
_»
LJ
CVJ
1
D
"i
>
n
j>
jj
K
-t
-v
(
c
V
UJ
>.
I-
^
UJ
a
>-
i
ro
n
I-
Al)
0
i?>
^ }
5
t£
U
1- g 1 "
o^ ^ 55 «:
ouj fE 1 :
^UJ i
!
fe ^ o
^~ m K 2
S m ^ '
y "( Z U- £
^GJ a o
t ^ o ;
^ir Jj o; z
^ to ^
§ LJJ Q_ U_
^~ LU ^5 U
_ j
O S
t- »
< ₯
LL.
1 w
UJ "
Q |_S
~ S
Z S
O fe
H S
Z5 "
"* 5 o
Z
IB -s
Z S
uJ , m s
Q^o R
O , w
< «, -J
« 2
« <
^
< u
^ tfl
< «
0 «
Z ~
UJ -
S
o
3 u-
3 X
i m
3 ro
r>
-------
tr
o
o
6
o
o
6
o
\-
D
a.
HI
CO
UJ
o
O
O
01
O
UJ
o:
en
o
UJ
*
il
CD
E
&
134
-------
C/J
2
OL.
O
LL.
ff.
O
O
O
P~
1
Z)
Q.
"Z.
i i N
LU ,
CVJ £t
_l Ul
^ l"~
0 >
0 ^
o
1 ^
It 1
1 x
ll 1
o <
2 °
O UJ
"Z. Q.
^ 1-
S
" K
r- ' ' '
^_ i r- ' ,
g is ! ! ' '
' t2
i *
a"~\ " i i
C/5 O f:
Q _5 0 1 ' I '
£ ' '
S
Q. ^ £
tri ^" "
Q- J3 £ !
O ^ 5
s!
s
^ *""* S i
1 ^* o I
m x*> * .
O ^ "
H^ a>
*-» 10
f- \ '
Si
2 *"* S
0 2 E
2 i. rt
In
£
^ 3 ?
K-l ^> *
I ^ «
2 S S
₯
5
2 -7 ;
tt ^ ?
0 S. o
S
CD
rf ^J «
i ^V rt
o 3- s
3
S
^ 'j "
'/s5 °
H 7»O S
Q s
R
S
S
O j
«j R
~
s
Q- ~
UJ ^ 2
2 O^ S
o * -
£
£ ii
2 $
e o
~ i
LJ =
0. e
< TE""
^ " *-
o » s
* o
*> <
H UJ
- 1
j
| !
j
i
i
\
\
|
, l
t
i
c
f
-x. -J to
0
r a
<
b
3
o
s
1 J
SI
s
IT
s
to
t _
£ q
i S ij.
e , o>
3 E O
E » £
" = S
c "~ K «"
- < » <
E t- « £1
S < » |_
(0- <
jrF ^
[ui- £
135
-------
CO
5
Q;
o
CD
Z
<
H
^
Q
r- T
Q. <
Z g
<
(V1
Q
UJ i
fVJ 1-
1 >
£t
§ Q
0 Z
u <
1 ("1
1 ^
1 (-J
Ul _J
Q
O H
2 z
o
CD Q.
Z
|_
Q;
V- W
t CL
_ ?~
I 1
< <
2 B
0 ^
Q
s
^
UJ
o:
1
i
o
^
o
z
Q. \^/
n i
(?)
5
CC
o
u.
ro
CO
z.
o
<\J
CO
o
^~
o
o
*"~
Q \
m "
2
o
> LJ J
Q >-C3
CO o1
0
*:
O.
5
h-
_|
». §
° <
Sf tr
O Q *
Ll- ^
l-|_
>o 2^
° Q;^
l~
o
LJ
SUJ
X
Q-
£5 Q
o ^<
cc 2
^
UJ
t-
^
Q
*
O
w
ro
£
ID
t
10
£
S
10
~
=
e
«
»
>"
»
m
»
n
-
u
_J
1
O-
1
Q
_J
1
Z
O
Q.
1
1
j
u
_l
f-
Q.
1
Q
_l
1-
Z
O
CL
i
<
«c
Q
Z
111
|g
₯
Si
s
kl
I
₯
t
£
S
w
5
~
^
S
»
ID
>-
U)
m
»
«i
OJ
-
_
O
(C (-
« §
u. i- e «
. UJ
u. m
- K
< § S
10 I ^
_r co n
°. _i C "-
If) CB
u. z
. UJ
X
< ^
~ m 1 *
< § S 5
z $ * S.
01 * 5 E
u. .= t-
136
-------
CO
5
o:
O
o
2
Q
o
o
^
1-
^
0 ru
h- '
^_) -J
Q. ^
2 5
Q;
Q
UJ £
CM h-
^ ^
"^ r^t
O~z_
1
1 Q
-J ^
UJ O
Q -1
O
^ r-
^
C3 o
2 Q.
\
ID
O <
o:
> uj
H Q.
' |
*5
^ ^
Q
Q
3>
^
UJ
en
1-
co
CO
^
o
2
"^ (w)
^
Q. "S
UJ
£
cc.
u.
"S
to
O
c
'-J
1 \
C/5 O
O ^
Q_ '"J
1 \
o: ^
o 5.
z
i ^
o §
2 ^~
' ^
o S£
2 B.
2 --J
' \
2 5-
2 3
S O
S.
< "3
i ^
o ^3-
_i
_j ^5
O
-^
UJ Q
O -1
UJ
Q.
X
2
<
Q
D
S
S
t
10
S
S
E
~
=
S
ai
>
i-
o>
10
>
T"
II
_J
H-
Q.
CM
1
0 -
_l
H
Z
O
[L
|
!
*
^-
-s.^
i
I
1
^n
1
*
1
4
i
1
i
I
u -«,^
j
^^
i
1
1
j
-~-,
I
1 ^
II
_l
1-
0.
EM
1
D
_J
Z
O
3.
[
1
I
!
S
(V
S
c
I
t
IS
£
H
E
!S
z
S
< <7>
O
t-
< £
(0
< f
Q m
UJ -
ir.
O
if
o>
0
in
X
"I
' CM
<
*
<
CM
<
^
O
U.
137
-------
tn
ij_
ry
u_
._
^
_
Q
O
O
g
O 2
H
/y
-1 UJ
UJ <
O uj
O o:
u? - Ld
t ^
_j (_
§ S
0 <
Q
*"C
^
4^T
UJ
UL
UJ
ff\
< (5)
Q.
LJ -S
\5/
i
£
*~ ^ s5
o o< -
S H
x o it
z
) SnjS^
^ _I>^O
a
UJ
a.
^
H
g
cc
en
K-
r-
1?
"
S
!L
^
^
o
1C
!c
0
s
U)
s
tr>
(D
S
z
w
Si
in
£
S
ui
5
ro
K
K
O
?
«r
s
v
tfi
J
;
»
;
0
S
S
S
s
s
s
Wl
?!
s
J>
m
CD
ru
S
S
~
1C
s
H
o
2
ID
t
*
£
S
5
~
=
9
O)
CO
K
10
«
»
fO
CO
-
"
2
<
a
<
»_
<
o
5
<
0
I
I
\
1
!
i
1
u
2
4
a
-------
tr\
u..
,,
o
2
^"
\-
Q
h-
I i
Q. >
2 ^
^
UJ §
CVI O
1 CD
2
^
O UJ
, «
1 1
_l £
UJ 1
Q ?;
0 Q
5
CD
2 f>
I
ID w
o °-
* ' >-
n" !_
<-£ t
S- <
Li H
^
_J Q
^
Z)
O
2
UJ
l_
1 1
CO
o
2
X,
ft Vi^/
III t
/CT\
VgJ/
S
U-
5
ttT ^_j
£ o'^
Zi 2°
o' ^
O
O
Z
<
to
z.
0
o
OJ
z.
o
o
z_
0
o
, .
_)
Q^
CD ^
^,
QJ^
>LiJ '
O^rn
CO^^
-o
Q
Q. rr
)5
ti |J-
^^
UJ
0.
V
H
<
^Z
Q
§
-
3
h
?
t
!£
in
2
12
~
z
o
o>
«
N
Kl
0
^
10
u
1
>-
IE
4
Q
Z
=
O
m
5
4
UJ
(E
|_
in
z
»
o
o
ro
t
*
n
2
5
~
Z
0
O)
CO
r-
UJ
10
V
n
<*
-
CVJ
^_
a:
o
z.
Z)
o
CD
^
2
1 1 1
cc
h-
'Z.
^
o
Q
<
fO
UJ
CL
h-
<
1
Q
o
Cft
o
c
"^
'
CO fo
Q 2
Q_
,
O rr,
_. ^
O .
"Z. ^~
! 1
O 5
Z --
z --~
1 x^
CVI CO
0 "5.
z
, H
TO C5
I 5
z
.
^- 1
o ^
°~ 5
o ^_
^_ ^
^ -J
_J "*
1^
UJ
Q.
2
^
o
«
s
s
Si
£
CD
r-
tc
IO
J
to
~
-
O
e»
E)
r~
(D
10
»
m
01
-
CVI
1
>
nr
<
a
z
3
o
CD
S
<
UJ
ec.
i-
en
z
$
O
a
<
K)
<
t-
<
a
z
UJ
0
01
E>
t
S
*?
2
12
~
r
0
CB
CD
1^
^^
<
S
a:
u.
139
-------
(T
O
Q
O
O
Q.
UJ
CM
O
I
UJ
Q
O
o
oc
o
LU
DC.
\-
O)
<
O
<
o
_J
o
o
I
IJ
o
o
o
Q
LU
UJ
z
LU
|_ £E
UJ (f)
S 00
< CL
m
m
30.
lit
UJ '
$
CO
LU
CH
D
in
UJ
z
o
o
_jO
OH-
UJQ*
CC
J3 Jt-
z '-1
£ <
UJ
Q.
^
1-
Q
g
2
GO
t
40
E
S
2
~
=
o
o)
CO
1^
J.
<1
»
Kl
tu
-
>-
O
o
-i
o
H
<
Z
_J
O
-1
<
O
o
_l
*s>
M "1
~*^,
1
_^
^-
--1
^ ^
^*-
^
"^.
_«
X"
**,
-^
-H
>*
>-
u
o
_l
o
^-
<
S
_i
0
_i
<
o
o
-1
S
Si
S
a
i
00
S
£
E
»
2
~
=
S
01
o>
t-
10
10
»
K)
"5
J
0
c
.2
£
5
I91
w
c
*
o
o
3
o>
0
o>
0
i>
E
10
o
o
"3 f
s i
§ E
^. « i-
°3 ^
0 -o
o" * co
Si 5
" 1 i
X .B X
O 0 >
« S
S
* z
OD L. i -
tc.
<
Q.
o.
11
II
T °-
D J£
* s
o i
140
-------
CO
S
cc
O
U.
Q
o
o
Q
I-
3
CL
LJ
(M
Z)
O
LU
O
O
O
QL
UJ
cr
CO
<
o
z
X
2
LJ
O
I
o
UJ
CC
Q
O
m
a
<
o
o
.
1-
^
h-
g
O)
2
£
t
£
£
s
E
2!
=
« X
o. O
cc CC
i-
u
«i Z
»
ft O
IM UJ
- 03
X
u
K
h-
o
_l
1
1
'
1
i
1
I
1
1
1
! 1
|
1
1
i
1
I
X
u
CC
2
o
UI
ffi
R
K
s
er
£
t
w
r
J
D
ff
=
X S
O »
CC o
r-
10
n
H- «
O >
_J ~
Q. -
3
B.
0
.Sot
* tt^
"-.is
. x c g
O a
£ P
- 0 &
S -5 *"
- O * to
" -S »
.. *>
t E
10 -e i
- -, E c
X c o
^* u
" * 5 i
5 Z |_
ceoo
O ui _i
u. m a
o.
o
o
I
u
cr
O
a.
x
u
(£
141
r us GOVERNMENT PRINTING OFFICE 1985- 559-111/20668
------- |