&EPA
United States
Environmental Protection
Agency
Environmental Research
Laboratory
Athens GA 30605
EPA-600/9-80-01 5
April 1980
Research and Development
Users Manual for
Hydrological
Simulation Program
Fortran (HSPF)
-------
RESEARCH REPORTING SERIES
Research reports of the Office of Research and Development, U.S. Environmental
Protection Agency, have been grouped into nine series. These nine broad cate-
gories were established to facilitate further development and application of en-
vironmental technology. Elimination of traditional grouping was consciously
planned to foster technology transfer and a maximum interface in related fields.
The nine series are:
1. Environmental Health Effects Research
2. Environmental Protection Technology
3. Ecological Research
4. Environmental Monitoring
5. Socioeconomic Environmental Studies
6. Scientific and Technical Assessment Reports (STAR)
7. Interagency Energy-Environment Research and Development
8. "Special" Reports
9. Miscellaneous Reports
This document is available to the public through the National Technical Informa-
tion Service, Springfield, Virginia 22161.
-------
EPA-600/9-80-015
April 1980
USERS MANUAL FOR HYDROLOGICAL
SIMULATION PROGRAM - FORTRAN (HSPF)
by
Robert C. Johanson
John C. Imhoff
Harley H. Davis, Jr.
Hydrocomp Incorporated
Mountain View, California 940UO
Grant No. R804971-01
Project Officer
Thomas 0. Barnwell
Technology Development and
Applications Branch
Environmental Research Laboratory
Athens, Georgia 30605
ENVIRONMENTAL RESEARCH LABORATORY
OFFICE OF RESEARCH AND DEVELOPMENT
U.S. ENVIRONMENTAL PROTECTION AGENCY
ATHENS, GEORGIA 30605
-------
DISCLAIMER
This report has been reviewed by the Environmental Research Laboratory, U.S.
Environmental Protection Agency, Athens, Georgia and approved for
publication. Approval does not signify that the contents necessarily
reflect the views and policies of the U.S. Environmental Protection Agency,
nor does mention of trade names or commercial products constitute
endorsement or recommendation for use.
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 analytical tools based on greater know-
ledge 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 Technology Development and
Applications Branch develops management or engineering tools to help pollu-
tion control officials achieve water quality goals through watershed manage-
ment.
The development and application of mathematical models to simulate the
movement of pollutants through a watershed and thus to anticipate environ-
mental problems has been the subject of intensive EPA research for several
years. The most recent advance in this modeling approach is the Hydrological
Simulation Program - FORTRAN (HSPF), which uses digital computers to simulate
hydrology and water quality in natural and man-made water systems. HSPF is
designed for easy application to most watersheds using existing meteorologic
ahd hydrologic data. Although data requirements are extensive and running
costs are significant, HSPF is thought to be the most accurate and appropri-
ate management tool presently available for the continuous simulation of
hydrology and water quality in watersheds.
David W. Duttweiler
Director
Environmental Research Laboratory
Athens, Georgia
iii
-------
ABSTRACT
The Hydrological Simulation Program - FORTRAN (HSPF) is a set of compu-
ter codes that can simulate the hydrologic, and associated water quality,
processes on pervious and impervious land surfaces and in streams and well-
mixed impoundments. The manual discusses the modular structure of the system,
the principles of structured programming technology, and the use of these
principles in the construction of the HSPF software. In addition to a pictor-
ial representation of how each of the 500 subprograms fits into the system,
the manual presents a detailed discussion of the algorithms used to simulate
various water quantity and quality processes. Data useful to those who need
to install, maintain, or alter the system or who wish to examine its structure
in greater detail are also presented.
The report was submitted in fulfillment of Grant No. R804971-01 by
Hydrocomp, Inc., under the sponsorship of the U.S. Environmental Protection
Agency. The report covers the period November 1, 1976, to November 30, 1978,
and work was completed as of January 16, 1980.
iv
-------
CONTENTS
Foreword ill
Abstract iv
Part
A Introduction 1
B General Principles 8
C Standards and Conventions 23
D Visual Table of Contents 45
E Functional Description 130
F Format for the Users Control Input 320
Appendices
I Glossary of Terms 635
II Sample Runs 642
III Program NEWTSS 656
IV Guide to the Programmers Supplement 674
-------
Introduction
PART A INTRODUCTION
CONTENTS
1.0 Purpose and Scope of the HSPF Software . 2
2.0 Requirements for HSPF 4
3.0 PUrpose and Organization of this Document 5
U.O Definition of Terms * 6
5.0 Notice of User Responsibility 6
6.0 Acknowledgments 6
-------
Introduction
1.0 PURPOSE AND SCOPE OF THE HSPF SOFTWARE
The use of models which simulate continuously the quantity/quality processes
occurring in the hydrological cycle is increasing rapidly. Recently there
has been a proliferation in the variety of models and in the range of
processes they simulate. This has been a mixed blessing to a user. To get
the benefits of simulation, he has to select a model from a bewildering
array and then spend much effort amassing and manipulating the huge
quantities of data which the model requires. If he wishes to couple two or
more subprocess models to simulate a complete process, he often encounters
further difficulties. The underlying assumptions and/or structures of the
subprocess models may make them somewhat incompatible. More frequently, the
data structures are so different that coupling requires extensive data
conversion work.
One reason for these problems is that the boom in modeling work has not
included enough work on the development of good model structures. That is,
very few software packages for water resource modeling are built on a
systematic framework in which a variety of process modules can fit.
With HSPF we have attempted to overcome these problems as far as possible.
HSPF consists of a set of modules arranged in a hierarchical structure,
which permit the continuous simulation of a comprehensive range of
hydrologic and water quality processes. Our experience with sophisticated
models indicates that much of the human effort is associated with data
management. This fact, often overlooked by model builders, means that a
successful comprehensive model must include a sound data management
component. Otherwise, the user may become so entangled in data manipulation
that his progress on the simulation work itself is drastically retarded.
Consequently, the HSPF software is planned around a time series management
system operating on direct access principles. The simulation modules draw
input from a Time Series Store and are capable of writing output to it.
Because these transfers require very few instructions from the user, the
problems referred to above are minimized.
The system is designed so that the various simulation and utility modules
can be invoked conveniently, either individually or in tandem. A top down
approach emphazing structured design has been followed. First, the overall
framework and the Time Series Management System were designed. Then, work
progressed down the structure from the highest, most general level to the
lowest, most detailed one. Every level was planned before the code was
written. Uniform data structures., logic figures, and programming
conventions were used throughout. Modules were separated according to
function so that, as much as possible, they contained only those activities
which are unique to them. Structured design has made the system relatively
easy to extend, so that users can add their own modules with relatively
little disruption of the existing code.
-------
Introduction
Now, a note on the initial contents of the system. Presently, it includes
modules which can handle almost all the functions which are available in the
following existing models:
(1) HSP (LIBRARY, UTILITY, LANDS, CHANNEL, QUALITY)
(2) ARM
(3) NFS
The HSPF software is not merely a translation of the above models, but a new
system with a framework designed to accomodate a variety of simulation
modules; the modules described above are the initial contents. Many
extensions have been made to the above models in the course of restructuring
them into the HSPF system.
It is hoped that HSPF will become a valuable tool for water resource
planners. Because it is more comprehensive than most existing systems, it
should permit more effective planning. More specifically, the package can
benefit the user in the following ways:
(1) The time-series-oriented direct access data system and its
associated modules can serve as a convenient means of inputting,
organizing, and updating the large files needed for continuous
simulation.
(2) The unified user oriented structure of the model makes it
relatively simple to operate. The user can select those modules
and options that he wishes to execute in one run, and the system
will ensure that the correct sets of code are invoked and that
internal and external transfers of data are handled. This is
achieved with a minimum of manual intervention. Input of control
information is simplified because a consistent system is used for
this data for all the modules.
(3) Because the system has been carefully planned using modern
top-down programming techniques, it is relatively easy to modify
and extend. The use of uniform programming standards and
conventions has assisted in this respect.
(1) Since the code is written almost entirely in ANSI standard
Fortran, implementation on a wide variety of machines is possible.
-------
Introduction
2.0 REQUIREMENTS FOR HSPF
In awarding the grant for development of HSPF, the EPA set the following
requirements:
(1) It must manage and perform deterministic simulation of a variety
of aquatic processes which occur on and under land surfaces and in
channels and reservoirs.
(2) It must readily accommodate alternate or additional simulation
modules.
(3) It must permit easy operation of several modules in series, and
thus be capable of feeding output from any operation to subsequent
operations.
(4) It must be in ANSI Fortran with minor specified extensions.
With the concurrence of the EPA, we expanded on these requirements:
(1) It must have a totally new design. Existing modules should not
merely be translated, but should be fitted into a new framework.
(2) It must be designed from the top down, using some of the new
improved programming techniques., such as Structured Design and
Structured Programming.
(3) Duplication of blocks of code which perform similar or identical
functions should be avoided.
(4) The user's control input must have a logically consistent
structure throughout the package.
(5) Uniform standards and practices must be followed throughout the
design, development and documentation of the system.
(6) It must have a conveniently operated disk-based Time Series Store
built on the principle of direct access.
(7) The design must be geared to implementation on larger models of
the current generation of "minicomputers." It must be compatible
with Operating Systems which share machine core space using either
the virtual memory approach or a conventional overlay technique.
-------
Introduction
3.0 PURPOSE AND ORGANIZATION OF THIS DOCUMENT
This report contains all the documentation of the HSPF system. It is
designed to:
(1) introduce new users to the principles and concepts on which the
system is founded
(2) describe the technical foundations of the algorithms in the
various application (simulation) modules
(3) describe the input which the user supplies to run the system
To meet these needs and, at the same time, to produce a document which is
reasonably easy to use, we have divided this report into several distinct
parts, each with its own organization and table of contents.
Part A (this one) contains introductory material.
Part B outlines the general principles on which the HSPF system is based.
This includes a discussion of the "world view" which our simulation modules
embody. A firm grasp of this material is necessary before the detailed
material can be properly understood.
Part C explains the standards and conventions which we employed to develop
the HSPF software. It discusses the advantages of Structured Programming
and describes how we implemented this technique in standard Fortran.
Part D is a visual table of contents. It displays pictorially the entire
hierarchy of subprograms, starting with the MAIN program and proceeding down
the "tree" to the most detailed level of the system. It gives the name,
number and a brief description of each of the 500 subprograms.
Part E documents the function of each part of the software. The
organization of this part follows the layout of the software itself. The
relationship between, and the functions of, the various modules are
described, starting at the highest most general level and proceeding down to
the lowest most detailed level following the numbering system in the visual
table of contents. The algorithms used to simulate the quantity and quality
processes which occur in the real world are described in this part.
Part F describes the User's Control Input; that is, the information which
the user must provide in order to run HSPF.
Material which might obscure the structure of this document if it were
included in the body of the report appears in Appendices. These include a
glossary of terms and sample runs.
-------
Introduction
M.O DEFINITION OF TERMS
In this document, terms which have a special meaning in HSPF, are enclosed
in quotes the first time they occur. Usually an explanation follows
immediately. A glossary of terms will be found in Appendix I.
5.0 NOTICE OF USER RESPONSIBILITY
This product has been carefully developed. Although the work included
testing of the software, the ultimate responsibility for its use and for
ensuring correctness of the results obtained, rests with the user.
The EPA and the developers of this software make no warranty of-any kind
with regard to this software and associated documentation, including, but
not limited to, the implied warranties of merchantability and fitness for a
particular purpose. They shall not be liable for errors or for incidental
or consequential damages in connection with the furnishing, performance or
use of this material.
While we intend to correct any errors which users report, we are not obliged
to do so. We reserve the right to make a reasonable charge for work which
is performed for a specific user at his request.
6.0 ACKNOWLEDGMENTS
This work was sponsored by the Environmental Reasearch Laboratory in Athens,
Georgia. David Duttweiler is the laboratory director and Robert Swank the
head of the Technology Development and Applications Branch, which supervised
the project. The Project Officers were Jim Falco and, later, Tom Barnwell.
Many people at Hydrocomp Inc. were involved in the project . The
role of each of the of the major contributors is summarized below. Persons
are listed approximately in order of time spent on the project, but this
does not necessarily reflect the relative value of their contributions.
Robert Johanson was Project Manager. He was responsible for project
coordination, development of the standards and practices and much of the
basic design work. He supervised closely the development of the application
modules and wrote the SNOW and PWATER sections of the PERLND module, and the
HYDR section of the RCHRES module. He was also responsible for the Run
Interpreter.
John Imhoff worked on the RCHRES module. He analyzed the HSP QUALITY code,
performed the detailed design of the new module and wrote the code and
documentation for it. He also coordinated the production of the functional
descriptions (Part E) of all the application modules.
Harley Davis designed and coded most sections of the PERLND module and all
of the IMPLND module. He also wrote the functional descriptions of those
modules.
-------
Introduction
Janice Walters worked on the Time Series Management System (TSMS),
particularly the TSGET and TSSMGR groups of subprograms.
Gerard Lum had the tedious, but critical, task of translating most of the
pseudo code into ANSI standard Fortran.
Delbert Franz participated in the overall design of the system. He also
supervised the work on the TSMS, and produced most of the code for the TSPUT
group of subprograms.
Tamas Eger did the detailed design and coding of most of the TSSMGR
subroutine group. He also wrote the stand-alone program NEWTSS.
Carolyn Karnos had the onerous task of assembling the Fortran subprograms
into the program file and contending with the "bugs" which arose in this
process.
Jean-Jaques Heler participated in the design of the system, especially in
shaping the concepts discussed in Part B. He also finished and tested the
subprogram groups TSSMGR and NEWTSS after Tamas Eger had to return to
Hungary.
Norman Crawford participated in the initial design of the system. As a
result, many of his ideas for the HSP-II system found their way into HSPF.
Eleanor Bassler worked on the design and organization of the Users Control
Input.
Alan Schiffer translated many of the RCHRES subprograms from pseudo to
Fortran code.
Diana Allred edited much of the text in this document, especially Part E.
Donna Mitchell also contributed to the editing effort.
Tony Donigian reviewed the text in Part E for technical accuracy.
Nancy Sharpe spent many hours at the copying machine, producing the working
documentation from which this document evolved.
Jack Kittle assisted in assembling the code into the program file, and also
set up the system for arranging data in the numerous versions of the COMMON
block.
Kevin Gartner translated some pseudo code into Fortran.
-------
General Principles
PART B GENERAL PRINCIPLES
CONTENTS
1.0 View of the Real World
1.1 General Concepts
1.2 Nodes, Zones, and Elements
1.3 Processing Units and Networks
2.0 Software Structure- ....................... 1
2.1 Concept of an "Operation" ................. 14
2.2 Time Series Storage . ................... 16
2.3 Times Series Management for an Operation ......... 17
2.U HSPF Software Hierarchy .................. 17
3.0 Structure of a Job ..................... . . 20
3.1 Elements of a Job ................... . . 20
3.2 Groups of Operations ..... ; ............. 20
FIGURES
Number Page
1-1 Nodes, zones and elements 10
1-2 Directed and non-directed graphs 12
1-3 Single- and multi-element processing units 13
2-1 Logical structure of the internal scratch pad 15
2-2 Activities involved in an operation 13
2-3 Overview of HSPF software 19
3-1 Schematic of data flow and storage in a single run 21
3-2 Extract from typical User's Control Input, showing how
grouping of operations is specified 22
-------
General Principles
1.0 VIEW OF THE REAL WORLD
1.1 General Concepts
To design a comprehensive simulation system, one must have a consistent
means of representing the prototype; in our case, the real world. We view
it as a set of constituents which move through a fixed environment and
interact with each other. Water is one constituent; others are sediment,
chemicals, etc. The motions and interactions are called processes.
1.2 Nodes, Zones, and Elements
The prototype is a continuum of constituents and processes. Simulation of
such a system on a digital computer requires representation in a discrete
fashion. In general, we do this by subdividing the prototype into
"elements" which consist of "nodes" and "zones."
A node corresponds to a point in space. Therefore, a particular value of a
spatially variable function can be associated with it, for example, channel
flow rate and/or flow cross sectional area. A zone corresponds to a finite
portion of the real world. It is usually associated with the integral of a
spatially variable quantity, for example, storage in a channel reach. The
zone the smallest unit into which we subdivide the world. The relationship
between zonal and nodal values is similar to that between the definite
integral of a function and its values at the limits of integration.
An element is a collection of nodes and/or zones. Figure 1-1 illustrates
these concepts. We simulate the response of the land phase of the
hydrological cycle using elements called "segments." A segment is a portion
of the land assumed to have areally uniform properties. A segment of land
with a pervious surface is called a "Pervious Land-segment" (PLS).
Constituents in a PLS are represented as resident in a set of zones (Fig.
1-1a). A PLS has no nodes. As a further example, consider our formulation
of channel routing. We model a channel reach as a one dimensional element
consisting of a single zone situated between two nodes (Fig. 1-1b). We
simulate the flow rate and depth at the nodes; the zone is associated with
storage.
The conventions of the finite element technique also fall within the scope
of these concepts. Figure 1-1c shows a two dimensional finite element used
in the simulation of an estuary. Three nodes define the boundaries of the
triangular element. A fourth node, situated inside, may be viewed as
subdividing the element into three zones. This last type of element is not
presently used in any HSPF module, but is included in this discussion to
show the generality provided by HSPF. The system can accommodate a wide
variety of simulation modules.
-------
General Principles
There are no fixed rules governing the grouping of zones and nodes to form
elements. The model builder must decide what grouping is reasonable and
meaningful, based on his view of the real world processes being simulated.
In the foregoing material we presented some elements used in HSP and other
systems. In general, it is convenient to define elements so that a large
portion of the real world can be represented by a collection of conceptually
identical elements. In this way, a single parameter structure can be
defined which applies to every element in the group. Thus, each element is
a variation on the basic theme. It is then meaningful to speak of an
"element type." For example, elements of type "PLS" all embody the same
arrangement of nodes and are represented by sets of parameters with
identical structure. Variations between segments are represented only by
variations in the values of parameters. The same applies to any other
element, such as a Reach, layered lake or a triangular finite element.
As illustrated in the above discussion, nodes are often used to define the
boundaries of zones and elements. A zone, characterized by storage,
receives inflows and disperses outflows; these are called "fluxes." Note
that if the nodal values of a field variable are known, it is often possible
to compute the zonal values (storages). The reverse process does not work.
1.3 Processing Units and Networks
To simulate a prototype we must handle the processes occurring within the
elements and the transfer of information and constituents between them. The
simulation of large prototypes is made convenient by designing a single
"application module" for a given type of element or element group, and
applying it repetitively to all similar members in the system. For example,
we may use the RCHRES module to simulate all the reaches in a watershed
using storage routing. This approach is most efficient computationally if
one element or group of elements, called a "processing unit" (PU), is
simulated for an extended period of time before switching to the next one.
To permit this, we must be able to define a processing sequence such that
all information required by any PU comes from sources external to the system
or from PU's already simulated. This can only happen if the PU's and their
connecting fluxes form one or more networks which are "directed graphs." In
a directed graph there are no bidirectional paths and no cycles. Figure 1-2
shows some directed and nondirected graphs.
The requirement that PU's form directed graphs provides the rule for
grouping elements into PU's. Any elements interacting with each other via
loops or bidirectional fluxes must be grouped into a single PU because none
of them can be simulated apart from the others.
Thus, we can have both single element and multielement PU's. A PLS is an
example of the former and a channel network simulated using the full
equations of flow exemplifies the latter (Fig. 1-3). A multielement PU is
also known as a "feedback region." The collection of PU's which are
simulated in a given run is called a "network."
11
-------
Cb)
CO
Cd)
Directed
Graphs,
unit", w'i"H*v-f
-------
PLS
PUs
(simu
-full ea
mo-Hfon)
. 1-3
-------
General Principles
The processes which occur within a PU are represented mathematically in an
"application model." The corresponding computer code is called an
"application module" or "simulation module."
2.0 SOFTWARE STRUCTURE
2.1 Concept of an "Operation"
A great variety of activities are performed by HSPF; for example, input a
time series to the Time Series Store, find the cross correlation coefficient
for two time series, or simulate the processes in a land segment. They all
incorporate two or more of the following functions: get a set of time
series, operate on the set of input time series to produce other time
series, and output the resulting time series. This applies both to
application modules (already discussed) and to "utility modules," which
perform operations ancillary or incidental to simulation. Thus, a
simulation run may be viewed as a set of "operations" performed in sequence.
All operations have the following structure:
SUPERVISE
OPERATIONS
(subroutine OSUPER)
i i i
GET TIME OPERATE PUT TIME
SERIES (utility SERIES
(subroutine group or (subroutine group
TSGET) application TSPUT)
module)
The OPERATE function is the central activity in the operation. This work is
done by an "operating module" (OM) and its subordinate subprograms. They
operate for a specified time on a given set of input time series and produce
a specified set of output time series, under control of the "operations
supervisor" (OSUPER). All of the pieces of time series involved in this
internal operation have the same interval and duration. They are therefore
viewed as written on an "internal scratch pad" (INPAD), resident in the core
of the machine (Fig. 2-1). The operating module receives the scratch pad
with some rows filled with input and, after its work is done, returns
control to the supervisor with another set of rows filled with output. The
operating module may overwrite an input row with its own output. The
computing module being executed, together with the options being invoked,
will determine the number of rows required in the INPAD. For example,
simulation of the hydraulic behavior of a stream requires relatively few
time series (eg. inflow, depth and outflow) but the inclusion of water
14
-------
General Principles
quality simulation adds many more time series to the list. Now, the total
quantity of machine space available for storage of time series is also fixed
(specified in a COMMON block) by the options in effect; this is the size
("area") of the INPAD. Since both the size (N*M) and number of rows (M) in
the INPAD are known, the "width" (no. of intervals.N) can be found. The
corresponding physical time is called the "internal scratch pad span
(INSPAN)."
Y\O.
I
4
r\o*>.
to
- - N
3
A
M
'. -Fhe^e '{& one -H we
per row.
15
-------
General Principles
The "get time series" function prepares the input time series. This work is
done by a subroutine group called TSGET. It obtains the correct piece of a
time series from the appropriate file, aggregates or disaggregates it to the
correct time interval, multiplies the values by a user specified constant
(if required), and places the data in the required row of the internal
scratch pad. Subroutine group TSPUT performs the reverse set of operations.
TSGET and TSPUT are sometimes bypassed if a required time series is already
in the INPAD when the operation is started, or if the output is being passed
to the next operation via the internal scratch pad.
Modules TSGET and TSPUT are part of the "time series management system"
(TSMS).
2.2 Time Series Storage
The time series used and produced by an operation can reside in three types
of storage.
(1) The Time Series Store (TSS)
This is the principal library for medium-long term storage of time
series. As far as the machine operating system is concerned, it
consists of a single large direct access file on a disc. HSPF
subdivides this space into many datasets containing time series. Each
is logically self-*contained but may be physically scattered through the
store. A directory keeps track of datasets and their attributes.
Before time series are written to the Time Series Store, the store must
be initialized and its directory created. This is done by executing
the separate program NEWTSS, which is documented in Appendix III.
(2) Sequential Files
These are files with a constant logical record length which are
resident on or routed to punched cards, a magnetic tape, a sequential
disk file or a line printer. Time series received from agencies such
as the National Weather Service are typically on sequential files
(cards or tape).
(3) Internal Scratch Pad (INPAD)
If two or more operations performed in sequence use the same internal
time step, time series may be passed between them via the INPAD.
Successive operations may simply pick up the data written by the
previous ones, without any external (disc) transfer taking place. This
is typically done when time series representing the flow.of water (and
constituents) are routed from, one stream reach to the one next
downstream.
16
-------
General Principles
2.3 Time Series Management For An Operation
Any operation involves a subset of the activities shown in Fig. 2-2. The
operating module expects a certain set of time series in the INPAD. The
operations supervisor, acting under user control, ensures that the
appropriate input time series are loaded from whichever source has been
selected, and informs the computing module of the rows in the INPAD where it
will find its input. Similar arrangements hold for output of time series.
2.M HSPF Software Hierarchy
The hierarchy of functions in HSPF is shown in Fig. 2-3. Some explanatory
notes follow.
The "Run Interpreter" is the group of subprograms which reads and interprets
the "Users Control Input." It sets up internal information instructing the
system regarding the sequence of operations to be performed. It stores the
initial conditions and the parameters for each operation in the appropriate
file on disc and creates an instruction file which will ensure that time
series are correctly passed between operations, where necessary.
The "TSS management" modules are those used to create, modify, or remove
data sets in the time series store.
The "Operations Supervisor" is a subroutine which acts on information
provided by the Run Interpreter, invoking the appropriate "application" or
"utility" modules. It provides them with the correct values for parameters
and state variables by reading the files created by the Run Interpreter.
Operating modules are either "application modules" or "utility modules."
They perform the operations which make up a run. Each time one of those
modules is called, an operation is performed for a period corresponding to
the span of the internal scratch pad (INSPAN). The Operations Supervisor
ensures that the correct module is invoked.
"Service subprograms" perform tasks such as reading from and writing to time
series storage areas, adding T minutes to a given date and time, to get a
new date and time, etc.
The "Time Series Management System" (TSMS) consists of all the modules which
are only concerned with manipulation of time series or the files used to
store time series. It includes the TSS management functions, and TSGET and
TSPUT.
17
-------
Supervisor
SUB-ROUTlME. CAUL
TIME SERIES TRAMSFER FATH
ACTIVITIES IMVOLVEP irs4 AN opERA-no^4
18
-------
Run
MAIN
Service
r
-r-3-s
module
I
TSPUT
I application and tfHl't-hy modole-s.
Time Series 'Store.
CVerview of H5FT=
-------
General Principles
3.0 STRUCTURE OF A JOB
3.1 Elements of a Job
A "JOB" is the work performed by HSPF in response to a complete set of Users
Control Input. It consists of one or more "RUNs" and/or "Time Series Store
Management" activities. A RUN is a set of operations which can be performed
serially, and which all cover the same period of time (span). The operations
are performed in a sequence specified in the Users Control Input. To avoid
having to store large quantities of intermediate data on disc, operations
may be collected in a group in which they share a common INPAD (INGRP).
3.2 Groups Of Operations
In most runs, time series have to be passed between operations. As
described in Section 2.2, each operation can communicate with three
different time series storage areas: the TSS, the INPAD, and sequential
files. This is illustrated in Fig. 3-1.
Potentially, any time series required by or output by any operation can be
stored in the TSS or a sequential file. The user simply specifies the exact
origin or destination for the time series, and the HSPF system moves the
data between that device and the appropriate row of the INPAD. This system
can also be used to transfer data between operations. However, it does
require that all transferred data be written to the TSS or a sequential
file. This may be very cumbersome and/or inefficient and it is better to
transfer data via the INPAD, where possible.
To transfer data via the INPAD, operations must share the same pad. This
means that all time series placed in the pad have the same time interval and
span. This requirement provides a logical basis for grouping operations;
those sharing a common INPAD are called an INGRP (Fig. 3-1). The user
specifies the presence of groups in his "Users Control Input (UCI)." A
typical sequence of input is shown in Fig. 3-2.
The user also indicates (directly or indirectly) in his control input the
source and disposition of all time series required by or output by an
operation. If he indicates that a time series must be passed to another
operation then the system assumes that the transfer will be made via the
scratch pad. If they are not in the same INGRP there is an error. Without a
common INPAD, the data must go via the TSS. The structure of the Users
Control Input is documented in Part F.
20
-------
OPN1
OPN4
INPAD
OPN5
OPW11
OPN11
OPN14
INPAO
TSS Time series store
OPH Operation
IN6RP lnt?rnal scratch pod
qroup
INPAD Internal ^cratch pad
FI6 3-1 6£rtEKUn"l£ OF »OA FLOW/
FOR A SINGLE. RUN
21
-------
General Principles
The sequence of events in a run is as follows (refer to Fig. 3-1).
(a) Operation 1 is performed until its output rows in the INPAD are filled.
(b) Data are transferred from those rows to other time series storage
areas, as required. If any of these data are not required by other
operations in INGRP1, their INPAD rows are available for reuse by other
operations in INGRP1.
(c) Steps (a) and (b) are repeated for each operation in INGRP1.
(d) Steps (a), (b), and (c) are repeated, if necessary, until the run span
is complete.
(e) The INPAD is reconfigured and work on operations 5 through 11 proceeds
as in steps (a-d) above. The step repeats until all INGRP's have been
handled. The run is now complete.
Note that reconfiguration of a scratch pad implies that its contents will be
overwritten.
OPN SEQUENCE
INGRP INDELT = 00:30
COPY 1
PERLND 1
END INGRP
PERLND 2 INDELT = 00:30
PERLND 3 INDELT = 00:20
INGRP INDELT = 00:30
COPY 2
RCHRES 1
RCHRES 3
RCHRES 5
RCHRES 20
RCHRES 22
RCHRES 23
RCHRES 7
RCHRES 8
RCHRES 50
RCHRES 100
RCHRES 200
END INGRP
INGRP INDELT = 00:10
ANALYS 1
PLTGEN 1
END INGRP
END OPN SEQUENCE
Fig. 3-2 Extract from typical Users Control Input,
showing how grouping of operations is specified
22
-------
Standards & Conventions
PART C STANDARDS AND CONVENTIONS
CONTENTS
Page
1 . Structured Programming Technology ................. 24
1.1 Introduction ......................... 24
1.2 Structure Charts ....................... 25
1.3 HIPO Diagrams ........................ 27
1.4 Pseudo Code ......................... 27
t
2. Terminology ............................ 29
3. Fortran Language ......................... 30
3.1 Use of ANSI Fortran ..................... 30
3.2 Allowable Extensions to ANSI Fortran ............. 30
4. Fortran Conventions ........................ 31
4.1 The Need for Conventions ................... 31
4.2 Conventions used to Implement the Standard Structure Figures . 31
4.3 Other Conventions ...................... 34
5. Conventions Used* in Functional Description
6. Method of Documenting Data Structures
6. 1 Structure of Data in Core
6.2 Structure of Data on Disc Files
7. Method of Handling Diagnostic Messages
8. References
23
-------
Standards & Conventions
1.0 STRUCTURED PROGRAMMING TECHNOLOGY
1.1 Introduction
Developing, modifying, or even trying to understand, a large computer
program can be a very frustrating activity, as many people who have been
involved in those kinds of work can testify. Typical problems encountered
are:
(1) It is hard to see the relationship between various parts of the
program. The underlying design, if there was any, has become
obscured as the program evolved to its current state.
(2) The logic of the program is very contorted. The flow chart
resembles a bowl of spaghetti.
(3) "Bugs" are difficult to locate.
(4) A seemingly trivial alteration to the program sets off a chain
reaction. The resulting problems take a long time to solve.
Why do these things happen? Must long programs necessarily suffer from
these ills?
There is no quick answer which covers all aspects of the problem. However,
many of these troubles are directly attributable to program structure, or
the lack of it. Typically, when the program was designed, an unsuitable
structure was chosen or no thought was given to the matter. The reason for
this situation is not hard to find. Most engineers are not taught how to
design programs; they are only taught how to write program code. There is a
great difference between these two activities. It is like the difference
between writing a novel and supplying short answers to a series of
questions.
The importance of structure emerged as computer scientists investigated the
reasons for these problems. As this work went on, a set of new disciplines
evolved, which emphasized structure in the design, implementation, and
documentation of software. Collectively, they are referred to as
"Structured Programming Technology." The individual techniques are:
(1) Structured Design
(2) Top Down Programming
(3) Structured Programming
(4) Hierarchy plus Input, Process, Output (HIPO) Diagrams
(5) use of a Development Support Library
It has been found that the use of these methods, eit-her singly or in
combination, leads to computer programs and documentation which are easier
to understand than those produced using conventional methods. Also, the
software has a unified underlying structure and it is easier to maintain.
24
-------
Standards & Conventions
Because of its obvious suitability, we made extensive use of Structured
Programming Technology on this project. To familiarize us with the
techniques, we acquired the "Textbook and Workbook on Structured
Programming," and "HIPO - a Design Aid and Documentation Technique" (IBM
1974). These books contain an independent study program which gives the
reader a thorough grounding in the underlying principles and practical
application of the techniques. Because the Textbook, Workbook, and HIPO
Manual are comprehensive and readable, we used them as reference works and
kept as close as possible to the concepts and terminology which they
contain. This should become apparent on reading further through this
document. Through the remainder of this manual, those documents will simply
be referred to as the "Textbook," "Workbook," and "HIPO Manual,"
respectively.
1.2 Structure Charts
A structure chart (SO is a diagram which shows the functions performed by
the various modules in a program and their hierarchy. The use of SC's in
Structured Program Design, including the methods used to compile them, is
described in Section 6 of the Textbook. We suggest that the reader review
that material before proceeding.
The HIPO method of documentation involves the use of a "visual table of
contents" (see HIPO Manual), which is very similar to an SC. In our
documentation system, we combine these two functions in a single diagram,
which we still call an SC, e.g. SC 4.0. The visual table of contents for
HSPF is given in Part D of this document.
To compile structure charts, we use a combination of the methods given in
the Textbook and in the HIPO Manual, together with some ideas and
conventions of our own:
(1) Typically, the SC for a large program extends over several pages.
We use the term SC either to describe the entire chart or to refer
to that part of it shown on a single page.
(2) Each box in an SC can be viewed in two distinct ways:
(a) as representing a single subprogram (the one which calls all
the subprograms which perform subordinate functions)
(b) as representing the entire group of subprograms mentioned in
(a)
In our documentation system we use this dual meaning to an
advantage. Details are given later.
25
-------
Standards & Conventions
(3) Each box contains a description of the action which it and all its
subordinate boxes collectively represent. It follows that the
description given in any box must be implied in the descriptions
given in all boxes in the chain connecting it to the MAIN program
box, e.g. the description for box 4.1.1 in SC 4.1 is implied in
(or is part of) the descriptions for boxes 4.1, 4.0, and 1.0.
Because it describes the action performed by one or more
subprograms, a description starts with the imperative form of a
verb, e.g. compute, find, get (HIPO Manual, p. 63). In an SC, the
description appears in lower case.
(4) The numbering system is as described in the HIPO Manual. Numbers
are placed inside the boxes, near the lower right hand corner. We
have introduced one extension to the numbering system because the
function "perform an operation" refers to any of the operating
modules in the system, which are all conceptually similar and of
equal rank. Therefore, we assign them all the number 4.2 and
distinguish them only by a subscript. For example, 4.2(3) refers
to the Reach/Mixed Reservoir Application Module, and 4.2(3).1
refers to Section 1 of that module (Hydraulics).
(5) There is a name next to each box, outside the top left hand
corner. It identifies the Fortran subprogram represented by the
box. Note that this name applies to a single box, not an entire
subprogram group. It is written in upper case.
(6) An SC bears the number of the highest level box which it depicts,
e.g. SC 4.0 shows subprogram 4.0 and its subordinates.
(7) Where an SC continues across more than one page we place "flags"
at points where the chart has been broken. They are labeled with
the chart number on which the information is continued. For
example , see the stubs at the bottom of SC 1.0 in Part D.
The entire set of SCs for the HSPF system appear together in the Visual
Table of Contents (Part D). It consists of the full set of Structure Charts,
arranged in numerical order. Its purpose is to enable a user to:
(1) see how the following material (the functional description in Part
E) is arranged
(2) get an overview of the en.tire system
(3) locate a given subprogram in the system
(4) see a subprogram in its surroundings, so its relationship to the
rest of the system can be better understood
(5) design an overlay structure, if HSPF is being installed on a
computer which requires it.
26
-------
Standards & Conventions
1.3 HIPO Diagrams
A HIPO diagram depicts graphically the function of a subprogram and its
subordinate subprograms. A HIPO Diagram bears the same number as the SC box
which it illustrates. Thus, for example, HIPO diagram 4.0 describes the
functions performed by the entire subprogram group 4 from the perspective of
the head subroutine (OSUPER).
Although some HIPO diagrams were prepared when HSPF was designed, we did not
find them very useful, so they are not included in this document. As we
become more familiar with Structured Program Design, we might develop a
better appreciation of their value.
1.4 Pseudo Code
Pseudo code is an English-like noncompilable language which describes
program logic. It is discussed in Section 4 of the Textbook. Its chief
advantage, compared to Fortran code, is that it encourages the writing of
programs with good structure which are easily read and understood. It does
not contain statement numbers or GO TO statements. We developed our own
version, almost identical to that used in the Textbook. The entire HSPF
system was written in pseudo code (Johanson, et al. 1979) before being
translated into Fortran. We observed the following conventions:
(1) Only valid structure figures are used. These include the
SEQUENCE, IFTHENELSE, WHILEDO, DOUNTIL, and CASE figures.
(2) Successive levels of indentation are displaced from each other by
two columns.
(3) Fortran identifiers and subprogram names are in upper case, e.g.
RCHTAB.
(4) Other identifiers are in lower case with initial capital letter,
e.g. Time-parms. This applies both to names of groups of
variables (e.g. Date-time) and to names of any individual
variables which have not been reduced to a Fortran identifier,
e.g. Land-segment.
(5) Keywords are in upper case, e.g. DOUNTIL, BEGIN, IF, THEN.
(6) Other text is in lower case (no initial capitals).
(7) "CALL" means a subprogram is called.
(8) "INCLUDE" means a named segment of code will be physically
inserted.
(9) In a DOUNTIL it is assumed that the initial value of the index (if
used) is 1 and that it is incremented by 1, unless other values
are previously assigned and/or specified.
27
-------
Standards & Conventions
(10) Character string constants are enclosed in single quotes, e.g.
'TSSM1.
(11) Two-valued flags have the value ON and OFF or 0 and 1.
(12) Comments are enclosed in brackets.
(13) Where we refer to all the elements of a vector or array we often
use the PL/1-like notation involving *. For example,
(a) we use
OM(») = FACT*OS(»)
instead of
DOUNTIL N = NOFLO
OM(N) = FACT*OS(M)
ENDDO
(b) we use
BEGIN SETTPT (RCHTAB(»), )
instead of
BEGIN SETTPKRCHTAB, )
to denote the entire vector RCHTAB
(14) Argument lists are enclosed in parentheses, with input arguments
first, input/output arguments next and output arguments last.
These groups are separated by 'bb1 or placed on separate lines as
is done in our Fortran.
(15) We often use a cross reference to save time and trouble, e.g.
.INCLUDE COMMON/SCRTCH/ version HYDR2
CALL AUXIL (argument list as in AUXIL)
(16) We use the following symbols for logical operators:
=, not=, >, <, >=, <=
(17) We follow our Fortran coding conventions, where applicable
(Section 4).
(18) Every program unit has a comment immediately following the BEGIN
statement, outlining its purpose (as in Fortran).
(19) We kept program units short, as the Textbook recommends. We made
free use of "CALL" to achieve this.
28
-------
Standards & Conventions
(20) Individual statements in pseudo code may correspond to individual
Fortran statements or to groups of Fortran statements, for
example:
(a) TABPT = TABPT - 1 (single Ftn stmt)
(b) Warn that extrapolation will take place (group of Ftn stmts)
We attempted to make our pseudo code sufficienty detailed that
inplementation in Fortran was simple. This helped eliminate
errors in the translation process.
2.0 TERMINOLOGY
Successful software development requires both clear thinking and
communication of ideas. If a member of a programming team is to function
effectively, he must have a clear picture of the problems he is to solve and
he must be able to express his solutions precisely. When team members
interact they must understand each other clearly. And the final product and
its associated documentation must be consistent and unambiguous. This
implies that a common set of technical terms and definitions should be used
throughout the project.
Unfortunately, in the software field there is considerable confusion over
the use of terms. The principal reason is that computer science is still a
young discipline. New terms are constantly appearing. The development of
Structured Programming Technology has spawned a whole set of jargon. There
is also confusion over the meaning of terms which have been in use for some
time, e.g. "module" and "segment." Many have poorly defined meanings and,
as a result, mean different things to different people. Obviously we need
to take the trouble to define our terms.
To minimize this kind of problem we developed a glossary of terms pertaining
to HSPF. It is located in Appendix I.
29
-------
Standards & Conventions
3.0 FORTRAN LANGUAGE
3.1 Use of ANSI Fortran
There are many different "dialects" of the Fortran language. Almost every
machine and operating system has its own variation on the basic theme. One
of the goals for this project was to produce software which would be
portable. That is, it should be capable of running on a wide variety of
machines. To do this we had to employ a common subset of the many dialects.
The American National Standards Institute has prepared a specification for
Fortran. Most manufacturers now implement all of its provisions, plus their
own extensions. We used ANSI Fortran together with a few specified
extensions.
3.2 Allowable Extensions to ANSI Fortran
We decided to extend beyond ANSI Fortran because, when applied to our type
of work, it does have some deficiencies which prohibit certain actions or,
at least, make them awkward. For an extension to be approved it had to meet
two criteria:
(1) There had to be a compelling reason for its use in this project.
(2) It had to be a feature incorporated in the Fortran languages of
all the following series of machines:
(a) IBM 360 and 370
(b) UNIVAC 1100
(c) HP3000
(d) PDP11
We figured that if these four systems all implemented a feature,
the chances were good that it was widely available.
The following extensions have been adopted:
(1) Direct access input/output
(2) Literals can be enclosed in single quotes e.g. 'ABC1 instead of
3HABC
(3) "Half-word" and "full-word" integers. Although ANSI Fortran does
not include two kinds of integers, differing in length, we did
make use of this feature. Unfortunately, it is not available in
the UNIVAC 1100 (and some other systems).
Adaptation of the code to a system which does not include extension (3)
would be quite costly; dispensing with extension (1) would be very costly.
30
-------
Standards & Conventions
4.0 FORTRAN CONVENTIONS
4.1 The Need for Conventions
In addition to defining a Fortran dialect for this software, we adopted a
set of conventions for its use. This was done because conventions promote
good programming style and, thus, help to produce consistent code. Usually,
there are many many ways of programming the solution to a given problem.
Although they may all use the same Fortran dialect and produce the same
results, one method may be vastly superior to another because it has a
clearer or better structure and implementation. While good computer
programming cannot be reduced to a set of rules, the use of appropriate
conventions can greatly assist in producing code which has a readily
observable and logical structure.
Serious communication problems can arise when several people work on the
same set of software. We wanted to produce consistent software, so it was
necessary to adopt conventions. However, this did involve a trade-off. Too
few conventions would have resulted in code which lacked cohesion; too many
would have caused programmers to feel stifled and restricted. We tried to
strike a reasonable balance.
4.2 Conventions Used to Implement the Standard Structure Figures
The Fortran language does not easily lend itself to structured programming
because:
(1) It is not block structured.
(2) It has an IF statement which allows conditional execution of only
a single statement on the true condition.
(3) It has a DO statement which allows for index condition testing
only.
Therefore, to write structured programs using Fortran, it is necessary to
use programming conventions. In this section the conventions for simulating
the following standard structure figures will be outlined:
(1) IFTHENELSE
(2) WHILEDO (called by some DOWHILE)
(3) DOUNTIL
(4) CASE
Because of the lack of block structure, it is necessary to use statement
labels in the simulation of the standard figures. These labels are given as
lowercase letters in the examples below.
31
-------
Standards & Conventions
4.2.1 IFTHENELSE
The IFTHENELSE figure tests a single predicate "p" to determine which of the
two function blocks (of code) F1 or F2 will be performed. The convention
for implementing the IFTHENELSE is:
IF (NOT.p) GO TO a
Code for F1
GO TO b
a CONTINUE
Code for F2
b CONTINUE
If the ELSE is not used, the figure reduces to the following:
IF (NOT.p) GO TO a
Code for F1
a CONTINUE
%
Statements within the clauses are indented two columns. The CONTINUE
statements terminate each clause and are coded in line with the IF. The
following may be used if the code for F1 is a single statement:
IF (p) Code for F1
4.2.2 WHILEDO
We call this figure WHILEDO and not DOWHILE (as in the IBM Manual) because
the name reminds one that the test is made before the block is "done." The
block will not be executed at all if the condition is initially false.
The WHILEDO is coded as follows:
C WHILEDO (optional additional explanatory comment)
a IF (.NOT.(p)) GO TO b
Code for F
GO TO a
b CONTINUE
The negation of the predicate "p" is used for clarity. In actual fact the
condition is generally given without the NOT. Statements within the figure
are indented two columns. The CONTINUE clause terminates the figure and is
coded in line with the IF.
32
-------
Standards & Conventions
U.2.3 DOUNTIL
The DOUNTIL figure provides essentially the same loop capability as the
WHILEDO, differing from it in two respects:
(a) The test of the predicate p is reversed. The WHILEDO terminates
when p is false; the DOUNTIL terminates when p is true.
(b) p is tested after each execution of F so that F is always executed
at least once.
This figure can be implemented in two ways:
Method 1
C DOUNTIL (optional explanatory comment)
a CONTINUE
Code for F
IF (.NOT.(p)) GO TO a
Statements within the figure are indented two columns. This method is
used if an index variable is not needed or if the predicate p is
complex.
Method 2
The standard Fortran DO is really a DOUNTIL based on an index. This
form has a cleaner appearance and will probably execute faster. It
should probably be used even if the index is only a counter.
DO a (index specification)
Code for F
a CONTINUE
Statements within the figure are indented two columns with respect to
the DO and CONTINUE. Note that, in terms of ANSI Fortran, the value of
the index is undefined after completion of the loop.
33
-------
Standards & Conventions
4.2.4 CASE
The CASE form is used to select one of a. set of functions to be performed
depending on the value of a variable aa with m possible values. The figure
is implemented as:
C CASENTRY I
k GO T0(a, b, ..., n), I
a CONTINUE
Code for case 1
GO TO t
b CONTINUE
Code for case 2
GO TO t
n CONTINUE
Code for case n
t CONTINUE
If it is possible for I to acquire an invalid value, the following code can
be added to check for that condition:
IF (I.GE.I.AND.I.LE.m) GO TO k
C Invalid value of I
code to handle invalid value goes here
Statements within each case are indented two columns. The CONTINUE
statements separate each case and are coded in line with the computed GO TO.
4.3 Other Conventions
4.3.1 Layout of program units
The Fortran language makes certain stipulations as to the order of
statements. We expanded on this, to improve the readability, clarity, and
consistency of our program units. The order of things is:
(1) A comment card containing the number of the program unit.
(2) The subprogram header line, including the list of dummy arguments.
(3) A comment of 1-10 lines which describes the general purpose of the
program unit.
(4) Type statements for the dummy arguments, if this is a subprogram.
34
-------
Standards & Conventions
(5) Any COMMON statements and their associated type and EQUIVALENCE
statements. The COMMON blocks are all named and occur in
alphabetical order.
e.g. COMMON/I0/INP,OUT,ERRMESS
INTEGER INP,OUT
REAL ERRMESS
COMMON/FARMS/A1,ALTER,COAT
INTEGER COAT
REAL A1,ALTER
(6) Local variable declarations and any associated EQUIVALENCE
statements.
(7) DATA statements, for local arrays and variables
(8) Statement function definitions
(9) FORMAT statements, ordered by statement number. Numbers of READ
formats start at 1000 and increase by 10; numbers of WRITE formats
start at 2000 and increase by 10.
(10) The executable code. This consists of the standard structure
figures and code segments implemented in Fortran using the
conventions described earlier.
4.3.2 Control structure of a program unit
Code segments may be composed of any of the Fortran statements. However,
the use of GOTO, IF, and DO statements is restricted to the standard
structure figures. Control flow enters a code segment at the top and leaves
at the bottom. No other entry points or exit points are permitted.
4.3.3 Type statements
The types and dimensions of all variables and arrays used in a program unit
are declared; none are typed by default. This is done using the following
statements, in the order given:
»
LOGICAL
INTEGER*2
INTEGER*^
REAL
DOUBLE PRECISION
35
-------
Standards & Conventions
Hollerith (character) data are stored as REAL variables with four characters
per storage location or as INTEGER*2 variables with two characters per
storage location. Note that this is designed to fit those machines on which
REAL variables occupy four bytes per storage location and integers normally
occupy 2 bytes. It is, to some extent, a machine dependent aspect of the
software.
The DIMENSION and IMPLICIT statements are prohibited, and dimensioning
information is never put in a COMMON statement. The identifiers in each
type statement are arranged in alphabetical order, except those which apply
to variables in COMMON, which are arranged in their physical sequence.
4.3.4 DATA statements
The DATA statement is used to assign initial values to a set of variables.
If it is used to initialize local variables in a subprogram, the compiler
ensures that this is done at the start of the run, not each time the
subprogram is called. This operation effectively makes the variable
"permanent." In a "stack" machine, such as the HP3000, this means that the
affected data items occupy space in the data segment throughout the
execution of the program, not just when the subprogram is executing. To
minimize the loss of working space in the data segment, we avoided the DATA
statement as far as possible.
4.3.5 EQUIVALENCE statements
Equivalence statements are used where it is necessary to do something
esoteric, such as overlaying sets of data which are not required by the
program at the same time, or the assignment of a single name to a group of
variables, which is done often in our data structures (Johanson, et al.
1979). The Chief Programmers regulated this type of work.
4.3.6 Identifiers
An identifier is the name of an array, variable, or subprogram. One of the
principles of Structured Programming is that "intelligent" identifiers be
used. That is, the characters in the identifier should convey to the reader
something about the data which they represent. In a large program there are
many entities which have to be named, preferably in such a way that the
relationships between them are indicated. For example, in a reach we might
have water, sediment size 1, sediment size 2, nitrate, DO, etc. Furthermore,
we might need to distinguish between the inflowing, stored, and outflowing
quantities of these constituents or we might need identifiers which indicate
the accumulated outflows of these constituents over the current day, month,
and year. Choosing identifiers for hundreds of entities, such that the
names are descriptive and fit into a systematic framework, is not a simple
task. Therefore, the assignment of these identifiers did not proceed on an
ad hoc basis, but a system was devised for naming the variables in each
major module. The variables used by each application module are listed in
the Programmer's Supplement (Johanson, et al. 1979).
36
-------
Standards & Conventions
4.3.7 Labeling of Statements
We label only those statements which need them. This improves the
efficiency of the object code. Executable statements are labeled in
ascending order, starting with 10 and increasing by 10, to allow later
insertion of additional labels without upsetting the order. Format
statements are collected at the head of the program unit. "Read" formats
come first starting at label 1000 and increasing by 10. "Write" formats
come next, starting at 2000 and increasing by 10.
Statement labels start in column 1 and are justified on the left.
4.3.8 Layout of executable statements
The executable parts of a program unit consist of the standard structure
figures arranged in a hierarchy which expresses the logic of the program.
The level in the hierarchy of any piece of code is reflected in its level of
indentation from the left margin. Long, complicated statements are avoided.
4.3.9 Continuation of a statement to more than one line
Because of the large number of continuation lines permitted by most
compilers, the continuation indicator is used to indicate the order of these
lines. We number them 1, 2, 3, etc. For declaration, CALL or subprogram
header statements, the final character on a continued line is a comma, to
eliminate doubt should the following line be missing. Similarly, for
assignment statements, the last character on a continued line is an
operator. Never is an identifier or constant broken at a line boundary. The
text on a continued line commences at the same level of Indentation as the
previous line (apart from the character in column 6, of course).
4.3.10 Mixed mode arithmetic
ANSI Fortran does not permit mixed mode arithmetic. We endeavored to
observe this requirement and use the "intrinsic" functions (see ANSI Manual)
to convert data to a uniform type where necessary. Mixed mode assignment
statements are permissible (see ANSI Manual).
37
-------
Standards & Conventions
4.3.11 Range and precision of numbers
In Fortran there are two types of numbers; INTEGER and REAL. Because these
numbers are represented in any machine by a finite number of bits, there is
a limit to the range of values which can be represented. With real numbers
there is also a limit to the precision with which numbers can be
.epresented. If a situation is encountered where a REAL variable cannot
represent a number with sufficient precision, the DOUBLE PRECISION option
can be used. The ANSI makes no stipulations as to the range and precision
-f variables used in Fortran; these factors vary from one type of machine to
another. We designed our program so that range and precision capabilities
of all of the following computer systems are met:
IBM 360 and 370
UNIVAC 1100
PDF 11
HP3000/II
Appropriate specifications are:
Type
DOUBLE PRECISION
REAL
INTEGER (2 byte)
Range Precision
From To (decimal places)
10»»-38
10**-38
-32767
10**38
10*»38
32767
16.5
6.5
exact
Note: In the case of DOUBLE PRECISION and REAL variables, negative and
positive ranges are the same. Only the positive range is indicated
above.
4.3.12 Use of comments
Comments are used to:
(1) Describe the purpose of a program unit (1-10 lines long).
(2) Separate and describe the function of logical blocks of code (code
segments).
(3) Clarify logic, where necessary, and to illuminate subtle points.
38
-------
Standards & Conventions
They are not used to provide a running commentary on the code or to explain
what is already obvious; e.g. C ASSIGN I TO Z. A well structured program
should not need a host of comments.
Comments appear immediately above the code to which they apply. They were
written at the same time as the executable code; these are usually better
than comments which are inserted afterwards.
Although good comments do nothing for the computer, care was taken over them
because they do help the person who has to read the code at a later stage.
Composing a good comment also helps the programmer to review the function of
that piece of code and, hence, to produce better code.
The text of a comment is indented to the same column as that of the code
which surrounds it, so that it does not impair the visual picture of program
structure conveyed by indentation.
Examples:
(1) c »»»»«»**»«»»««»«»»»»»****
C
C WARNING! This may not
C work on XXXX installation
C
c «*»««»»««*««»»«*»*«*«•*»»
(2) A = TEMP«A
C ACCA will be used in the 'Totals' Block
ACCA = ACCA + A
4.3.13 Transfer of data between program units
In Fortran, data can be transferred between program units in two ways: by
using argument lists or COMMON blocks. Argument lists are used whenever
feasible because they specify exactly those items of data passed between the
calling and the called program units. COMMON blocks, however, tend to grow
as software is developed because they usually serve several program units.
Typically, there are more data items in these blocks than any given pair of
program units requires. Thus, it is difficult to tell exactly which items
are being passed between a pair of program units. The development and
maintenance of programs becomes more difficult if COMMON blocks are used
where argument lists would have sufficed.
A further advantage of argument lists is that they permit arrays in the
called subprogram to have variable dimensions. Implementation of this
feature is discussed in the following section.
39
-------
Standards & Conventions
4.3.14 Argument lists
In ANSI Fortran, an item in a list of actual arguments can be any of the
following:
(1) a Hollerith constant
(2) a variable name
(3) an array element name
(4) an array name
(5) any other expression
(6) the name of an external procedure
There is little need to use types (5) and (6) and they are avoided where
possible.
Arguments can serve either or both of two functions: to transfer information
to the called subprogram (input mode) or to retrieve information from it
(output mode). Although Fortran does not require that arguments be
separated according to their function, we believe it is helpful to do so. We
arrange arguments in the following order:
(input,bbinput/output,bboutput)
or
(input,bb
input/output,bb
output)
bb indicates two consecutive blanks.
Examples:
CALL JUNK (JIN1,JIN2,bbFLAG1,FLAG2,bbJOUT1)
input input/output output
CALL WORKER (bbbbbbbbVALUE)
output
CALL LONG(INP1,INP2,etc (input)
2 MOD1,MOD2,etc (input/output)
3 OUT1,OUT2,etc) (output)
40
-------
Standards & Conventions
An argument list can be used to vary the dimensions of arrays in the called
subprogram. In ANSI Fortran this is done by including the name of the
affected array, together with scalar variables which indicate its
dimensions, in the argument list.
Example:
SUBROUTINE DOIT (A.MAXI,MAXJ)
INTEGER MAXI.MAXJ
REAL A(MAXI.MAXJ)
etc.
We make use of this feature when we have reason to believe that an array
will have its dimensions changed as the model is applied to different
watersheds or as it is implemented on different computer installations. In
this way only the highest subprogram in which the array occurs will need to
be recompiled when the program is reconfigured. This construction is also
sometimes used in a subprogram which is called by several other program
units. Thus, the dimensions of dummy arrays in the subprogram automatically
agree with the dimensions of the real arrays in the calling program units.
U.S.15 COMMON blocks
Only labeled COMMON is used. The COMMON statement contains a list -of
variable and array names. It does not contain dimension information; this
is put in the associated type statement.
Each COMMON block is immediately followed by statements which declare the
type and dimensions of all the variables in the block. The layout of the
type statements is described in Section 4.3.3.
Example:
COMMON/ABAC/HEAVY,LIGHT,LIGHTR.K1
INTEGER K1
REAL LIGHT(20),LIGHTR(10,2)
DOUBLE PRECISION HEAVYC10)
Because COMMON blocks usually serve several program units and because we
wanted to control their use carefully, the layout of data in COMMON was
supervised by the Chief Programmer.
4.3.16 General comments on programming style
Good programmers write code for the benefit of other people, not just for
the computer. They avoid the use of "clever" constructs which, although
efficient in execution, are obscure.
In general, if we were faced with a conflict between machine efficiency and
clarity of the code, our chosen solution favored the latter goal.
41
-------
Standards & Conventions
5.0 CONVENTIONS USED IN FUNCTIONAL DESCRIPTION
The primary purpose of the Functional Description (Part E) is:
(1) to describe the functions performed by the various subprograms in
more detail than can be achieved in the Structure Charts
(2) to explain the technical algorithms and equations which the code
implements.
Subprograms are described in numerical order in the text. This system
provides a logical progression for the descriptions. General comments
regarding a group of subprograms can be made when the "top" subprogram is
described, while details specific to an individual subordinate subprogram
can be deferred until that part is described. For example, a general
description of the PEHLND module (Section 4.2(1)) is followed by more
detailed descriptions of its twelve sections, ATEMP (Section 4.2(1).1)
through TRACER (Section 4.2(1).12).
6.0 METHOD OF DOCUMENTING DATA STRUCTURES
6.1 Structure of Data in Core
The way in which we arrange the variables used in our programs is important.
We structure them, as far as possible, using techniques like those used in
Structured Program Design. We try to group data items that logically belong
together.
Most of the variables in an Operating Module are contained in the Operation
Status Vector (OSV). The OSVs for the application modules are shown in the
Programmer's Supplement (Johanson, et al. 1979). The format used to document
a data structure is similar to that used to declare a "structure" in PL/1.
We do this because the technique is logical and convenient, not because of
language considerations.
6.2 Structure of Data on Disk Files
The HSPF system makes use of several different classes of disk-based data
files:
(1) The Time Series Store (TSS) is described in Section 2.0, Part E
and in Appendix III.
(2) The instruction files (OSUPFL, TSGETF, TSPUTF) and the OSVFL are
documented in the Programmer's Supplement.
42
-------
Standards & Conventions
(3) The information file (INFOFL), error message file (ERRFL) and
warning message file (WARNFL) are self documenting. One need only
list the file and read it to understand its contents.
7.0 METHOD OF HANDLING DIAGNOSTIC MESSAGES
HSPF makes use of two kinds of diagnostic message; error messages and
warnings. These messages are all stored on two files; ERRFL and WARNFL.
This system has at least two advantages:
(1) Because the messages are not embedded in the Fortran, they do not
normally occupy any core storage. This reduces the length of the
object code.
(2) The files are self documenting. They contain not only all the
messages, but other explanatory material. A user need only obtain
a line printer listing of the files to get an up-to-date copy of
this documentation.
Each message has been given a "maximum count". If the count for a message
reaches this value, HSPF informs the user of the fact. Then:
(1) If it is an error message,-HSPF quits.
(2) If it is a warning, HSPF continues but suppresses any future
printing of this message.
In addition to the above features, the Run Interpreter has been designed to:
(1) Stop if 20 errors of any kind have been detected. This gives the
user a fair number of messages to work on, but avoids producing
huge quantities of error messages, many of which may be spurious
(say, if the code could not recover from early error conditions).
(2) Stop at the end of its work if any errors have been detected by it.
Thus, HSPF will not enter any costly time loop if the Run
Interpreter has found any errors in the User's Control Input.
43
-------
8.0 REFERENCES
International Business Machines Inc. 197M- Structured Programming Textbook &
Workbook — Independent Study Program.
International Business Machines Inc. 1974- HIPO - A Design Aid and
Documentation Technique, Report GC20-1851-1. 130 pp.
American National Standards Institute. 1966. USA Standard Fortran, Standard
X3.9-1966. 36 pp.
Johanson, R.C., J.C. Imhoff and H.H. Davis, Jr. 1979- Programmer's
Supplement for the Hydrological Simulation Program - Fortran. This material
is on magnetic tape. See Appendix IV.
44
-------
Visual Table of Contents
PART D
VISUAL TABLE OF CONTENTS
General Comments 49
FIGURES
Structure Page
chart no.
1.0 Upper levels of HSPF system 49
1.2 Service subprograms available to the entire HSPF system . 50
2.0 Subroutine group TSSMGR 51
2.001 Service subprograms for subroutine group TSSMGR 52
3.0 The Run Interpreter 53
3.01 Service subprograms for the Run Interpreter 54
3.2 Subroutine group SEQBLK of the Run Interpreter 55
3-4 Subroutine group OPNBLK of the Run Interpreter 56
3.4.1 Service subprograms for processing input for
operating modules 57
3-4(1) Processing of input for the PERLND module 58
3.4(1).10 Processing of input for the PEST section of the
PERLND module 59
3.4(1). 11 Processing of input for the NITR section of the
PERLND module 60
3.4(1).12 Processing of input for the PHOS section of the
PERLND module 61
3.4(2) Processing of input for the IMPLND module 62
3.4(3) Processing of input for the RCHRES module 63
3.4(3).2 Processing of input for the HYDR section of the
RGHRES module 64
3.5 Processing of User's Control Input which deals
with time series 65
3.5.01 Service subprograms for the TIMSER subroutine group ... 66
3.5.2 Processing entries in the EXT SOURCES Block 67
3.5.2.2 Subroutine group EXTTS 68
3.5.2.3 Check and expand a reference to an Operation time series. 69
45
-------
Visual Table of Contents
3.5.3 Processing entries in the NETWORK-block 70
3.5.4 Processing entries in the EXT TARGETS Block 71
3.5.6 Subprograms involved in allocation & deallocation
of INPAD rows 12
3.5.8 Subroutine group TINSTR 15
3.5.8.01 Service routines for subroutine group TINSTR '4
3.5.8.2.3 Subroutine group PINIT ]?
4.0 Operations group of modules
4.01 Service subprograms for the Operations Supervisor .... ''
4.1 TSGET '*
4.1.01 Service routines for TSGET and TSPUT '*
4.1.1 The GETTSS section of the TSGET module °^
4.1.2 The GETSEQ section of the TSGET module °^
4.2(1) Pervious land-segment application module
4.2(1).2 The SNOW section of modules PERLND and IMPLND •"
4.2(1).3 The PWATER section of the PERLND application module ... °^
4.2(1).4 The SEDMNT section of the PERLND application module ... °5
4.2(1).7 The PQUAL section of the PERLND application module. .... °°
4.2(1).8 The MSTLAY section of the PERLND application module ... 87
4.2(1)=9 The PEST section of the PERLND application module .... 88
4.2(1).10 The NITR section of the PERLND application module .... 89
4.2(1).11 The PHOS section of the PERLND application module .... 90
4.2(1).12 The TRACER section of the PERLND application module ... 91
4.2(1).13 Subroutine group PPTOT 92
4.2(1).14 Subroutine group PBAROT 93
4.2(1).15 Subroutine group PPRINT 94
M.2(1).15.1 Subroutine group PERACC 95
4.2(1).15.2 Subroutine group PERPRT 96
4.2(1).15.3 Subroutine group PERRST 97
4.2(2) Impervious land-segment application module 98
4.2(2).3 The IWATER section of the MPLND application module ... 99
4.2(2).4 The SOLIDS section of the IMPLND application module . . . 100
4.2(2).6 The IQUAL section of the IMPLND application module. ... 101
4.2(2).? Subroutine group IPTOT 102
4.2(2).8 Subroutine group IBAROT .... 103
4.2(2).9 Subroutine group IPRINT 104
4.2(2),9.1 Subroutine group IMPACC 105
4.2(2).9.2 Subroutine group IMPPRT 106
4.2(2).9.3 Subroutine group IMPRST 107
4.2(3) Reach/mixed reservoir application module 108
4.2(3).1 The HYDR section of the RCHRES application module .... 109
4.2(3).4 The HTRCH section of the RCHRES application module. ... 110
4.2(3).5 The SED section of the RCHRES application module Ill
4.2(3).7 The RQUAL section of the RCHRES application module. ... 112
4.2(3).7.1 OXRX subroutine group. DO, BOD simulation 113
4.2(3).7.2 NUTRX subroutine group. Primary inorganic
N and P simulation 114
4.2(3).7.3 Simulate plankton populations and associated reactions. . 115
4.2(3).7.3.3 Simulate phytoplankton. . 116
4.2(3).7.3.5 Simulate benthic algae. ... 117
46
-------
Visual Table of Contents
4.2(3).7.4 Simulate pH and carbon species 118
4.2(3).8 Subroutine group RPTOT 119
4.2(3).9 Subroutine group RBAROT 120
4.2(3).10 The RPRINT section of the RCHRES application module ... 121
4.2(3).10,1 Subroutine group RCHACC 122
4.2(3).10.2 Subroutine group RCHPRT 123
4.2(3).10.3 Subroutine group RCHRST 124
4.2(13) The Display utility module 125
4.2(14) The Duration Analysis utility module I26
4.3 TSPUT 127
4.3.1 The PUTTSS section of module TSPUT 128
4,3.1.3 Subroutine group FILTSS I29
47
-------
Visual Table of Contents
GENERAL COMMENTS
The entire set of structure charts for the HSPF system appears on the
immediately following pages, forming a visual table of contents for the
software. For a discussion of the conventions followed in compiling the
charts, refer to Section 1.2 of Part C "Standards and Conventions".
Note that the following material gives a complete picture of the system. It
starts at the highest, most general, level and proceeds down to the lowest,
most detailed, levels. A user need not assimilate it all before using HSPF.
Initially, one may not wish to proceed more than two or three levels down
the structure tree. Later, having become more familiar with the system, one
may wish to explore it in more detail.
In general, a subprogram calls only those subprograms situated immediately
below it. However, "service subprograms" are an exception to this rule.
These routines perform elementary tasks, usually for more than one calling
subprogram. Therefore, they are arranged in groups; each group is situated
immediately above the subprograms which its members can serve. For example,
the subprograms in Structure chart 1.2 can be called by any subprogram in
HSPF; those in Structure chart 3.01 can be called by any subprogram in the
Run Interpreter (Structure chart 3.0). This arrangement makes the structure
charts more compact and readable and will assist programmers who need to
design overlay structures for the code, or partition it in some other way.
48
-------
USRRDR.
MAIN
preprocees
user's control
input
1.1
TSSMQR
manaqe
lime Series
1.0
\
Z.O
)
1
provide a
system for
operating on
time series
1.0
MTERP
1
1
service
subprograms
'-^ J
u >
OSUPER
interpret
run insts. in
User's
Control
Input
3.0
3
supervise
and peribrm
operation*
4.0
1
4.0 \
Structure chart 1.0 Upper leveb of H5PF
49
-------
OAYVM.
Accvec
linearly
interpolate
a value for
this day
I.Z.I
TRNVBC
accutn. the
elements of
on* vector
to another
I.I.Z
service
sub-proqs
servec
perform a
linear
transform.
on the
element* of
a vector
l.t.»
1.1
ceisus
set all the
elements of
a vector to
a specified
value
l.t.4
FAHREN
convert a
temp in
deqF
deqC
to
l.t.S
BRROR
convert a
tamp in
dea.C
to deqF
1.26
produce on
error
message
1.1.1
WARN
produce a
warning
U1
o
ire.
FNDKWO
locate a
qiven
keyword in
a keyword
library
l.t.9
DAYMON
find the no.
of day* in a
qlven year
and month
I Z 10
PUMPtR
dump a
section of
the Users
Control
Input
1.2 11
VW-MO_J
determine
whether a
qiven value
is in
VALLIB
1.2.12
eerosv
qet OSV
fronn
OSVPL
>MTO«V I
pot OSV
into
OSVFL
IZ.Il
l.t.14
OSVRD
qet o
record
OSVWR I
write a
record
I.Z 13.1
I.2.I4-.I
EXOA.TE
KBUPP
Convert
date/time
to external
•format"
I
Z 15
WBOFF
read Tse>
record
into
buffer
iaifc
write
-toTS
LPVEAR
ir
i>
i z .n
Set
ueap
-roip
year
lilft
MBXT
Computie.
time -from
sfart of
base year
I Z 11
PYPMON
tOTSS
buff*'"
IZ.10
OlFTIM
•find the no.
of day«
(n o qiven
year
month
I-Z.ZI
4iff betw
t-yyo ddte/
times
i.z.ia
Structure chart 1.2 tiervice subprograms available to the entire H5PF system
-------
T9SM6R
Ul
adds
datasets
to th« time
series store
reads the
user's input
stream
Structure chart ZO Sobroutine qroup TSSMGR-
-------
ECHO
echoes the
user input
streom
1.001
ROKVND
reads key-
word* from
the
information
file
Z.OOZ
CHK-STR.
service
routines
t.OOl
checks the
validity of Q
user Input
string
2.00J
PUT
write* an
array of a
given length
toTSS
Z.OO4-
1
CHKINT
check* the
validity of a
user input
number
2.00S
SET
reads an
array of a
given tenqth
from T5S
Z.006
Ul
N)
PRTML
prints a
dataset label
to the line
printer
2.007
i
»UTSPA
output space
information
about the
T5S
Z.000
PROCHK.
classifies a
chunk of
free space
according to
Its «iie
£.000.1
SETUP
sets up Q
library of key
words for
ADqUPOWE
$ SCRATCH
command*
2.009
i
KDMKT
reads the
next keyword
from the
user Input-
stream
roio
FNDRMO
finds the last
keyword at
o qiven level
in ihe u?er
input streom
2010.1
i
•MOCHK.
finds the first
fitted chunk.
of free
space on the
T5S
2.OII
Struchire chart 2.001 Service subprograms 1br Subroutine
qroup TSSMGR
-------
INTBRP
Interpret a
RUN Data Set
in the UCI
5.0
6LOBLK.
SBQftLX.
OJ
process the
GLOBAL
Block
3.1
r:—. 1
Service
Subprogram*
**^l I
5.oi
TAMLK.
process the
OPN 5EQ
Block
3.Z
OPMBLK
locate tobies
in the
F TABLES
Block
3.3
TIM5ER
process
input for
each
operating
module
3.4
process
blocks which
deal with
time series
vurite the
ppi Supervisor
instrflle
3.6
P5PCCL
process
Specidl
Acti on
instrs.
3.7
3.2
3.4
35
Structure chart 3.0 The Run Interpreter
-------
Ln
service
sub-proas for
INTERP
group
3.01
YRMtNft
STD/CTE
find the offset
of a date/ time
from the
*Wt of the
year
3.OI. 1
BWDATE
check validity
of a startinq
dote /time
$ concert
internal
format
to
3.01.2
LOCBLK
check validity
of an endinq
date/ time *
convert to
internal
format
3.01.3
OPNNO
locate ihe
blocks of text
in a given
ranqe of
text
3.01.5"
find OPMNO
of a qiven
opn_ id in
the OPNTA&
5.01.6
Structure chart 3.01 Service subproqrarrvs for the Run Interpreter
-------
S60BLK.
Ui
Ui
KIBVNOPN
a new
operation
3.t.l
NEWIN6
1
3.E/ZL
i
i
1
1
4EWWG
a new
INGROUP
3.Z.2
46VSJEY.G
!
3es
process the
OPN
SEQUENCE
Block
MevviexQ
3.Z
a new
EVGROUP
3.Z.?
i
c
1
2ND1NG
EHOBX6
end of
INGROUP
iwoexa 1
i_
3.Z.4
— '
3.2.5"j
end of
EV.6ROUP
3.E.5
Structure chart 3.a Subroutine qroup 5EQ.BLK of the Run Interpreter
-------
OPNftlK
Ui
process input
for each
operating
module
— —
3.4
— -
process input
for all
operations of
ihe same -type
3.4 (NO
I service i
-^>' subprograms
' -T~-J
$.4:1
process
genera (
Input
3.4-CN).!
process
input for
section 1
5.4CN).Z
process
input for
section K-l
3.4(W).K
This substructure is duplicated for each operating module (utility or application
module ) in the system.
The numbering
shown is that -for the Nth operatinq module..
I
Structure chart 3.4 Subroutine ajoup OPNBLK. of the RUKi Interpreter
-------
Cn
ITAJ&LE.
process doto
from a table
containing
integer
values
3.4.1.1
TA&RBC
find the record
(in a table)
containing
info, for a
specified
operation
3.4 U.\
service
subprograms
for processing
input to
operating
modules
RTABL.6
process data
from a table
containing
real /values
3.4.1.
process data
contained
in an
Ftable
3.4.1.3
TA&RBC
I
I
I
L_ J
Stroctore chart 3.4.1 Service 5ubproqram5 for processing input for operating modules.
-------
Ul
oo
PPERUM
process input
to PERLWD
module
3.4-(0
1
PPGEKI |
process the
general
Input
3.4(0.1
PPWTQS |
Section
PWTGAS
3.4 (O.T
>TRAC |
Section
TRkCER
3.4(0.13
1
,
PATEKAP |
Section
ATEMP
3.4(0.£
PPQUA.L. |
Section
PQUW-
3.4CO.6
PBRR5T I
<~ ~\
set f lux.
accumulators
to lero
4.2C17.15.3
PSK10W |
Section
SNOW
3.4(03
PM4TU |
Section
MSTLW
3.4 CO. 9
1
?PWATR |
Section
PWATER
3.4(0.4
PPE6T |
Section
PEST
3.4C0.1O
1
3A.I Y\ \f\ \
*r\.\i. \O y
\
WEDMT |
Section
SEDMNT
3.4(0.5
PNITR |
Section
N\TR
3.4(0. l\
1
2 A- t \"\ It \
5.*ni).ii y
i
*5TMP I
Section
PSTEMP
3.4(0.fc
PPHO* |
Section
PHOS
3.4(0. It
1
a .A. / i^ 19 x
3.^-CU.It. N
Structure chart 3.4CI). Processing of input for the PERLND module
-------
PPB5T
process
input for
Section
PEST
3.4-C0.10
SOLDAT
vo
qet data
on soil
properties
3.4-CO.io.l
process
parms reqd
for I*t order
kinetics
3.4CO.IO.Z.
process a
•table of
order parms
3.4C0.10.2.1
•SVALP
parms
forsinqle
valued
ads/des
procfi95 a
table of 5inqle-
valued
ads/des
parm*
NON.6VP
porms reqd
for non-
sinqle valued
ads/des
3.4-C 0.10.4.
ptpceee.
ini-Hal
process a
•table of
ads/de-5
parms
3.4-6). 10.4.
Structure chart 3.4CO.IO frocessinq of input for the PE^T section
of ttTe PERLKD module
-------
ON
O
50UDAT _[_
l~ "1
1 1
1 1
1 1
1 1
1 1
|_ 3/KO-KM J
PLNTPKft
process
plant uptake
parms
3-4X0. U.I
process
input for
'Section
KUTR
34 CO. 1 1
PIR6.TP
~1
1
5VALP I
i r ~i
1 1
i i
i i
1 1
1 1
i i
3.4-CO-'03
l_ ' _J
r
1
i
i
1
1
[_ 5.4-1
-1
-Structure chart 3.4CO.H Pfoce^ing of input -for-the
•section erf module
-------
PPHOS
50LDAT I
r
1
WJ4TPI4 I
~i r ~i
1 1 1
input for
-section
PH05
3.4C0.1Z.
EHRSTP
r ~i
WWP 1 sroRae 1
r ~i r "i
1 II 1
V •
I I
I I
I I
I I
I I
I I
I I
' 3.4CO.UM I
[_ _ _3ACO.Ua_ J
[_ 5.4 CO .10.3- J
-3trviciure chart 3.4-CO.IZ Proce^sinq of input -for 1he PH06 median of
module
-------
PIMPLN
K)
PIWATR
i
Section
IWWTER
3.4(Z).4
FMOEKl
process
general input
P5OLID
Section
SOLIDS
tt)5
i
i
process
input to the
IMPLND
module
3.4(2)
?ATEMP
... 1
Section
KTEMP
3.40XZ
L_ . 1
PIWTQS
Section
IWTQAS
i
i
!
PSNONW
1
Section
SNOW
3.4(0.3
PIQUAL
Section
\QUAO-
IMPR5T |
i set flu* i
accumulators
to zero
1 1
Structure chart 3.4(2") Processing of input for the IMPLND module
-------
ON
PRCHRE
process input
to the
RCHRES
module
3.4 C3)
PRG6M |
process
general input
3.4(311
c
>HVDR |
Section
HYDR
1
i
POXRX
(
7ADCAU |
Section
ADCALC
3.4(313
PFSTO
1
Section
FSTORD
3.4(317
1
Section
3.4t5).8.l
l
l
PCONS
1
Section
COMS
3.4C314
PRQUAL |
Sections
C*RX -thru
PHCAfcB
RNUTK*. |
Section
NUTRY
l
l
]
PPLANK.
PHTRCW | P5EO |
Section
HTRCH
3.4(31?
Section
SED
3.4(31.6
KCHRST 1
-1- - — |
set flux
accumulators
torero
I
Section
PLANK
3.4(3>ft3
1
PPUCM*.
Section
PHCA.RB
3.4(3>.8.4
Structure chart 3.4(5). Processing of input for the RCHRES module
-------
PHY&R
DISCH
compute
initial
discharge
rates
3.4(3>.2.l
process
Input for
Section
HYDR
3.4(3)Z
FMOROW
~l
AOXIL 1
r
!
1
[_ 4.2(3). 1.3
Structure chart 3.4(3X2
Proce*sinq of input for the HYDR 5ection
of the RCHRES module
-------
TIMSBR
Ui
f ' I
process service |
blocks which subprograms i
deal with \; [
timeserie* X< [
3.5 | 3.5.01 i
" T
3.5.01 y
CHKTi* QPFINT SRCBLK WETBLK TAR8LK SAME M.LOC OSVPPT TINSTR I
look at the initialize process process the process the chain entries allocate
T53 TS&ET/PUT EKT NETWORK. &X with same rows in
descriptor instr. files SOURCBS Block TARGETS VQlue INRW>
Block Block
3.5.1.1 3.5.1.1 J.'j.fc 3.5.3 35.4 3.5.5 15.6
assiqn
values to
flatj-ptri.
inOSVs
35.7
qenarate
final
F56ET/PUT
instructionf
3.5.0
| OSVPRO
3.5.z \ 3.5.3 y m y 3.5.6 y
process a
group of
entries in
WORKFL
3.5.71
3,5.6 ^>
Structure chart 15 Processinq of User's Control Input which deals with time series
-------
subprograms
3.5.OI
WKDMP1
dump time
series linkage
info from
WORKFL
3.5.01.1
dump
primitive
TSGET/T5PUT
instrociions
from WORKFL
35.0 l.fc
Structure chart 3.5.OI Service subprograms for ihe TIMSER subroutine group.
-------
TOPTMO
EKTTS
check
multiple
tarqet
reference
3.5.2.1
check and
expand refs.
to ev-t time
series
3.5.2.2
SRCBLK
process EXT
SOURCES
Block
5.5.
OPMTS
check and
expand a
ret to an
opn time
series
3.5. E.3
3.512.5
PAIRS
match source
and tarqet
time series
3.5.ZA
CHAIN
chain entries
inOPN
SEQUENCE
order
55.Z.5
Structure chart 3.5.2 Processing entries in the EXT SOURCES Block
-------
6-XTT&
check/
expand
reference*
to external
•hme series
TSSDS
seoos
check/
expand
references
toTSS
dcrhiset-
35-Z.2.I
check/
expand
references
to a
sequential
file
DSCHK.
PROMCKA
check
validity of
dataset
reference
member
3.9.2.2.1.2
R6UFP
read a
record -from
T5S
I
i :__j
Structure chart- 35.Z.2 Subroutine group EXITS
68
-------
OPMT&
chk £ expand
a ref to an
Opn. time
series
NA6MTS
proc«s
-------
KIETBV.K
o
TOPTNO 1
1 1
1 1
1 1
1 1
process
NETWORK
Block
3.5.3
1
OPKITS PAIRS CHAIN I
I ' 1 1 1 1
II II
II II
II II
•"1
-
3.5.5.4 j
Structure chart 35.3 Processing entries in the NETWORK Block
-------
process EXT
TARGETS
Block
BXTTS <
r 1
1
1
I
1 '
1 35.2.2 |
3PKJTS J
35Z3
WMRS
r -
1 1
I 1
1 ,
1 1
1 35.2.4- 1
CH/UN
r "
1
I
1
1
1
.1.
35.2.6
chart 35.4 Proceesinq errtries in-the EYT TARGETS Block
-------
ALLOC
TARGET
process
entries,
considering
opn os a
target
35.6.1
allocate
INPAD
rows
3.5.6
SOURCE
process
entries,
considering
opn as a
source
ROW
—
allocate
first
available
row
3.5.6.1.1
*
rsvNS
— -
generate a
primitive
T56ET/TSPUT
instruction
3.5.6.U
CHANGS
— -
chanqe
subseq.
references
to a source
or target
3.5.6.1.3
release
rows which
hove been
so f laqqed
3.5.6.3
I
qroup
Structure chart 3.5.6 Subprograms involved in allocation $ deallocation of IWPAD rows
-------
•polrvte
CO
lRlNT
date/time
of-first
M.
T1NSTR
service
i 1
L- T-J
IT-SB.
fird
-Hme
PINSTR
3.5-8.01
V
PNtT
FUNT_ L
L 3.^.0.1.
—
TINOOrl
output a
TS6-E.T/PUT
3.6.6
~r/NS~TR
-------
Service
routines
SVv/ORP
Get o
w/ord
PVAU
•Prom
the TSS
data set
3S.9-OI
VROFP
PUT one
or
more.
values
into TSS
&&*«
Compute-
A rt • |_
-f-rom
- of
3.5-6.03
Structure, Chart 3.5.0-01 service routines.-for
-------
FLFRKl
ev.
Ul
TSB
-for
PIN IT
wilue*
r
L_3L&fttXlj
AlUNC
APP
S.S6.23.2.2.
ADJUB
Of
APP
3.g.g.2.3.Z.3
P8
51NTJ
L _ _?^ft-LJ _l
fnend-
PIMIT
-------
TSGET
cjrt reqd.
input- -Hme
series 4.1
supervise 4
perform oper-
4.0
I
i — •
j service
->l
L___4£!_J
perform an
operation
4.2(N)
E> i
I
Operating
Module, *adr|
Operating
Module,
4.2(Nl
Operating
ect. k
4.2 (NVk
•This substructure is duplicated -for each operating module
application module) in the ^stenn.
The numbering Astern shown is-that tor the M*h oper
-------
service
subprograms
&A.LCHK
perform
ADDT1M
a
material
balance
check
4.O\
SPECL ,
increment
the dote
and -time
4.OZ
*
perforrr\
Special
Actions
4.03
Structure chart 4.01 Service subproqrarrv* for the Operation* Supervisor
77
-------
bu-Pfcsr
I fWJ-
6ETTSS
T55
-4.1.1
F1LLWS
4-1.1.2.
-to
4.I.I-Z.I
BC^I
I
The
of
80
-------
oo
read
card
Aim
Hydro <£>cv\p
4.I.ZJ.Z
Structure
CHK.SE.GL
-f.1
service
routine
4.U-
MOCRP
read
HVdr
sem-
4UI.4
MDVAI-
4.12./.4-I
FIFCRP
15
4.1.2.U5
red A
\
5
4.1.2.1.6
The
or» of -the
-------
oo
NJ
Perform
computations
on a segment
of prrviouc
land
4.1(1)
1
ATEMP SNOW
Correct air
t*mp»rature
for elevation
difference
4.1 CO. 1
PWAJER
Simulate the
accumulation
and melting
of vtow and
ice
4.1(1). 2
4.zaY7>
SEDMNT
Simulate
water
budget
for pervioos
land segment
4.1 CO
I
m.3\
2
Produce
and
P5TEKAP
remove
sediment
4.1(0.
C0.4)
4
PWTGAS PQUAL
E<»tiniote
temperatures
4.1(l1.5
Ultimate simulate Duality
water coovtrtucflt?
temperature otiog simple
and dissolved 9*1 rela-liomhipf
c»ncentratior« tai*>> s»iitne«t <•
water yidJ
1
4.2dVV>
•"" 1
&.e>R\- CUE.VA\CAL SSCTIOKS
WVSTLM
estimate the
moisture < the
PEST
A
Simulate
pesttcide
behavior in
detail
VPTOT P8AROT
PVace point-
valued output
in INPAD
4.aciyn
4.7.(0.9
tt(,V»>
MiTR
ilats
nitrogen
betv
ivior in
detail
PPRINT
Place bar-
valued output
in
INOA.D
4.2(0.
14
4.1(1"). 10
|
4.a(0.,o>
Produce
pri
ntc
&
ou-t^ut
4.2(0.15
PWOS
Simulate
?ho
sohoros
behavior in
detail
4.-z(iYu
\
•P^CER
Simula-re, the
moM
ement of
a tracer
Cca
iKervafiNe")
4.a (ft. a.
4.2(,y)>
i
chart 4.1 (0 Ptrviout land-segment application module
-------
Simulate
and
ot
at\d ice
oo
CO
METEOR
EFFPRC
Estimate
meteorological
conditions
4.2(0-2.1
Determine
effect of
precipitation
on the pack
4.2(0-2-2-
Corn pad:
4.2(0-2-3
frum the
pack
4.16V2.4
COOL6R
hart
meH-
4.26)2.?
Simulate
loss of tieat
•Gram -me
pack
4.2(0.2.6
Calculate
4-.2(0.2.4.1
Warm 1lie
pocW if
possible
4.2 (A.2.7
UQUIP
Melt the paalc
wing any
remaining
Handle liquid
4.2(1) .2.
6M6LT
otcurence of
pack ice
4.2 6) .2-10
Wlelt -the pack
iftlno, heat
•fVom ihe
4.2(0.2-.! I
state
variables
wheh snow-
pack disap-
pears
chari 4.2 (l).2. Trie 5WOvg Action
mcxioles- PE£1ND and
-------
PWATER
Simulate \Nater,
for pervioos land
segment 4.xCO.3
ICEPT
3URFA.C
UIOME
HONE I 5WAT&g I EVAPT I
simulate
interception
4-.20Y3.I
oo
water available
for i
and runoff
\orter-flovN
Simulate
upper zone
behavior
4.2(05.4.
of
UIIMF
FROUTE.
Compute
uvf low IP
upper zone
lower zone
behOMior
qroond-
\Nate-r
behoviotr
Simulate
evapotwns -
fira-hon (6T1)
4.2O3.7.
ENICW
ET
-from
fcaseflc
4.2Ci:mi
trom
tion
-i.-JLC03.ia
rone
_L
uja-fer
41(03.74
tooer
zone
4.2CO.S7.<
determine
^ur-face
rimof-f
42(l>5.2.»4
ETUIS
upper zone
subrou-Hne
4.2(0.34.1
Structure Chart 4.2CH.3 The PWATER section 0*
application module
-------
SEDKAMT
fVoducc
and
sediment
4.1(0.4
DETACH
oo
01
Detach so\l
by rainfall
4.10V4.1
/Utermritw Subroutines
surface flow
using
method 1
4.1(0.4.2
ATTACH
by
4.1(0.
Simulate
He-attachment
of de-tached
sediment
4.1(0.4.4
Structure chart 4.1 OY4 the S6WANT SediOH off
application module
-------
PQUAL
Simulate
constituents
osing simple-
relationship?
with -sediment
and water yield
QU/M.SD
oo
Remote by
association
with
sediment
Accumulate,
and remove by
a cDmtaO"t unit
rate and t>y
flow
OU&LIP
Simulate by
concentration
in interflow
Simulate by
concer>tra-tion
in
round water
4.1 ro -7-4
Structure chart 4.2 dVl The ^>Q\JM section
PERLWD applicntion module
the
-------
NASTIA»Y
estimate mois-
ture storaqes
$ *dute fluxes
Q% fractions of
stored material
4.2 (0-6
TOWAY
for
topsoil layers
subsoil lasers,
chart 4.^ CO.ft The MSTLKf section
-------
9E3T
simulate
pesticide
behavior in
detail
AsGRGET
required
aqri-chem
section*
1
series
^ i
ncal
0-4-1
SDFRAC I
calculate trie
fraction of
the surface
layer that i«
eroding
4/z.ayq.z.
SBDMOV, [
move the
chemical on/
with sediment
4.ZC0.9.3
TOPMOV |
move solutes
in the -topsail
(surface £
upper lasers)
t update
storages
4.ZC0.9.-*
PfiTR%N |
perform
reactions on
pesticide
4,1 C 0-9.5
dU&MOV
move solutes
in the sub-
surface lovers
(lower i qroond-
water) 4- update
5ii>ra9e5
4.2(0-9.6
00
oo
PIRORD
SV
NONSV
adsorb/desorb
using first
order
kinectic*
a4«orb/desorb
usinq the
einqle value
Freundlich
method
ITER
1
adsorb/desorb
usino the non-
sinqle value
Freundlich
method
ITER
iterate for
adsorbed $
solution value
on Freundlich
isotherm
I
I
degrade
chemical
Structure chart 4.Z(0.9 The PEST section
of the PERLND application module
-------
MITR
nitrogen
behavior in
detail
i
see PEST
«ecV\on
^&ftC*D^k£ i
9W* *^*^™« i
l|
see PEST
secHon
co
[ -see PEST
1 section
TOPVIOV
see PEST
section
perform
reactions on
ni-Vroqen
foctns
3UBMoy_|
seciton
PIKORP i
i •-
!«*ePE5T
section
SV
1
see PEST
sedlon
^-iM-ii
BR I
PEST
section
.3.61
Structure chart 4.2 (I). IO The NITR
eecKon of ihe PERLND application module
-------
PH08
simulate
phosphorus
behavior in
detail
AGR8KT 1 9DFR>>C 1 8EDMOV
___ir i__!p
J-,
TOPMOV 1
r
, *ee PEST 1 see PEST , see PEST i ' see PEST
. section <
1 i
i i
i I
1 1
1 1
i |
section section
1 section
1
1
1
1
[ 4.XCO-^-IJ j_ 4.Z£lV*l.lJ L 4.1CO-«I.3_| |_ 4,Z(l).q.4J
PHOWtt*
perform
reactions
SU5WVOV }
on
phosphorus
forms
MI.I
i
1 see PEST
1 section
1
1
1
i
L 4.2_ \_
see PEST
section
i_ _____ t_ __ _j
I sec PEST
i section
I
1TER. j
j~»ee PEST
1 section
Structure chart 42(0.11 The PROS section
of the PERLND application module
-------
TRACER
simulate
movement
of a tracer
4-.ZO).IZ.
1
TOPMOV 1 SUdMOV 1
r- -•• -j [- -•-
see PEST I see PEST
1 section 1 section
1 I
: i
4..Z(I>*.*
4..iO)q.t>
Structure chart 4:2(1). 12 The TRACER section
of the PERLKJD application module
91
-------
vo
to
Section
ATEMP
4-.Z f 0.13.1
5K10W
PPTOT
place point-
valued time
n
PWATER
4-.2CO-B.4-
P5TEVAP
•Section
PQUM-
MSTIJPT
Section
M.ITR5T
W05PT
Section
WITR.
TRACPT
Section
PHOS
5eciion
5trudure chart 4:2 C0.6 5ubrovriine qroup PPTOT
-------
P6AROT
place bar-
tt
in
VO
U)
5K10W
•section
-5BDVAMT
?QOM_
lection
4.2 (i>. 14-. 8
Action
*>Vrudure chart 4
5ibrou-t1ne
?BA»Ror
-------
PERACC
accumulate
fluy.es for
printout
4-.ZCO.IS.I
4-.ZCD.I5.I
PPR1NT
perform
materials
balances fc
produce
printout
4.ZCD.15
PERPRT
mat
balances
write to
printer
PERRST
reset
$ state
variables
Structure chart 4.1 (I). 15" Subroutine qroup PPR.INT
94
-------
PERA.CC
accumulate
fluxes for
printout
*. I
SNOACC
accumulate
ftuxet-for
section
SNOW
IS II
PWAKCC |
accumuVatg
flukes for
section
PWKTER
SDACC |
accumulate
flukes for
section
SEDMNT
PWQACC |
accumulate
fluxes, -for
section
PWTOA.5
PQACC
PSTA.CC
accumulate-
fluxes for
section
PQUAL
NtTACC
accumulate
fluxes for
section
PEST
4..ZCO.I5.1.6
PHOACC
accumulate
fluxes for
section
NITR
TRACC
accumulate
fluxes for
section
PHOS
accumulate
fluxes for
section
TRACER
4-.2.CD.I5.I.9
Structure chart 4.Z(O.I5.\ Subroutine group PERACC
-------
PERPRT
VD
AIRPRT I
for module
section
ATENAP
«(«.«*.,
PQPRT |
tor module
section
44C11.IW.7
SNOPRT |
for module
section
SNOW
4.ZU).»*.Z.Z
MSTPRT |
for module
section
MSTLAV
4.ZID.I9.Z.8
PWAPRT
perform
materials
balances t
write to
printer
i
1
for module
section
PWATER
4-.ZO)
.1*2.3
PESPRT
for module
section
PEST
*.zcn
-15.Z.9
SEOPRT |
tor module
section
SEDMNT
4.ZOVI5.Z.4
NITPRT ]
for module
section
N1TR
4.ZCD.I52.IO
PSTPRT |
for module
section
PSTEMP
**cn.ww
PHOPRT |
for module
section
PHOS
4..ZO).K.Z,,
PWTPRT |
for module
section
PWTGAiS
4.Z(O.I9.Z.fr
TRAPRT 1
for module
sechon
TRACER
4-.zco.,^,z
Structure chart4.Z(O.I5.Z Subroutine qroup PERPRT
-------
PERRST
reset flux
t state
variables
4.1(11.17.3
VO
SNORST
reset for
section
SNOW
44 <». 1*3.1
PVNARST
react for
section
PWATEfc.
4.1(1). 15. S.Z
SOWST |
reset for
section
SEDMNT
4.t(l).IT.».3
PWQRST
reset for
section
PWT6AS
4.ZCD.IV.3.4
PQRST
PSTR*T
KIITRST
PHORST
TRRST
reset for
section
reset for
section
PEST
4.2. (0.15.16
reset for
section
NITK
section
PHOS
reset for
section
TRACER,
Structure chart 4.Z(t)tl5.3 Subroutine group PERRST
-------
Perform
computations
on a segment
o4 impervious-
land
4.1 (21
ATEMP | SNOW 1 IWWER
1 11 i
1 II 1
| (see module 1 1 ( see module |
| PERLKID 1 1 1 PeRLNDl |
1 II 1
SOUDS
vwater budget
for Impervious
land segment
4.2. (U^
I
\
4.1 (1\5>
Accumulate
and removfr
-solids
4.1^.4-
1
\
IWTGAS
Estimate
water
temperatures
and di solved
ga« cones
4.1 Cl\*
>
vo
00
IQUAL
1PTOT
5imulate quality
constituents
u-sing simple
relationships urHt
wlids and /or
water v/ield
4.1 (zVfc
4.Z<
»>>
IfeAJROT
Place point-
vialued output
in INP&D
4.1(2.^.1
42.
Ct).7>
IPRINT
Place
bar-valued
Output
In IUP/XD
4.2
\
.(2)8
Produce
printed
output
>
4.2 (
».^>
Structure chart 4.1 (tt Impervious land-5»gment application module
-------
INNATER
RETN
Simulate
moisture
retention
4.2/2Y3.I
Simulate
vwater budget
•for impervious
land segment
4.2 C2Y3
t)et ermine
how much
of the mois-
ture supply
runs off
ENRfeTM
from
retention
storage
•structure chart 4.2 fl\3 The
of the IWWLNO application module
99
-------
SOLIDS
Accumulate
and
solids
subroutines
50SL&2
Wash, off
cohd?
method I
ACCUfA
using
method 2.
Acaimulate ^
remoMc solids
independentlv
of runoff
4.1 (1Y4.3
chart 4.1f2).4 The SOLIDS section of
application module
100
-------
1QUM
Simulate
vwashof-f of
qualitw
relationship*
u>i-Vh a ter
-4.im.fe
\MASHSD
WASUDF
Simulate \xy
association
Accumulate
tenwe by a
coirstant uoti
rate and by
a/trlaad
4.2(zU.z
chart 4.lf2^.t The
of -the \V/\PUKU^ abortion module
101
-------
tPTOT
AIRTPT
r
time series in
INPAD
Section
ATEMP
4:2(0- ai
SOUD IP
1
SNOWPT
r
~1
Section
SNO
w
4.2 (O.I5.X
iwa*iP
Section
50U
£>5
4-2C2T7.2
^ecti
l\NT<
IWATIP
^ec-Hon
IW/KTP^-
4.U2T7,
1QAU.1P
on
&M>
42(Z).7.3
*fcticn
IQUAL
*.B17.4
5truc-iure chart 4/.2C2>.T Subroutine qrajp
102
-------
JtolE>
"5ection
50LI05
4-.2C2").8.2
i
place bar-
valued time
series In
WPAD
SNOK/Pfr
r
1>ec-
^MC
__±l
W<&Sl&
- -— — —
Hon
)\W
CIV14-.2
IVNT6A5
!
IWATIB
WATER
4-.ZC2Y8. 1
QAL-IB
Section
\QUAL
•5truc-K«e chart 4-.
qroup
103
-------
IMPA.CC
accumulate
flukes -for
printout
IPRtNT
perform
materials
balances f
produce
printout
_L
IMPPRT
perform mat.
balances 4
write to
printer
IMPRST I
reset flu*
$ state
variables
Structure chart 4.1 (2).q Subroutine qroup I PRINT
104
-------
IMPACC
o
Oi
SMOfeCC | IVNAKCC j
1" ~1
(see
subroutine
qroup
PPR1NT)
I 4.zoYis.i.i!
accumulate
fluv.es for
section
IWKTER
4.^c^•).q.l.^
accumulate
fluy.es for
printout
i
5V.O/V.CC
accumulate
f lu*es ibr
section
SOLIDS
4.Z(2>.
,,.3
IWQACC
accumulate
fluxes for
section
IWTQAS
4.2(2).
,.1.4-
IQA.CC 1
accumulate
fluxes for
section
4.2CZ).
-------
WPPRT
perform
materials
balances
write to
printer
AIRPRT 1 SNOPRT | JWAPRT
i i i — * — i
| for module j | for module j
section section
| ATEMP i j SNOW |
1 II 1
1 II 1
[ 4.zcn.i^.z.i J 1 4.tcn.i*t.2 1
for module
section
IWA.TER
^ *9 /^>\ Q *
^r, fc \~ /«T« t
2.3
SLDPRT
for module
section
SOLIDS
4.1C2V,.
2.4
nwePRT
for module
section
IWT6AS
««v,
2.6
IQPRT 1
for module
section
IQUAL
4.ZC2>.,.2.fe
Structure chart 4.Z(Z).°,.Z Subroutine qroup IMPPRT
-------
o
00
T\w.nis\_.a
rWfcrm
computations
for a reach or
mhced reservoir
4.2(3)
1
HVDR
ADCALC
Simulate
hydraulic
behavior
4.2 (3). 1
CONS
Prepare to
Simulate
advection of
entrained
constituents
4.1(3). 1
4.1 C3Y1 \
FSTORD
RaUAL
Simulate
behavior of
constituents
undergoing
first order
decay
4.2(5). fc
i
ADVECT
HTRCU
Simulate
behavior of
conservative
constituents
4.2 (n%
ADVECT
Simulate advection
of Constituent
totally entrained
in water
RPTOT
Simulate
behavior of
constituents
involved in
biochemical
traftHbrmations
i
i
i
4.21^.7 y
SEt>
Simulate heat
exchange, and
water temp.
4.2ft).4
Simulate
behavior of
inorganic
sediment
4.1 BY 4 \
R6AROT
Put current
values of point
valued time
series in
INPAD
4.2. BV8
42X3)
4.2(?
« >
RPR\MT
rlrf current
values of
bar- valued -Hme
scries in
INPAD
4.2 ft). 9
.9 )>
4-2(3
Produce
printed
outfot
4.1 ft). 10
* y
4.ZC3).IO \
Structure chart 4.1 W Reach/mixed reservoir
application module.
-------
HYDR
Simulate hydraulic
behavior
ROUTE
NOROUT
Calculate outflows
usina, hydraulic
routinq
Calculate outflows
without usinq
hydraulic routing
4.* (»V M
SOLVE
Find outflow demands
corr*sp.tDQ specified
row in RCHTAB
4.1 (9). 1.1.Z.
Solve routing
eqns. used in
case I
4.1 (*).). I.*
EM AND I
I _________ )
Compute values
of auxiliary state
variables
FNDROVN
Swrch RCHTA& for
row* which*bracket"
a qiven VOL
' 4.2 (3). I.1.1
Structure chart 4.Z.tt).l The H^DR section of
ihe RCHRES application module
-------
HTECH
Simulate
h&at
and water
temp
4-.lM.4-
JWECT_ |__.
ADWECT
r
i i
Correct air
temp for
elevation
difference
structure chart 4.2 ttV4 The UTRCH section
of the RCH^tc, application module
110
-------
SEP
behm/ior of
inorganic
sediment
WA&HLD
Simulate
behavior of
wa-shload
4.1«VJ.l
SA.UDID
Simulate
teha\/ior of
eand/grxwel
ADVECT
i
I
1
SINK
~~i
i
i
i
!
_J
Calcula-h?
quantity of
material
settling out of
control volume
4.1 OVV. 1. 1
chart 4.Z 6)5 The SCO section
of the RCMRK afpliartion module
ill
-------
EQUAL
Simulate
of
involved in
biochemical
-transformations
OVRX
Simulate
primary DO,
EOD balances
42 6U!
4/lftYl.l \
KiUTR*
Determine
primary
inorganic, nitrogen
and phosphorous
balances
4.1 (3V7.2.
i«
• )
PLANK
Simulate
behavior of
plankton poplns.
and associated
reactions
4.1
4*m
\ v *»
7-3
.T3
•>
PHCARB
Simulate pH
and carbon
species
4.1 ('
iVl.4
4.2(5V7.4 \
chart 4.2 (^.1 The RQUA.L
of the RCHfcES application (nodule
-------
OXRX
Simulate
primary DO,
BOD balances
ADX/ECT -DlNK O*ftEN
r i r -1- i
L - U J
OKRCA
Simulate
bcnVhal
otojde
:n
demand and
benthal
of BOD
relcace
rtn.i.i
Calculate
o*ygett
reaeration
4.2ft\l.l.2
60DOEC
•
Calculate
ROD deca>/
4.lft\T».l
5Vructure chart 4.1 (^.1.1 OXRX
subroutine group. bO, BOD simulation
-------
deVermvne
primary
inorganic.
nitrogen and
priosphorooc
balance?
4.1 (V).-i.i
AJDN/EOT
f
B6NTH
Simulate
benthal release
of constttuents
4AWVT
n
MITRIF
Simulate
Nitrification
process
4.1 1^. 1.1
.1
DEHIT
Simulate
denitrification
process
^V.l W.T2
%
I
ieCBAL
Pterfofm
material?
balance for
trans-forma-Hoo
•from organic -to
inorganic material
4.1 ftVj.1.4-
Structure chart 4-.l(?VT2 KlUTRX suba)irtine grow
vnor^aniL Kl and P simula-Hon
-------
PLANK
simulate
plankton
populations
and associated
reactions
4.ZC3X7.3
AOVPLK
ADVECT
advecr
plankton
4.2 C 3X7. 3.1
1
SINK
1
1
1
UTRCH
1
1
1
1
calculate
Jiqht-related
information
for alqal
simulation
4.2C3X7.3.2
j 4.2(3X3.1 |
I j
PHYRX
simulate
phytopiankton
4.2(3X7.3.3
4.tO).7.
>
> ,
j 4.2C3X7
i
ZORX
simulate
phytopiankton
4.2(3X7.3.4
u
i
BALRY.
simulate
benthlc
olqae
uan*.*y
Structure chart 4.2(3).7.3 5»mulaie plankton populations and
associated reactions
-------
AiL&RO
G&OCUVC
Calculate unit
gcouAYt anA
respiration
rates for
algae
4r.l(V>.7. 3.3.1
Gftedc. utaetaer
Simulate
tAiybpUunH^D
n
4.2 ft). 7.3.?
1
~T
| WVDTH
compote! growth
rate demand
•too much of any
material
-4.7 ft") n. 3.3.1
1
1
|
1
1
1
Calculate
fhvtoplflnlctor
death
A.a(O.7.?.?.
r
| OR06AU
1
^
i
|
I
I
i
ViUTfcUP
VevHbnm ma^rials
balance -for -Vans-
•ftjrmafion -firom
liVinA 'b dead
organic matarial
4..1(3).7.V3.4
f
?ferforrT>ma^cak
tolance-fcrirdn?-
•forma-tton -ftom
inorganic -to
organic material
4-.2(V)."r.3.3.?
I CJIst/L/V >-\
I aRouo B i
chflri-4.2^.73.?. Simulate
-------
[„_.
6roop k in structure
chart 4.1(^.7.5.3
Simulate berrthic
algae
BAIDTH
Calculate bent hie
alqae deaih
r
1.
6rouo d in
I structure chart
| 4.2(?). 7-7.3
i
Sfrocture diari 4.20^.7.45
algae
-------
PHCARB
4.2(3)7.4
ADVJCT
BEWTH __
4.Z(3)3.l
PHCALC
Calculate
PH
4.Z(3).7.4.l
Structure chart 4.Z(3).7.4 Simulate pH and carbon species
118
-------
N>
O
HXDR
RSARPT
place bar-
valued time
56,rie5 in
cows
4.2(^.2
4:2(3-^.3
FSTRC>
•Section
P5TORD
UTRE.
4ZCO-1.7
PHB.&
4:2 C^-
chart- 4.2 &).°( -Subrouiine
-------
RCHA£C
accumulate
flukes for
printout
RPRtNT
produce
printed
output
4/z.w.io
I
RCHPRT I
perform
materials
balances t
write to
printer
4AC3UO.2.
*(3>>Qy 4to?IQV
RCHRSTJ
resert flu*
e^tatie
variables
4.2(3). 10.%
.10.3\
Structure chart-4.Z(3).IO The RPR I NT section
of the RCHRES application module
121
-------
ro
RCH^CC
fluxes for
printout
HYOACC
for module
section
HYDR
4.ZCf).l O.I.I
CONA.CC
for module
section
CONS
4ACW.
10.1*
HT/XCC
for module
section
HTRCH
4.Z.C*)
.10.1.*
SEOACC
for module
station
3ED
«W
.10.1.*
FSTACC |
for modul<
section
FSTORD
4,.m
>.IO.I.S
OXA.CC
NUTACC
PLKACC
PHA.CC
•for subroutine
group
OXRX
4.t(*).lO.I.6
for subroutine
group
NUTR^
«f».IO.I.T
for subroutine
qroup
PLAN^
4.1 0). 10- 1.8
for subroutine
group
PHCflRB
««mau*
Structure chart4.Z(3).lO.\ Subroutine qroup RCHA.CC
-------
RCHPRT
HYDPR.T |
•for module
section
VNDR
4«Cf>.IO«.l
CONPRT
•for module
section
CONS
4.ZW.JO.Z/Z.
perform
materials
balances and/or
printer
4.2(3). 10.2
i
HTPRT |
for module
section
HTRCH
«CWftt..
SeOPRT |
for module
section
SED
4.2 C37. 10-2.4
FSTPRT |
for module
section
FSTORD
*«»*»
ho
u>
OXPRT I
•for subroutine
group
OXRX
4:2 W. 10.2 .6
NUTPRT I
for subroutine
group
NUTRX
4.2C?). 10.2.7
PLKPRT
•for subroutine
qroup
PLANK
PHPRT I
Ibr subroutine
group
PHCAR&
4.1 «). 10. 21
Structure chert 4.1 W.IO.^ Subroutine qroup RCHPRT
-------
HYPRST I
•for module
section
HY&R
ho
RCHRST
and state
variables
X
CONRST
HTRST
for module
section
CONS
for module
section
HTRCH
SEOR&T
PSTRST
for module
section
SED
for module
section
FSTORD
OV.RST
NUTR6T
•for subroutine
qroup
4-2 (». 10.3.
•for subroutine
9foup
NUTRX
PUK.R6T I
for subroutine
group
PL/XNK.
for subroutine
group
PHCflRB
Structure chart 4.2 (5). 103 Subroutine qroup RCHRSJ
-------
PISPLY
Display a
time, series
AG®
to a larger
(service
TRANS
5HDISP
Transform
data to
Produce
Short - spar\
Produce
annual
4. 2(^3). 3
Struct
rucure c
K
-------
DURANL
Pe-rfornrN
diuraVio
OURDMP
Update
accumulator
4.204V
4.2(14^
DURPUT
Final
processing
4.2-642
DUROUT
set of
-hables
Structure chart 4.2(\4) TKe Duration Analysis u
mo
dule
126
-------
PUTT65
-to
Some service rotrfinc* u&«d bu
~
4.?> T-SPUT
127
-------
•for REPLACE.
to
VO
m/tfd /
-1.5.1 A.I.
T^E*
state currant
-fr^rwe irv
TiSS
-4.3. J. 3
-for
. M
UPPLAb
XVWT
vH-
43.1.3-4
"BUf I NT_ _[
L-Utl'LJ
NEWT5&
new
-------
Functional Description
PART E FUNCTIONAL DESCRIPTION
CONTENTS
General Comments
Descriptions:
1.0 MAIN Program 134
2.0 Manage the Time Series Store (Module TSSMGR) 134
3.0 Interpret a Run Data Set in the User's Control Input
(Module INTERP) 135
4.0 Supervise and Perform Operations (Module OSUPER) 137
4.03 Perform Special Actions (Subroutine SPECL) 138
4.1 Get Required Input Time Series (Module TSGET) 139
4.2 Perform an Operation 139
4.2(1) Simulate a Pervious Land Segment
(Module PERLND) 142
4.2(1).1 Correct Air Temperature for Elevation
Difference (Section ATEMP) 142
4.2(1).2 Simulate Accumulation and Melting of
Snow and Ice (Section SNOW) 143
4.2(1).3 Simulate Water Budget for a Pervious
Land Segment (Section PWATER) 157
4.2(1).4 Simulate Sediment Production and
Removal (Section SEDMNT) 176
4.2(1).5 Estimate Soil Temperatures
(Section PSTEMP) 182
4.2(1).6 Estimate Water Temperature and Dissolved
Gas Concentrations (Section PWTGAS) . . . 133
4.2(1).? Simulate Quality Constituents Using
Simple Relationships with Sediment and
Water Yield (Section PQUAL) 184
4.2(1).8 Estimate Moisture Content of Soil Layers
and Fractional Fluxes of Solutes
(Section MSTLAY) 192
4.2(1).9 Simulate Pesticide Behavior in Detail
(Section PEST) 195
4.2(1).10 Simulate Nitrogen Behavior in Detail
(Section NITR) 207
4.2(1).11 Simulate Phosphorous Behavior in Detail
(Section PROS) 203
4.2(1).12 Simulate Movement of a Tracer
(Section TRACER) 205
130
-------
Functional Description
4.2(2) Simulate an Impervious Land Segment
(Module IMPLND) 207
4.2(2).3 Simulate the Water Budget for an Imper-
vious Land Segment (Section IWATER) ... 207
4.2(2).4 Simulate Accumulation and Removal of
Solids (Section SOLIDS) 210
4.2(2).5 Estimate Water Temperature and Dissolved
Gas Concentrations (Section IWTGAS) . . . 214
4.2(2).6 Simulate Washoff of Quality Constituents
Using Simple Relationships with Solids
and Water Yield (Section IQUAL) 215
4.2(3) Simulate a Free-flowing Reach or Mixed Reservoir
(Module RCHRES) 219
4.2(3).1 Simulate Hydraulic Behavior
(Section HYDR) 219
4.2(3).2 Prepare to Simulate Advection of Fully
Entrained Constituents (Section ADCALC) . 238
4.2(3).3 Simulate Conservative Constituents
(Section CONS) 240
4.2(3).4 Simulate Heat Exchange and Water
Temperature (Section HTRCH) 244
4.2(3).5 Simulate Inorganic Sediment (Section SED) 249
4.2(3).6 Simulate First Order Decay Constituents
(Section FSTORD) 254
4.2(3).? Simulate Constituents Involved in
Biochemical Transformations 256
4.2(3).7.1 Simulate Primary DO and BOD
Balances (Section OXRX) . . . 258
4.2(3).7.2 Simulate Primary Inorganic
Nitrogen and Phosphorous
Balances (Section NUTRX) . . 266
4.2(3).7.3 Simulate Plankton Populations
and Associated Reactions
(Section PLANK) 273
4.2(3).7.4 Simulate pH, Carbon Dioxide,
Total Inorganic Carbon, and
Alkalinity (Section PHCARB) . 296
4.2(11) Copy Time Series (Utility Module COPY) 302
4.2(12) Prepare Time Series for Display on a Plotter
(Module PLTGEN) 302
4.2(13) Display Time Series* in a Convenient Tabular
Format (Utility Module DISPLY) 304
4.2(14) Perform Duration Analysis on a Time Series
(Utility Module DURANL) 309
4.2(15) Generate a Time Series from One or Two Other
Time Series (Utility Module GENER) 316
4.3 TSPUT Module 317
References 317
131
-------
Functional Description
FIGURES
Number Page
3.0-1 Functions and data transfers involved in the operations
portion of HSPF 136
4.2(1).2-1 Snow accumulation and melt processes 144
4.2(1).2-2 Flow diagram of water movement, storages and phase
changes modeled in the SNOW section of the PERLND and
IMPLND Application Modules 146
4.2(1).3-1 Hydrologic cycle 158
4.2(1).3-2 Flow diagram of water movement and storages modeled in
the PWATER section of the PERLND Application Module . . . 160
4.2(1).3-3 Determination of infiltration and interflow inflow ... 164
4.2(1).3-4 Fraction of the potential direct runoff retained by the
upper zone (FRAC) as a function of the upper zone soil
moisture ratio (UZRAT) 167
4.2(1).3-5 Fraction of infiltration plus percolation entering
lower zone storage 172
4:2(1).3-6 Potential and actual evapotranspiration from the lower
zone 175
4.2(1).4-1 Erosion processes 177
4.2(1).4-2 Flow diagram for SEDMNT section of PERLND Application
Module 178
4.2(1).7-1 Flow diagram for PQUAL section of PERLND Application
Module 186
4.2(1).8-1 Flow diagram of the transport of moisture and solutes,
as estimated in the MSTLAY section of the PERLND
Application Module ; . . . 193
4.2(1).9-1 Flow diagram showing modeled movement of chemicals in
solution 197
4.2(1).9-2 Freundlich isotherm calculations 200
4.2(1).10-1 Flow diagram for nitrogen reactions 204
4.2(1).11-1 Flow diagram for phosphorus reactions 206
4.2(2)-1 Impervious land segment processes 208
4.2(2).3-1 Hydrologic processes 209
4.2(2).4-1 Flow diagram of the SOLIDS section of the IMPLND
Application Module 212
4.2(2).6-1 Flow diagram for IQUAL section of IMPLND Application
Module 216
4.2(3)-1 Flow of materials through a RCHRES 220
4.2(3)-1-1 Flow diagram for the HYDR Section of the RCHRES
Application Module 222
4.2(3).1-2 Graphical representation of the equations used to compute
outflow rates and volume 225
4.2(3).1-3 Typical RCHRES configurations and the method used to
represent geometric and hydraulic properties 226
4.2(3).1-4 Graphical representation of the work performed by
subroutine ROUTE 229
132
-------
Functional Description
4.2(3).1-5 Graphical representation of the work performed by
subroutine NOROUT 234
4.2(3).1-6 Illustration of quantities involved in calculation of
depth 236
4.2(3).2-1 Determination of weighting factors for advection
calculations 240
4.2(3).3-1 Flow diagram for conservative constituents in the CONS
section of the RCHRES Application Module 240
4.2(3).4-1 Flow diagram for HTRCH section of the RCHRES Application
Module 245
4.2(3).5-1 Flow diagram for washload in the SED section of the
RCHRES Application Module 25°
4.2(3).5-2 Flow diagram for sandload in the SED section of the
RCHRES Application Module 251
4.2(3).6-1 Flow diagram for first order decay constituents in the
FSTORD section of the RCHRES Application Module 255
4.2(3)-7.1-1 Flow diagram for dissolved oxygen in the OXRX subroutine
group of the RCHRES Application Module 259
4.2(3).7.1-2 Flow diagram for biochemical oxygen demand in the OXRX
subroutine group of the RCHRES Application Module .... 260
4.2(3).7.2-1 Flow diagram for inorganic nitrogen in the NUTRX
subroutine group of the RCHRES Application Module .... 267
4.2(3).7.2-2 Flow diagram for ortho-phosphate in the NUTRX group of
the RCHRES Application Module 268
4.2(3).7.3-1 Flow diagram for phytoplankton in the PLANK section of
the RCHRES Application Module 274
4.2(3).7.3-2 Flow diagram for dead refractory organics in the PLANK
section of the RCHRES Application Module 274
4.2(3).7.3-3 Flow diagram for zooplankton in the PLANK section of the
RCHRES Application Module 275
4.2(3).7.3-4 Flow diagram for benthic algae in the PLANK section of
the RCHRES Application Module 275
4.2(3).7.3-5 Relationship of parameters for special advection of
plankton 277
4.2(3).7.3-6 Zooplankton assimilation efficiency 259
4.2(3).7-4-1 Flow diagram of inorganic carbon,In the PHCARB section
of the RCHRES Application Module 297
4.2(13)-1 Sample Short-Span Display (First Type) 306
4.2(13)-2 Sample Short-Span Display (Second Type) 307
4.2(13)-3 Sample Long-Span (Annual) Display 308
4.2(14)-1 Definition of Terms Used in Duration Analysis Module. . . 310
4.2(14)-2 Sample Duration Analysis Printout 312
133
-------
Functional Description
GENERAL COMMENTS
For a discussion on how this part of the documentation is organized, refer
to Section 5 in Part C "Standards and Conventions.",
The subprograms are discussed in numerical order, following the numbering
system used on the structure charts (Part D). When studying the
descriptions which follow, you will find it helpful to refer to the
appropriate structure chart and to the pseudo code in the Programmer's
Supplement (Johanson, et al. 1979).
1.0 MAIN Program
The MAIN program stands at the head of the structure (Structure Chart 1.0)
and calls, directly or indirectly, all the other modules in the system. The
functions performed are:
(1) Preprocess the Users Control Input (UCI). Subroutine USRRDR transfers
the UCI to a direct access file, appends a value at the end of each
record which points to the next non-comment record, and recognizes data
set headings and delimiters: RUN, END RUN, TSSM, END TSSM.
(2) Call TSSMGR if a TSSM data set has been found.
(3) If a RUN data set has been found, call INTERP to interpret it and then
call OSUPER to supervise execution of it.
2.0 Manage the Time Series Store (Module TSSMGR)
General Description of Module TSSMGR
This module maintains a user's Time Series Store (TSS) and performs some
housekeeping chores associated with the data sets in it. From the point of
view of the computer's operating system, the TSS is a single file (which may
be very large). However, the HSPF software can store many distinct time
series in this file. This permits a user easily to keep track of the
various time series with which he is dealing. Furthermore, he need only
refer to a single disc file for all his time series input and output needs,
no matter how many time series are involved. This simplifies communication
with the computer system.
Time series are arranged within the TSS in one or more "data sets". The
contents of each data set and its physical location in the TSS are recorded
in a directory located at the start of the TSS. The primary function of
module TSSMGR is to maintain this directory, from input supplied by the
user. He can add a new data set label to the directory, update certain parts
of the label (such as the SECURITY option), scratch a data set label (and,
thus, the data set contents), extend the space allocated to a data set, or
134
-------
Module TSSMGR
show the contents of one or all of the labels in the directory. The
commands used to achieve this are documented in Part F, Section 2.
The design of the TSS is based on our experience with HSPX and HSPII.
Extensions to the HSPX time series management system include:
(1) The storage of data in compressed form. Disk space is saved by
improving the method of encoding values which occur in a sequence, such
as a string of zeros.
(2) The inclusion of multimembered data sets in the TSS. This permits
storage of several time series in the same data set, which is useful if
one needs to save several related functions of time. For example, if a
user needs to save the discharge of water and associated constituents
from a stream reach, he needs to create and refer to only one data set
(with several members). Up to 20 time series may be stored in a single
data set.
Further details of the organization of a TSS are given in the documentation
for program NEWTSS (Appendix III).
Once the label for a dataset has been created and space reserved for it in
the TSS, time series data can be stored in the data set. This is done by an
operating module (Section 4.0); either a utility module (eg. COPY) or an
application module (eg. PERLND).
3.0 Interpret a RUN Data Set in the User's Control Input (Module INTERP)
General Description of Module INTERP
This module, known as the Run Interpreter, translates a RUN data set in the
User's Control Input (documented in Section 4 of Part F) into many
elementary instructions, for later use by other parts of the system, when
the time series are operated on. To do this, the Run Interpreter performs
such tasks as:
(1) Check and augment the data supplied by the user.
(2) Decide which time series will be required and produced by each
operation, based on the user's data and built-in tables which contain
information on the various operations.
(3) Allocate INPAD rows to the various time series.
(4) Read the control data, parameters, and initial conditions supplied for
each opertion, convert them to internal units, and supply default
values where required.
The output of the Run Interpreter is stored in disk-based files containing
instructions to be read by the Operations Supervisor, TSGET and TSPUT
(Figure 3.0-1). The instruction files are:
135
-------
Run Interpreter
(1) The Operations Supervisor Instruction File (OSUPFL). This file contains
instructions which the Operations Supervisor reads to manage the
operations in a run. This includes information on:
(a) the configuration of the scratch pads (time intervals and widths)
(b) the configuration of the EXGROUPs and INGROUPs (number of
EXGROUPs, number INGROUPs in each EXGROUP, operations in each
INGROUP, etc.) (EXGROUPs have not yet been implemented)
(2) The Operation Status Vector File (OSVFL). The operations in a run are
interrupted every time an INPAD span is completed (Part B, Section
3.2) To save machine core space, the system is designed so that the
various operations all use the same area of core. This requires that
upon interruption, all information necessary to restart an operation be
stored in a disk file. The data, called the "Operation Status Vector"
(OSV), reside in a string of contiguous locations in core and have a
structure specified in the Programmer's Supplement (Johanson, et al.
1979). The disc file OSVFL contains an exact copy of the OSV for each
operation. It is used to restore the OSV in core when the operation is
resumed after interruption.
(3) The Input Time Series Instruction File (TSGETF) and the Output Time
Series Instruction File (TSPUTF). These files contain instructions
which govern the transfer of pieces of time series into and out of the
INPAD, respectively. Each instruction enables module TSGET to retrieve
a specified piece of time series from one of the source volumes (Figure
3.0-1), transform it to the interval and form required for the INPAD,
and insert it in the desired row of the INPAD. In the case of TSPUTF,
the sequence is the reverse of that just described.
Each operation has its own set of instructions in TSGETF and TSPUTF
which are read whenever modules TSGET and TSPUT are called upon to
service that operation (every INSPAN).
CO The Special Action Instruction File (SPACFL). Each record of this file
contains a single special action instruction, which specifies the
action required to be taken in a given operation at a specific time,
e.g. report operation state, modify a state variable, (not yet
implemented)
The structures of these files are documented in the Programmer's Supplement.
4.0 Supervise and Perform Operations (module OSUPER)
Function of Operations Group
The Operations group of modules handles all the manipulations of time series
and thus, performs most of the work in a run. Subroutine OSUPER controls
the group. It performs some of the tasks itself, but it invokes subordinate
modules to do other tasks.
137
-------
Operations Group
General Description of Subroutine OSUPER
The primary tasks of subroutine OSUPER are to ensure that the various
operations in the run are called in the correct sequence and that the
associated time series and OSVs are input and/or output at the
required junctures (see Part B, Section 3.2). OSUPER uses a nest of
DO-loops to control the sequencing. The instruction file OSUPFL
specifies the ranges of the loops and supplies information ("keys")
which enable OSUPER, TSGET and TSPUT to correctly access the other
instruction files. OSUPER reads an instruction from disc each time an
operation starts a new INSPAN. Using this information, it then:
(1) calls TSGET, to supply the required input time series (using
TSGKST, TSGKND)
(2) reads the OSV from disc (using keys OSVKST, OSVKND)
(3) calls the operating module (OMCODE indicates which one is to be
called)
When the INSPAN is over, OSUPER:
(1) writes the OSV to disc (using keys OSVKST, OSVKND)
(2) calls TSPUT, to output time series (using keys TSPKST, TSPKND)
4,03 Perform Special Actions (Subroutine SPECL)
HSPF permits the user to perform certain "Special Actions" during the
course of a run. A special action instruction specifies the following:
1. The operation on which the action is to be performed (eg. PERLND
10)
2. The date/time at which the action is to be taken.
3. The type and location, within COMMON block SCRTCH, of the data
item to be updated.
U. The action to be performed. Two choices are available:
a) Reset the variable to a specified value
b) Increment the variable by a specified value
The special action facility is used to accomodate things such as:
138
-------
Operations Group
1. Human intervention in a watershed. Events such as plowing,
cultivation, fertilizer and pesticide application, and harvesting
are simulated in this way.
2. Changes to parameters. For example, a user may wish to alter the
value of a parameter for which 12 monthly values cannot be
supplied. He can do this by specifying a special action for that
variable. He could reset the parameter to its original value by
specifying another special action, to be taken at a later time.
At present, Special Actions can only be performed on variables in the
PERLND module. The input is documented in Section 4.10 of Part F.
4.1 Get Required Input Time Series (module TSGET)
The task of this module is to insert in the INPAD all input time
series required by an operation. OSUPER calls it (once) each time an
operation is to commence an INSPAN, passing to it the keys of the
first and last records in TSGETF which must be read and acted upon.
Each instruction causes a row of the INPAD to be filled. TSGET can
draw its input time series from any of the following source "volumes":
TSS, sequential file and INPAD (Figure 3.0-1).
TSGET will, if necessary, automatically transform the time interval
and "kind" (Appendix III) of the time series, as it is brought from
the source location to the INPAD (target). TSGET can also perform a
linear transformation on the values in a time series; for example, if
the source contains temperatures in degrees C and the INPAD needs them
in degrees F.
4.2 Perform an Operation
Function of an Operating Module
An operating module (OM) is at the center of every operation (Part B,
Section 2.1). When the Operations Supervisor calls an OM the time
series which it requires are already in the INPAD. The task of the OM
is to operate on these input time series. The results of this work
are:
(1) updated state variables. The OM constantly updates any state
variables. These are located in the OSV. Thus, when the OM
returns control to the Operations Supervisor, which copies the
OSV to disc, the latest values of all state variables are
automatically preserved.
(2) printed output. The OM accumulates values, formats them and
routes these data to the line printer.
139
-------
Operations Group
(3) output time series. The OM places these in the INPAD but is not
concerned with their ultimate disposition; this is handled by
module TSPUT.
Note that all time series simultaneously present in an INPAD have the
same constant time interval. This implies that, internally, all time
series involved in an operation have the same time interval.
Externally, the time series may have differing time intervals. Part
of the function of modules TSGET and TSPUT is to convert time series
from external to internal time intervals and vice versa.
Sub-divisions in an Operating Module
An operating module may be divided into several distinct "sections,"
each of which may be selectively activated in a given run, under the
user's control, e.g. the Pervious Land-segment module (PERLND)
contains twelve sections, the first being air temperature correction,
the last tracer (conservative substance) simulation (Structure Chart
4.2(1)). The operating procedure is as follows: in each time interval
of the INSPAN, the operating module calls each of its active sections
in the order in which they are built into the code (the sequence can
not be altered by the user). When the INSPAN has been covered the
operating module returns control to OSUPER which determines the next
action to be taken. This procedure implies that an operating module
must be arranged so that a section is called after any others from
which it requires information. For example, in the Pervious
Land-segment module (Structure Chart 4.2(1)), the sediment calculation
section may use data computed by the snow and water balance-sections
but not by Sections 5 through 10. This kind of information flow is
called an inter-section data transfer (ISDT).
Partitioning of an Operation
A user may activate one group of module sections in an initial run and
other groups in subsequent runs. Thus, he may "partition" an
operation. For example, he may wish to calibrate the hydraulic
response of a set of river reaches before moving on to simulate the
behavior of constituents contained in the water. If this type of work
involves ISDT's between the sections handled in different runs, it
follows that:
(1) The time series involved in the ISDT's must be stored between
runs, probably in the TSS.
(2) In the second run the system will expect the user to specify external
sources for all these time series.
140
-------
Operations Group
Some users will be confused by the rules for partitioning operations, but
our experience indicates this will be outweighed by the flexibility which it
brings.
Numbering of Operating Modules
In principle, there is no limit to the number of operating modules which the
system can accommodate. Ultimately, we expect a large number of modules
ranging from very simple utility modules (eg. COPY) to very complex
simulation algorithms (eg. PERLND). Although the size and complexity of the
modules vary greatly, they all are, logically, of equal rank (Figure 2-3,
Part B). The adopted numbering system reflects this. Every operating
module is identified by the number 4.2 (Structure Chart 4.0) and is
distinguished from the others only by a subscript. For example, the
Pervious Land-segment Module is 4.2(1) and the Reach/Mixed Reservoir Module
4.2(3).
Inserting Additional Operating Modules
A user may insert additional modules. To do this he must:
(1) Write or adapt his operating module. This includes restructuring the
data into an OSV which conforms with the requirements of the HSPF
system. (This task may be time consuming).
(2) Add a section of code to the Run Interpreter to interpret the UCI for
the new module.
(3) Add data to the information file (INFOFL) and, if necessary, to the
warning and error message files.
(4) Make minor changes to subroutines OPNBLK and OSUPER.
Types of Operating Modules
There are two types of operating modules; utility modules and application
modules. Utility modules perform any operations involving time series which
are essentially auxiliary to application operations, e.g. input time series
data from cards to the TSS using COPY, multiply two time series together to
obtain a third one, plot several time series on the same graph. The utility
modules perform many of the functions which were previously part of HSP
LIBRARY or HSP UTILITY. They are given numbers starting with 4.2(11).
Application (simulation) modules represent processes, or groups of
processes, which occur in the real world. They have been allocated numbers
4.2(1) through 4.2(10) although, at present, only three application modules
are written.
141
-------
Module PERLND
4.2(1) Simulate a Pervious Land Segment (Module PERLND)
A land segment is a subdivision of the simulated watershed. The boundaries
are established according to the user's needs, but generally, a segment is
defined as an area with similar hydrologic characteristics. For modeling
purposes water, sediment, and water quality constituents leaving the
watershed move laterally to a downslope segment or to a reach/reservoir. A
segment of land which has the capacity to allow enough infiltration to
influence the water budget is considered pervious. In HSPF, PERLND is the
module that simulates the water quality and quantity processes which occur
on a pervious land segment.
The primary module sections in PERLND simulate snow accumulation and melt
(Section SNOW), the water budget (section PWATER), sediment produced by land
surface erosion (section SEDMNT), and water quality constituents by various
methods (section PQUAL and the agri-chemical sections). Other sections
perform the auxiliary functions of correcting air temperature (section
ATEMP) for use in snowmelt and soil temperature calculations, producing soil
temperatures (section PSTEMP) for estimating the outflow temperatures and
influencing reaction rates in the agri-chemical sections, and determining
outflow temperatures which influence the solubility of oxygen and carbon
dioxide. Structure Chart 4.2(1) shows these sections and their
relationships to each other and to PPTOT, PBAROT, and PPRINT. The last three
sections manipulate the data produced. Section PPTOT places state varables
(point values) and PBAROT places flux variables which are actually averages
over the interval (bar values) into the INPAD. PPRINT produces the printable
results in the quantity and frequency that the use specifies. The sections
in the structure chart are executed from left to right.
4.2(1).1 Correct Air Temperature for Elevation Difference (Section ATEMP of
Modules PERLND and IMPLND)
Purpose
The purpose of ATEMP is to modify the input air temperature to represent the
mean air temperature over the land segment. This module section is part of
both PERLND or IMPLND. Air temperature correction is needed when the
elevation of the land segment is significantly different than the elevation
at the temperature gage. If no correction for elevation is needed, this
module section can be skipped.
Method
The lapse rate for air temperature is dependent upon precipitation during
the time interval. If precipitation occurs, a wet lapse rate of 0.0035
degrees F per foot difference in elevation is assumed. Otherwise, a dry
lapse rate, that varies with the time of day, is used. A table of 24 hourly
142
-------
Module Section ATEMP
dry lapse rates varying between 0.0035 to 0.005 is built into the system.
The corrected air temperature is:
AIRTMP = GATMP - LAPS*ELDAT (1)
where:
AIRTMP = corrected air temperature in degrees F
GATMP = air temperature at gage in degrees F
LAPS = lapse rate in degrees F/ft
ELDAT = elevation difference between the land segment and the
gage in ft
4.2(1).2 Simulate Accumulation and Melting of Snow and Ice (section SNOW of
modules PERLND and IMPLND)
Purpose
SNOW deals with the runoff derived from the fall, accumulation and melt of
snow. This is a necessary part of any complete hydrologic package since
much of the runoff, especially in the northern half of the United States, is
derived from snow conditions.
Approach
Figure 4.2(1).2-1 illustrates the processes involved in snow accumulation
and melt on a land segment. The algorithms used are based on the work by
the Corps of Engineers (1956), Anderson and Crawford (1964), and Anderson
(1968). Empirical relationships are employed when physical ones are not
well known. The snow algorithms use meteorologic data to determine whether
precipitation is rain or snow, to simulate an energy balance for the
snowpack, and to determine the effect of the heat fluxes on the snowpack.
Five meteorologic time series are required by SNOW for each land segment
simulated. They are:
precipitation
air temperature
solar radiation
dewpoint
wind velocity
A value from each of these time series is input to SNOW at the start of each
simulation interval. However, some of the meteorological time series are
only used intermittently for calculating rates, such as in the calculation
of the potential rate of evaporation from the snowpack.
Air temperature is used to determine when snow is falling. Once snow begins
to accumulate on the ground, the snowpack accumulation and melt calculations
take place. Five sources of heat which influence the melting of the
snowpack are simulated:
143
-------
+ + + +
+ + +
y -• .//>'.
R AIN/SN6W /' ','•
DETERMINATION' AIR
+ ,*+ / '•-' /TEMPERATURE
,'.. / DEWPOINT
/' .<'•/ TEMPERATURE
WIND
/ ' /.!
' \ ,
COWOINS^YiON \ \
HEA.T EXCHANGE
GROUNDMELT
LIQUID STORAGE
IN SNOWPACK
LAND SURFACE
AREAL EXTENT
OF SNOW COVER
Figure 4.2(1).2-1 Snow accumulation and melt processes
-------
Module Section SNOW
1. net radiation heat (RADHT), both longwave and shortwave
2. convection of sensible heat from the air (CONVHT)
3. latent heat transfer by condensation of moist air on the
snowpack (CONDHT)
4. heat from rain, sensible heat from rain falling
(RNSHT) and latent heat from rain freezing on the snowpack
5. conduction of heat from the underlying ground to the snowpack
(GMELTR)
Other heat exchange processes such as latent heat from evaporation are
considered less significant and are not simulated.
The energy calculations for RADHT, CONVHT, and CONDHT are performed by
subroutine HEXCHR while GMELTR is calculated in subroutine GMELT. Latent
heat from rain freezing is considered in subroutine WARMUP. RNSHT is
computed in the parent subroutine SNOW. For uniformity and accounting,
energy values are calculated in terms of the water equivalent which they
could melt. It takes 202.4 calories per square cm on the surface to melt
one inch water equivalent of snow at 32 degrees F. All the sources of heat
including RNSHT are considered to be positive (incoming to the pack) or
zero, except RADHT which can also be negative (leaving the pack).
Net incoming heat from the atmosphere (the sum of RADHT, CONVHT, CONDHT, and
RNSHT) is used to warm the snowpack. The snowpack can be further warmed by
the latent heat released upon rain freezing. Any excess heat above that
required to warm the snowpack to 32 degrees F is used to melt the pack.
Likewise, net loss of heat is used to cool the snowpack producing a negative
heat storage. Futhermore, incoming heat from the ground melts the snowpack
from the bottom independent of the atmospheric heat sources except that the
rate depends on the temperature of the snowpack.
Figure 4.2(1).2.2 gives a schematic view of the moisture related processes
modeled in section SNOW. Precipitation may fall as rain or snow on the
snowpack or the ground. Evaporation only occurs from the frozen portion of
the pack (PACKF). The frozen portion of the pack is composed of snow and
ice. The ice portion of PACKF is considered to be in the lower part of the
snowpack, so it is the first to melt when heat is conducted from the ground.
Similarly, the snow portion of PACKF is the first to melt when atmospheric
heat increases. Melted PACKF and rain falling on the snowpack produce the
water portion of the total snowpack which may overflow the capacity of the
pack. The water yield and rain on the bare ground becomes input to module
section PWATER or IWATER. These moisture related processes as well as the
heat exchange processes are discussed later in more detail.
Heat transfer from incoming rain (RNSHT) to the snowpack is calculated in
the parent subroutine SNOW (Section 4.2(1).2). The following physically
based equation is used:
RNSHT = (AIRTMP - 32.0)*RAINF/W.O (2)
145
-------
/SNOWE\
I evaporation >
I from PACKF/
/ PREC \
V precipitation/
o
SNOWF
snowfall
on
PACKF
rainfall
total snowpack
/ rainfal I
\ enter in
rainfall
on ground
MELT
atmospheric
heat melt m
ftomPMfffi
Csnow m
FREEZE
freeling
of WCKW
to
MELT
atmospheric
heat melt
irom PACKF
(wiovw
GMEUTR
qround
heat melt
from pack
Cice
first
PACKW
liquid
water
qround
heat melt
from pack
first!
frozen portion
WYIELD
water
yield
from
PACKW
Fiqure 4.2(0.2-^ Flow diaqram of water movement, ^toraqas and phase
modeled in -the 9MOW section of the PERLK1D and IMPLK1D Application Modules
146
-------
Module Section SNOW
where:
AIRTMP = temperature of the air in degrees F
RAINF = rainfall in inches
144.0 = factor to convert to equivalent depth of melt
32.0 = freezing point in degrees F
Other characteristics of the snowpack are also determined in the main
subroutine SNOW. The fraction of the land segment covered by the snowpack
is estimated by merely dividing the depth of the snowpack by a cover index
(COVINX) which is a function of the parameter COVIND and the history of the
pack as explained in subroutine EFFPRC. The temperature of the snowpack is
estimated by:
PAKTMP = 32.0 - NEGHTS/(0.00695*PACKF) (3)
where:
PAKTMP = mean temperature of the snowpack in degrees F
NEGHTS = negative heat storage in inches of water equivalent
PACKF = frozen contents of the snowpack in inches of water
equivalent
0.00695 = physically based conversion factor
4.2(1).2.1 Estimate Meteorological Conditions (subroutine METEOR)
Purpose
Subroutine METEOR estimates the effects of certain meteorological conditions
on specific snow-related processes by the use of empirical equations. It
determines whether precipitation is falling as snow or rain. The form of
precipitation is critical to the reliable simulation of runoff and snowmelt.
When snow is falling, the density is calculated in order to estimate the
depth of the new snowpack. The fraction of the sky which is clear is also
estimated for use in the radiation algorithms, and the gage dewpoint is
corrected if it is warmer than air temperature. ;
Method
The following expression is used to calculate hourly the effective air
temperature below which snowfall occurs:
SNOTMP = TSNOW + (AIRTMP - DEWTMP)*(0.12 + 0.008*AIRTMP) (4)
where:
SNOTMP = air temperature below which snowfall occurs in degrees F
TSNOW = parameter in degrees F
AIRTMP = air temperature in degrees F
DEWTMP = dewpoint in degrees F
SNOTMP is allowed to vary in this calculation by a maximum of one degree F
from TSNOW. When AIRTMP is equal to or greater than SNOTMP, precipitation is
assumed to be rain.
147
-------
Module Section SNOW
When snowfall occurs, its density is estimated as a function of air
temperature according to:
RDNSN = RDCSN + (AIRTMP/100.0)**2 (5)
where:
RDNSN = density of new snowfall (at zero degrees F or greater)
relative to liquid water
RDCSN = parameter designating density of new snow at an air temperature
of zero degrees F and lower, relative to liquid water
RDNSN is used in subroutine EFFPRC to calculate the new depth of the
snowpack resulting from the addition of the snow. This and all other snow
density terms are in water equivalent (inches) per depth of the snowpack
(inches).
The fraction of the sky which is clear (SKYCLR) is needed for the
calculation of the longwave back radiation to the snowpack from the clouds
(done in subroutine HEXCHR). SKYCLR is set to the minimum value of 0.15
when precipitation occurs. Otherwise, it is increased each simulation time
interval as follows:
SKYCLR = SKYCLR + (0.0004*DELT) (6)
where:
DELT = simulation time interval in min
SKYCLR increases until either it reaches unity or precipitation causes it to
be reset.
A gage dewpoint higher than air temperature is not physically possible and
will give erroneous results in the calculation of snowpack evaporation.
Therefore, dewpoint is set equal to the air temperature when this situation
occurs. Otherwise, the gage dewpoint is used.
4.2(1).2.2 Determine the Effect of Precipitation on the Pack (subroutine
EFFPRC)
Purpose
The purpose of this subroutine is to add the falling snow to the pack,
determine the amount of rain falling on the snowpack, and adjust the
snowpack dullness to take into account new snow.
Method
The amount of precipitation falling as snow or rain is determined in
subroutine METEOR. Subroutine EFFPRC accounts for the influence that
snowfall and rain have on the land segment. The subroutine begins by
increasing the snowpack depth by the amount of snow falling on the pack
divided by its density.
148
-------
Module Section SNOW
The fraction of the land segment which is a covered by the snowpack (SNOCOV)
is determined by re-evaluating the index to areal coverage (COVINX) When
the frozen contents of the pack (PACKF) exceeds the value of the parameter
describing the maximum PACKF required to insure complete areal coverage by
snow cover (COVIND), then COVINX is set equal to COVIND. Otherwise, COVINX
is equal to the largest previous value of PACKF. SNOCOV is PACKF/COVINX if
PACKF < COVINX. The amount of rain falling on the snowpack is that fraction
of the precipitation which falls as rain multiplied by the SNOCOV. Rain
falling on the snowpack will either freeze, adding to the frozen portion of
the pack and produce heat used to warm the pack (see subroutine WARMUP), or
it will increase the liquid water content of the pack (see subroutine
LIQUID). Any rain not falling on the pack is assumed to land on bare
ground.
When snowfall occurs, the index to the dullness of the snowpack (DULL) is
decreased by one thousand times the snowfall for that interval. However, if
one thousand times the snowfall is greater than the previous value for DULL,
then DU LL is set to zero to account for a new layer of perfectly
reflectable snow. Otherwise, when snowfall does not occur, DULL is increased
by one index unit per hour up to a maximum of 800. Since DULL is an
empirical term used as an index, it has no physical units. DULL is used to
determine the albedo of the snowpack which in turn is used in the shortwave
energy calculations in subroutine HEXCHR.
4.2(1).2.3 Compact the Pack (subroutine COMPAC)
Purpose
The addition of new snow will reduce the density as well as increase the
depth of the snowpack as in subroutine EFFPRC. The pack will tend to
compact with age until a maximum density is reached. The purpose of
subroutine COMPAC is to determine the rate of compaction and calculate the
actual change in the depth due to compaction.
Method
When the relative density is less than 55 percent compaction is assumed to
occur. The rate of compaction is computed according to the empirical
expression:
COMPCT s 1.0 - (0.00002*DELT60*PDEPTH*(0.55 - RDENPF)) (7)
where:
COMPCT = unit rate of compaction of the snowpack per interval
DELT60 = number of hours in an interval
PDEPTH = depth of the snowpack in inches of total snowpack
RDENPF = density of the pack relative to liquid water
The new value for PDEPTH is COMPCT times PDEPTH. PDEPTH is used to
calculate the relative density of the snowpack which affects the liquid
water holding capacity as determined in subroutine LIQUID.
149
-------
Module Section SNOW
4.2(1).2.4 Simulate Evaporation from the Pack (subroutine SNOWEV)
Purpose
The SNOWEV subroutine estimates evaporation from the snowpack (sublimation).
Method
Evaporation from the snowpack will occur only when the vapor pressure of the
air is less than that of the snow surface, that is, only when the air vapor
pressure is less than 6.108 mbar which is the maximum vapor pressure that
the thin surface film of air over the snowpack can attain. When this
condition is met the evaporation is computed by the empirical relationship:
SNOWEP = SNOEVP*0.0002*WINMOV*(SATVAP - VAP)*SNOCOV (8)
where:
SNOWEP = potential rate of evaporation from the frozen part of the
snowpack in inches of water equivalent/interval
SNOEVP = parameter used to adjust the calculation to field conditions
WINMOV = wind movement in miles/interval
SATVAP = saturated vapor pressure of the air at the current air
temperature in mbar
VAP = vapor pressure of the air at the current air temp, in mbar
SNOCOV = fraction of the land segment covered by the snowpack
The potential (SNOWEP) will be fulfilled if there is sufficient snowpack.
Otherwise, only the remaining pack will evaporate. For either case,
evaporation occurs only from the frozen content of the snowpack (PACKF).
4.2(1).2.5 Estimate Heat Exchange Rates (except ground melt and rain heat)
(subroutine HEXCHR)
Purpose
The purpose of this subroutine is to estimate the heat exchange from the
atmosphere due to condensation, convection, and radiation. All heat
exchanges are calculated in terms of equivalent depth of melted or frozen
water.
Method of Determining Heat Supplied by Condensation
Transfer of latent heat of condensation can be important when warm moist air
masses travel over the snowpack. Condensation occurs when the air is moist
enough to condense on the snowpack. That is, when the vapor pressure of the
air is greater than 6.108 mbar. This physical process is the opposite of
snow evaporation; the heat produced by it is calculated by another
empirical relationship:
150
-------
Module Section SNOW
CONDHT = 8.59*(VAP - 6.108)*CCFACT*0.00026*WINMOV (9)
where:
CONDHT = condensation heat flux to the snowpack in inches of water
equivalent/interval
VAP = vapor pressure of the air at the current air temp, in mbar
CCFACT = parameter used to correct melt values to field conditions
WINMOV = wind movement in miles/interval
CONDHT can only be positive or zero, that is, incoming to the pack.
Method of Determining Heat Supplied by Convection
Heat supplied by turbulent exchange with the atmosphere can occur only when
air temperatures are greater than freezing. This convection of heat is
calculated by the empirical expression:
CONVHT r (AIRTMP - 32.0)*(1.0 - 0.3*MELEV/10000.0)* (10)
CCFACT*0.00026*WINMOV
where:
CONVHT = convective heat flux to the snowpack in inches of water
equivalent/interval
AIRTMP = air temperature in degrees F
MELEV = mean elevation of the land segment above sea level in ft
In the simulation, CONVHT also can only be positive or zero, that is, only
incoming.
Method of Determining Heat Supplied by Radiation
Heat supplied by radiation is determined by:
RADHT = (SHORT + LONG)/203.2 (11)
where:
RADHT = radiation heat flux to the snowpack in inches of
water equivalent/interval
SHORT = net solar or shortwave radiation in langleys/interval
LONG = net terrestrial or longwave radiation in langleys/interval
The constant 203.2 is the number of langleys required to produce one inch of
melt from snow at 32 degrees F. RADHT can be either positive or negative,
that is, incoming or outgoing.
SHORT and LONG are calculated as follows. Solar radiation, a required time
series, is modified by the albedo and the effect of shading. The albedo or
reflectivity of the snowpack is a function of the dullness of the pack (see
subroutine EFFPRC for a discussion of DULL) and the season. The equation for
calculating albedo (ALBEDO) for the 6 summer months is:
151
-------
Module Section SNOW
ALBEDO = 0.80 - 0.10*(DULL/24.0)**0.5 (12)
The corresponding equation for the winter months is:
ALBEDO = 0.85 - 0.07*(DULL/24.0)**0.5 (13)
ALBEDO is allowed a minimum value of 0.45 for summer and 0.60 for winter.
The hemispheric location of the land segment is taken into account for
determini ng summer and winter in using the above equation. This is done
through the use of the latitude parameter which is positive for the northern
hemisphere.
Once the albedo of the pack is found then solar radiation (SHORT) is
modified according to the equation:
SHORT = SOLRAD*(1.0 - ALBEDO)*(1.0 - SHADE) (14)
where:
SOLRAD = solar radiation in langleys/interval
SHADE = parameter indicating the fraction of the land segment which
is shaded
Unlike shortwave radiation which is more commonly measured, longwave
radiation (LONG) is estimated from theoretical consideration of the emitting
properties of the snowpack and its environment. The following equations are
based on Stefan's law of black body radiation and are linear approximations
of curves in Plate 5-3, Figure 6 in Snow Hydrology (Corps of Engineers
1956). They vary only by the constants which depend on air temperature. For
air temperatures above freezing:
LONG = SHADE*0.26*RELTMP + (1.0 - SHADE)*(0.2*RELTMP - 6.6) (15)
And for air temperatures at freezing and below:
LONG = SHADE*0.20»RELTMP + (1.0 - SHADE)*(0.17*RELTMP - 6.6) (16)
where:
RELTMP = air temperature minus 32 in degrees F
6.6 = average back radiation lost from the snowpack in open
areas in langleys/hr
Since the constants in these equations were originally based on hourly time
steps, both calculated values are multiplied by DELTOO, the number of hours
per interval, so that they correspond to the simulation interval. In
addition, LONG is multiplied by the fraction of clear sky (SKYCLR) when it
is negative to account for back radiation from clouds.
152
-------
Module Section SNOW
4.2(1).2.6 Simulate Loss of Heat from Pack (subroutine COOLER)
Purpose
The purpose of this code is to cool the snowpack whenever it is warmer than
the ambient air and thus loses heat. This is accomplished by accumulating
negative heat storage which increases the capacity of the pack to later
absorb heat without melting as simulated in subroutine WARMUP.
Method
In every interval where there is heat loss to the atmosphere and the
temperature of the snowpack is greater than the air temperature, the
negative heat storage will increase; that is, the pack will cool. However,
there is a maximum negative heat storage. The maximum negative heat storage
that can exist at any time is found by assuming a linear temperature
distribution from the air temperature which is considered to be above the
pack to 32 degrees at the bottom of the snowpack. This maximum negative
heat storage is calculated hourly as follows:
MNEGHS = 0.00695*(PACKF/2.0)*(-RELTMP) (17)
where:
MNEGHS = maximum negative heat storage, inches of water equivalent
PACKF = water equivalent of the frozen contents of the snowpack,
inches
RELTMP = air temperature above freezing in degrees F
The accummulation of the negative heat storage is calculated hourly from the
following empirical relationship:
NEGHT = 0.0007*(PAKTMP - AIRTMP)*DELT60 (18)
where:
NEGHT = potential rate of cooling of the snowpack in inches of water
equivalent per interval
PAKTMP = mean temperature of the snowpack in degrees F
AIRTMP = air temperature in degrees F
DELT60 = number of hours per interval
NEGHT is added to the negative heat storage (NEGHTS) every interval except
when limited by MNEGHS. NEGHTS is used in the parent subroutine SNOW to
calculate the temperature of the snowpack and in subroutine WARMUP to
determine the extent that the pack must be warmed to reach 32 degrees F.
4.2(1).2.7 Warm the Snowpack if Possible (subroutine WARMUP)
Purpose
This subroutine warms the snowpack to as much as 32 degrees F when possible.
153
-------
Module Section SNOW
Method
When there is negative heat storage in the pack (see subroutine COOLER for a
discussion of NEGHTS), and there is net incoming energy as calculated in
previous subroutines, then NEGHTS will decrease resulting in a warmer
snowpack and possible melt.
The calculations in this subroutine are merely accounting. They decrease
NEGHTS to a minimum of zero by subtracting the net incoming heat. If any
negative heat storage remains, then the latent heat released by the freezing
of any incoming rain is added to the pack. Since NEGHTS and all other heat
variables are in units of inches of melt, the inches of rain falling on the
pack and freezing is subtracted from NEGHTS without any conversion.
4.2(1).2.8 Melt the Pack Using Any Remaining Heat (subroutine MELTER)
Purpose
MELTER simulates the actual melting of the pack with whatever incoming heat
remains. Any heat which was not used to heat the snowpack in subroutine
WARMUP can now be used to melt the snowpack.
Method
This subroutine is also merely an accounting subroutine. The net incoming
heat has already been calculated in terms of water equivalents of melt.
Hence, any remaining incoming heat is used directly to melt the snowpa ck
either partially or entirely depending on the size of the snowpack.
4.2(1).2.9 Handle Liquid Water in the Pack (subroutine.LIQUID)
Purpose
Subroutine LIQUID first determines the liquid storage capacity of the
snowpack. It then determines how much liquid water is available to fill the
storage capacity. Any liquid water above the capacity wil 1 leave the
snowpack unless it freezes (see subroutine ICING).
Method
The liquid water holding capacity of the snowpack can be at the maximum as
specified by the parameter MWATER, at zero, or somewhere in between
depending on the density of the pack: the less dense the snowpack the
greater the holding capacity. The following relationships define the
capacity:
for RDENPF > 0.91,
PACKWC = 0.0 (19)
for 0.6 < RDENPF < 0.91,
PACKWC = MWATER*(3.0 - 3.33*RDENPF) (20)
154
-------
Module Section SNOW
for RDENPF < 0.61,
PACKWC = MWATER (21)
where:
PACKWC = liquid water holding capacity of the snowpack in
in./in.
MWATER = parameter specifing the maximum liquid water content of
the snowpack in in./in.
RDENPF = density of the snowpack relative to liquid water
MWATER is a function of the mass of ice layers, the size, the shape, and
spacing of snow crystals and the degree of channelization and honeycombing
of the snowpack.
Once PACKWC is calculated, it is compared to the available liquid water in
the pack PWSUPY. PWSUPY is calculated by summing any storage remaining at
the start of the interval, any melt, and any rain that fell on the pack
which did not freeze. If PWSUPY is more than PACKWC, then water is yielded
to the land surface from the snowpack.
4.2(1).2.10 Simulate Occurrence of Ice in the Pack (subroutine ICING)
Purpose
The purpose of subroutine ICING is to simulate the possible freezing of
water which would otherwise leave the snowpack. This freezing in turn
produces ice or frozen ground at the bottom of the snowpack. In this
subroutine, the ice can be considered to be at the bottom of the pack or
frozen in the ground below the snow portion of the pack thus extending the
total pack into the soil. This subroutine may only be applicable in certain
areas; therefore, it is user optional.
Method
The freezing of the water yield of the snowpack depends on the capacity of
the environment to freeze it. Every day at approximately 6 a.m. the capacity
is reassessed. A new value is estimated in terms of inches of melt by
multiplying the Fahrenheit degrees of the air temperature below 32.0 by
0.01. This estimate is compared with the freezing capacity if any which
remains from the previous 24-hr period. If it is greater, then the new
estimated capacity replaces the old, else the old value remains as the
potential. Any water yield that occurs freezes and is added to the ice
portion of the snowpack until the capacity is met. Any subsequent water
yield is released from the snowpack.
155
-------
Module Section SNOW
4.2(1).2.11 Melt the Pack Using Heat from the Ground (subroutine GMELT)
Purpose
The purpose of the GMELT subroutine is to simulate the melt caused by heat
conducted from the surface underlying the snowpack. This ground heat melts
the pack only from below. Therefore, melt from this process is considered
independent of other previously calculated heat influences except for an
indirect effect via the temperature of the snowpack. Unlike the other melt
processes, ground heat melts the ice portion of the snowpack first since ice
is considered to be located at the lower depths of the pack.
Method
The potential rate of ground melt is calculated hourly as a function of
snowpack temperature (PAKTMP) and a lumped parameter (MGMELT). MGMELT is
the maximum rate of melt in water equivalent caused by heat from the ground
at a PAKTMP of 32 degrees F. MGMELT would depend upon the thermal
conductivity of the soil and the normal depth of soil freezing. The
potential ground melt is reduced below MGMELT by 3 percent for each degree
that PAKTMP is below 32 degrees F to a minimum of 19 percent of MGMELT at 5
degrees F or lower. As long as a snowpack is present, ground melt occurs at
this potential rate. ;
1.2(1).2.12 Reset State Variables When Snowpack Disappears (subroutine
NOPACK)
Purpose
This code resets the state variables (for example, SNOCOV) when the snowpack
completely disappears.
Method
The frozen contents of the snowpack required for complete areal cover of
snow (COVINX) is set to a tenth of the maximum value (COVIND). All other
variables are either set to zero or the "undefined" value of -1.0E38.
156
-------
Module Section PWATER
4.2(1).3 Simulate Water Budget for a Pervious Land Segment (Section PWATER
of Module PERLND)
Purpose
PWATER is used to calculate the components of the water budget, primarily to
predict the total runoff from a pervious area. PWATER is the key component
of module PERLND; subsequent major sections of PERLND (eg. SEDMNT) depend on
the outputs of this section.
Background
The hydrologic processes that are modeled by PWATER are illustrated in
Figure 4.2(1).3-1. The algorithms used to simulate these land related
processes are the product of over 15 yr of research and testing. They are
based on the original research for the LANDS subprogram of the Stanford
Watershed Model IV (Crawford and Linsley 1966). LANDS has been incorporated
into many models and used to successfully simulate the hydrologic responses
of widely varying watersheds. The equations used in module section PWATER
are nearly identical to the ones in the current version of LANDS in the PTR
Model (Crawford and Donigian 1973), HSP (Hydrocomp 1976), and the ARM and
NFS Models (Donigian and Crawford 1976 a,b). However, some changes have been
made to LANDS to make the algorithms internally more amenable to a range of
calculation time steps. Also, many of the parameter names have been changed
to make them more descriptive, and some can be input on a monthly basis to
allow for seasonal variations.
Data Requirements and Manipulation
The number of time series required by module section PWATER depends on
whether snow accumulation and melt are considered.
When such conditions are not considered, only potential evapotranspiration
and precipitation are required.
However, when snow conditions are considered, air temperature, rainfall,
snowcover, water yield, and ice content of the snowpack are also required.
Also, the evaporation data are adjusted when snow is considered. The input
evaporation values are reduced to account for the fraction of the land
segment covered by the snowpack (determined from the generated time series
for snow cover), with an allowance for the fraction of area covered by
coniferous forest which, it is assumed, can transpire through any snow
cover. Furthermore, PET is reduced to zero when air temperature is below the
parameter PETMIN. If air temperature is below PETMAX but above PETMIN, PET
will be reduced to 50% of the input value, unless the first adjustment
already reduced it to less than this amount.
The estimated potential evapotranspiration (PET) is used to calculate actual
ET in subroutine group EVAPT.
157
-------
Ul
oo
I I I I I I I I 1 ill I I
Precipitation
Interception
Evapotranspiration /
t /
• 'flfe-T^^
• V .rflo» •
Groundwater table
Figure 4.2(1).3-1 Hydrologic cycle
-------
Module Section PWATER
Approach
Figure 4.2(1).3-2 represents the fluxes and storages simulated in module
section PWATER. The time series SUPY representing moisture supplied to the
land segment includes rain, and when snow conditions are considered, rain
plus water from the snowpack. SUPY is then available for interception.
Interception storage is water retained by any storage above the overland
flow plane. For pervious areas, interception storage is mostly on
vegetation. Any overflow from interception storage is added to the
optionally supplied time series of surface external lateral inflow to
produce the total inflow into the surface detention storage.
Inflow to the surface detention storage is added to existing storage to make
up the water available for infiltration and runoff. Moisture which directly
infiltrates moves to the lower zone and groundwater storages. Other water
may go to the upper zone storage, may be routed as runoff from surface
detention or interflow storage, or may stay on the overland flow plane, from
which it r.uns off or infiltrates at a later time.
The processes of infiltration and overland flow interact and occur
simultaneously in nature. Surface conditions such as heavy turf on mild
slopes restrict the velocity of overland flow and reduce the total quantity
of runoff by allowing more time for infiltration. Increased soil moisture
due to prolonged infiltration will in time reduce the infiltration rate
producing more overland flow. Surface detention will modify flow. For
example, high intensity rainfall is attenuated by storage and the maximum
outflow rate is reduced. The water in the surface detention may also later
infiltrate reoccurring as interflow, or it can be contained in upper zone
storage.
Water infiltrating through the surface and percolating from the upper zone
storage to the lower zone storage may flow to active groundwater storage or
may be lost by deep percolation. Active groundwater eventually reappears as
baseflow, but deep percolation is considered lost from the simulated system.
Lateral external inflows to interflow and active groundwater storages are
also possible in section PWATER. One may wish to use this option if an
upslope land segment is significantly different to merit separating it from
a downslope land segment and no channel exists between them. This capability
was not included in the previous models.
Not only are flows important in the simulation of the water budget, but so
are storages. As stated, soil storage affects infiltration. The water
holding capacity of the two soil storages, upper zone and lower zone, in
module section PERLND is defined in terms of nominal capacities. Nominal,
rather than absolute capacities, serve the purpose of smoothing any abrupt
change that would occur if an absolute capacity is reached. Such capacities
permit a smooth transition in hydrologic performance as the water content
fluctuates.
159
-------
1
c
/TAET \ -|
^total actual I
inter- i
•^B
/"SUPY \
/precipor rain]
\4- snowpack 1
\water yi«ldy
i
ception
evapo- CEPS
V. ration y interception ET-evapotransplration
>. N^/ storage
fsURlTl
external
lateral
surface
\inflow / ,
\x /
i
51
sur
dete
sto
i
/ CEPO \
{ interception 1
\ outflow /
\ /
\ Inflow* / surface
N^ / outflow
r ^L f
IRS *^J
faoe .
ntton i IFWLI | s ^x
•^o0* external _liiA
lateral . f WO
interflow inlerflo*
;__..i , ______ •HITTIftlll
V input / ouniow
1 \X , IFW5 V /
A 1 " interflowa ^f ••
7 x^N. storage
/IFWI\
interflow
input
» * from
____ 1 «urfneet
Figure 4.ZC 0.3-2. Flow diagram of water movement and storages
modeled in the PWMER section of the PERLND Application Module.
160
-------
UZET
upper
zone
ET
6
L7ET
lower
rone
ET
o
AGWET
qround
water
ET
o
INFIL
VnflltratIon
/
/infiltration^
\percolation tej
\towerzones,
i •
6
UZI
upper
zone
inflow
UZS
upper zone
storage
LZi
lower
zone
inflow
AOWI
active
qround
water
irvfioiv
PIGWI
deep
ercdation/
lateral
qround
water
inflow^
lower zone
storage
AOWS
active
qround water
storaqe
AOWO
qround
water
outflow!
Fiqure 4.e CI )>*. (Continued)
161
-------
Storages also affect evapotranspiration loss. Evapotranspiration can be
simulated from interception storage, upper and lower zone storages, active
groundwater storage, and directly from baseflow.
Storages and flows can also be instrumental in the transformation and
movement of chemicals simulated in the agri-chemical module sections. Soil
moisture levels affect the adsorption and transformations of pesticides and
nutrients. Soil moisture contents may vary greatly over a land segment.
Therefore, a more detailed representation of the moisture contents and
fluxes may be needed to simulate the transport and reaction of agricultural
chemicals. Following the ARM Model, HSPF permits a segment to be further
divided into conceptual "blocks" which represent the areal variations of the
watershed in more detail. Further explanation of these conceptual areal
blocks or zones is given in the ARM Model report (Donigian and Crawford
1976a). ARM uses five blocks but HSPF allows the user to specify from one
to five.
The following subroutine descriptions will explain in more detail the
algorithms of the PWATER module .section. Further detail can be found in the
pseudo code and the reports cited above.
.2(1).3.1 Simulate Interception (subroutine ICEPT)
Purpose
The purpose of this code is to simulate the interception of moisture by
vegetal or other ground cover. Moisture is supplied by precipitation, or
under snow conditions, it is supplied by the rain not falling on the
snowpack plus the water yielded by the snowpack.
Method
The user may supply the interception capacity on a monthly basis to account
for seasonal variations, or may supply one value designating a fixed
capacity. The interception capacity parameter can be used to designate any
retention of moisture which does not infiltrate or reach the overland flow
plane. Typically for pervious areas this capacity represents storage on
grass blades, leaves, branches, trunks, and stems of vegetation.
Moisture exceeding the interception capacity overflows the storage and is
ready for either infiltration or runoff as determined by subroutine group
SURFAC*. Water held in interception storage is removed by evaporation; the
amount is determined in subroutine EVICEP.
162
-------
Module Section PWATER
4.2(1).3.2 Distribute the Water Available for Infiltration and Runoff
(subroutine SURFAC)
Purpose
Subroutine SURFAC determines what happens to the moisture on the surface of
the land. It may infiltrate, go to the upper zone storage or interflow
storage, remain in surface detention storage, or run off.
Method
The algorithms which simulate infiltration represent both the continuous
variation of infiltration rate with time as a function of soil moisture and
the areal variation of infiltration over the land segment. The equations
representing the dependence of infiltration on soil moisture are based on
the work of Philip (1957) and are derived in detail in the previously cited
reports.
The infiltration capacity, the maximum rate at which soil will accept
infiltration, is a function of both the fixed and variable characteristics
of the watershed. Fixed characteristics include primarily soil permeability
and land slopes, while variables are soil surface conditions and soil
moisture content. Fixed and variable characteristics vary spatially over
the land segment. A linear probability density function is used to account
for areal variation. Figure 4.2(1).3-3 represents the
infiltration/interflow/surface runoff distribution function of section
PWATER. Careful attention to this figure and Figure 4.2(1).3-2 will
facilitate understanding of subroutine SURFAC and the subordinate
subroutines DISPOS, DIVISN, UZINF, and PROUTE.
The infiltration distribution represented by Figure 4.2(1).3-3 is focused
around the two lines which separate the moisture available to the land
surface (MSUPY) into what infiltrates and what goes to interflow. A number
of the variables that are used to determine the location of lines I and II
are calculated in subroutine SURFAC. They are calculated by the following
relationships:
IBAR = (INFILT/(LZS/LZSN)**INFEXP)*INFFAC O>
IMAX = INFILD«IBAR
IMIN = IBAR - (IMAX - IBAR)
RATIO = INTFW*(2.0*»(LZS/LZSN))
where:
IBAR = mean infiltration capacity over the land segment in
in./interval
INFILT r infiltration parameter in in./interval
163
-------
5
MSUPY
"5
«
J
o
block number
lint II (interflow*
infiltration copocity)
•line I (infiltration capacity)
I I
50
% of &rea
too
potential surface
detention/runoff
potential
interflow inflow
>• potential direct runoff
Fiqure 4.2(I).3-J Determination of infiltration
and interflow inflow
164
-------
Module Section PWATER
LZS = lower zone storage in inches
LZSN = parameter for lower zone nominal storage in inches
INFEXP = exponent parameter greater than one
INFFAC = factor to account for frozen ground effects, if applicable
TMAX = maximum infiltration capacity in in./interval
INFILD = parameter giving the ratio of maximum to mean infiltration
capacity over the land segment
IMIN = minimum infiltration capacity in in./interval
RATIO = ratio of the ordinates of line II to line I
INTFW = interflow inflow parameter
The factor that reduces infiltration (and also upper zone percolation) to
account for the freezing of the ground surface (INFFAC) is 1.0 if icing is
not simulated. When icing occurs, the factor is 1.0 minus the water
equivalent of ice in of the snowpack to a minimum of 0.1.
The parameter INTFW can be input on a monthly basis to allow for variations
throughout the year.
When the land segment is separated into conceptual areal blocks as
designated by the vertical subdivisions diagrammed in Figure 4.2(1).3-3,
corresponding IMAX and IMIN values must be determined for each block:
IMNB = IMIN + (BLK - 1)*(IMAX - IMIN)/NBLKS (5)
IMXB = IMNB + (IMAX - IMIN)/NBLKS (6)
where:
IMNB = minimum infiltration capacity for block BLK in in./interval
BLK = block number
NBLKS = total number of blocks being simulated
IMXB = maximum infiltration capacity for block BLK in in./interval
4.2(1)3.2.1 Dispose of Moisture Supply
(subroutine DISPOS)
Purpose
Subroutine DISPOS determines what happens to the moisture supply (MSUPY) on
either an individual block of the land segment (if NBLKS is greater than 1)
or on the entire land segment (if NBLKS is equal to 1).
Method
This subroutine calls subordinate routines DIVISN, UZINF, and PROUTE. DIVISN
is called to determine how much of MSUPY falls above and below line I in
Figure 4.2(1).3-3. The quantity under this line is considered to be
165
-------
Module Section PWATER
infiltrated. The amount over the line but under the MSUPY line (the entire
shaded portion) is the potential direct runoff '(PDRO), which is the combined
increment to interflow, and upper zone storage plus the quantities which
will stay on the surface and run off. PDRO is subdivided by line II. The
ordinates of line II are found by multiplying the ordinates of line I by
RATIO (see subroutine SURFAC for definition). The quantity underneath both
line II and the MSUPY line but above line I is called potential interflow
inflow. This consists of actual interflow plus an increment to upper zone
storage. Any amount above line II but below the MSUPY (potential surface
detention/runoff) is that portion of the moisture supply which stays on the
surface and is available for overland flow routing, plus a further increment
to upper zone storage. The fractions of the potential interflow inflow and
potential surface detention/runoff which are combined to compose the upper
zone inflow are determined in subroutine UZINF.
1.2(1).3.2.1.2 Compute Inflow to Upper Zone (subroutines UZINF1 and UZINF2)
Purpose
The purpose of this code is to compute the inflow to the upper zone when
there is some potential direct runoff (PDRO). PDRO, which is determined in
subroutine DISPOS, will either enter the upper zone storage or be available
for either interflow or overland flow. This subroutine determines what
amount, if any, will go to the upper zone storage.
Method
The fraction of the potential direct runoff which is inflow to the upper
zone storage is a function of the ratio (UZRAT) of the storage to the
nominal capacity. Figure 4.2(1).3-4 diagrams this relationship. The
equations used to define this curve follow:
FRAC = 1 - (UZRAT/2)*(1/(4 - UZRAT))**(3 - UZRAT) (7)
for UZRAT less than or equal to two. For UZRAT greater than two,
FRAC = (0.5/(UZRAT-1))**(2*UZRAT-3) (8)
where:
FRAC = fraction of PDRO retained by the upper zone storage
UZRAT = UZS/UZSN
Since UZS and FRAC are dynamically affected by the inflow process it becomes
desirable when using particularly large time steps to integrate over the
interval to find the inflow to the upper zone. This is done in subroutine
UZINF1. The solution is simplified by assuming that inflow to and outflow
from the upper zone is handled separately. Considering inflow, the following
differential equation results:
166
-------
Thus
Module Section PWATER
d(UZS)/dt = (d(UZRAT)/dt)*UZSN = PDRO*FRAC
(9)
d(UZRAT)/FRAC = (PDRO/UZSN)*dt
Now taking the definite integral of both sides of the equation:
•UZRAT
INTGRL
•L
(10)
(11)
where:
t1 = time at start of interval
t2 = time at end of interval
|-i^^SaSS?rf«Srf?:a.
over
e
time step, due to its being filled (Figure 4.2(1).3-4).
CC
u.
l.OO
.60
.60
.40
.ZO
.00
X
\
1.0 2.0
UZRKT
Figure 4.2(1).3-4 Fraction of the potential direct runoff
retained by the upper zone (FRAC) as a function of the
upper zone soil moisture ratio (UZRAT)
167
-------
Module Section PWATER
M.2(1).3.2.1.3 Determine Surface Runoff (subroutine PROUTE)
Purpose
The purpose of subroutine PROUTE is to determine how much potential surface
detention runs off in one simulation interval.
Method of Routing
Overland flow is treated as a turbulent flow process. It is simulated using
the Chezy-Manning equation and an empirical expression which relates outflow
depth to detention storage. A more detailed explanation and derivation can
be found in the reports cited in the initial background discussion. The rate
of overland flow discharge is determined by the equations:
for SURSM < SURSE
SURO = DELT60*SRC*(SURSM*(1.0 + 0.6(SURSM/SURSE)**3)**1.67 (12)
for SURSM >= SURSE
SURO = DELT60*SRC*(SURSM*1.6)**1.67
where:
SURO = surface outflow in in./interval
DELT60 = DELT/60.0 (hr/interval)
SRC = routing .variable, described below
SURSM = mean surface detention storage over the time interval in
inches
SURSE = equilibrium surface detention storage (inches) for current
supply rate
DELT60 makes the equations applicable to a range of time steps (DELT). The
first equation represents the case where the overland flow rate is
increasing, and the second case where the surface is at equilibrium or
receding. Equilibrium surface detention storage is calculated by:
SURSE = DEC*SSUPR**0.6 (13)
where:
DEC = calculated routing variable, described below
SSUPR = rate of moisture supply to the overland flow surface
There are two optional ways of determining SSUPR and SURSM. One option -
the same method used in the prior models HSP.ARM and NPS - estimates SSUPR
by subtracting the surface storage at the start of the interval (SURS) from
the potential surface detention (PSUR) which was determined in subroutine
DISPOS. The units of SSUPR are inches per interval. SURSM is estimated as
the mean of SURS and PSUR. The other option estimates SSUPR by the same
method except that the result is divided by DELT60 to obtain a value with
168
-------
Module Section PWATER
/
units of inches per hour. SURSM is set equal to SURS. This option has not
been used in prior models, but is dimensionally consistent for any time
step.
The variables DEC and SRC are calculated daily in subroutine SURFAC, but
their equations will be given here since they pertain to routing. They are:
DEC = 0.00982*(NSUR*LSUR/SQRT(SLSUR))**0.6 (1H)
SRC = 1020.0*(SQRT(SLSUR)/(NSUR*LSUR)) (15)
where:
NSUR = Manning's n for the overland flow plane
LSUR = length of the overland flow plane in ft
SLSUR = slope of the overland flow plane in ft/ft
NSUR can be input on a monthly basis to allow for variations in roughness of
the overland flow plane throughout the year.
4.2(1).3.3 Simulate Interflow (subroutine INTFLW)
Purpose
Interflow can have an important influence on storm hydrographs particularly
when vertical percolation is retarded by a shallow, less permeable soil
layer. Additions to the interflow component are retained in storage or
routed as outflow from the land segment. Inflows to the interflow component
may occur from the surface or from upslope external lateral flows. The
purpose of this subroutine is to determine the amount of interflow and to
update the storage.
Method of Determining Interflow
The calculation of interflow outflow assumes a linear relationship to
storage. Thus outflow is a function of a recession parameter, inflow, and
storage. Moisture that remains will occupy interflow storage. Interflow
discharge is calculated by:
IFWO = (IFWK1«INFLO) + (IFWK2*IFWS) (16)
where:
IFWO = interflow outflow in in./interval
INFLO = inflow into interflow storage in in./interval
IFWS = interflow storage at the start of the interval in inches
IFWK1 and IFWK2 are variables determined by:
IFWK1 s 1.0 - (IFWK2/KIFW)
IFWK2 a 1.0 - EXP(-KIFW)
169
-------
Module Section PWATER
and
KIFW = -ALOG(IRC)«DELT60/2l|.0 (19)
where:
IRC = interflow recession parameter, per day
DELT60 = number of hr/interval
24.0 = number of hours per day
EXP = Fortran exponential function
ALOG = fortran natural logarithm function
When a pervious land segment is divided into more than one block, the
algorithms are applied separately to each block. IRC is the ratio of the
present rate of interflow outflow to the value 24 hours earlier, if there
was no inflow. IRC can be input on a monthly basis to allow for variations
in soil properties throughout the year.
.2(1).3.*» Simulate Upper Zone Behavior (subroutine UZONE)
Purpose
This subroutine and the subsidiary subroutine UZONES are used to calculate
the water percolating from the upper zone. Water not percolated remains in
upper zone storage available for evapotranspiration in subroutine ETUZON.
Method of Determining Percolation
The upper zone inflow calculated in DISPOS is first added to the upper zone
storage at the start of the interval to obtain the total water available for
percolation from the upper zone.
Percolation only occurs when UZRAT minus LZRAT is greater than 0.01. When
this happens, percolation from the upper zone storage is calculated by the
empirical expression:
PERC = 0.1*INFILT*INFFAC*UZSN*(UZRAT - LZRAT)**3 (20)
where:
PERC = percolation from the upper zone in in./interval
INFILT = infiltration parameter in in./interval
INFFAC = factor to account for frozen ground, if any,
UZSN = parameter for upper' zone nominal storage in inches
UZRAT = ratio of upper zone storage to UZSN
LZRAT = ratio of lower zone storage to lower zone
nominal storage (LZSN)
The upper zone nominal capacity can be input on a monthly basis to allow for
variations throughout the year. The monthly values are interpolated to
obtain daily values. When a pervious land segment is divided into more than
one block, the algorithm is applied separately to each block.
170
-------
Module Section PWATER
4.2(1).3.5 Simulate Lower Zone Behavior (subroutine LZONE)
Purpose
This subroutine determines the quantity of infiltrated and percolated water
which enters the lower zone. The infiltrated moisture supply is determined
in subroutine DISPOS. The percolated moisture from the upper zone is found
in subroutine UZONE.
Method
The fraction of the direct infiltration plus percolation that enters the
lower zone storage (LZS) is based on the lower zone storage ratio of
LZS/LZSN where LZSN is the lower zone nominal capacity. The inflowing
fraction is determined empirically by:
LZFRAC = 1.0 - LZRAT*(1.0/(1.0 + INDX))**INDX (21)
when LZRAT is less than 1.0, and by
LZFRAC = (1.0/O.O + INDX))**INDX (22)
when LZRAT is greater than 1.0. INDX is defined by:
INDX = 1.5*ABS(LZRAT - 1.0) + 1.0 (23)
where:
LZFRAC = fraction of infiltration plus percolation entering LZS
LZRAT = LZS/LZSN
ABS = function for determining absolute value
These relationships are plotted in Figure 4.2(1).3-5. The fraction of the
moisture supply remaining after the surface, upper zone, and lower zone
components are subtracted is added to the groundwater storages.
1.2(1).3.6 Simulate Groundwater Behavior (subroutine GWATER)
Purpose
The purpose of this subroutine is to determine the amount of the inflow to
groundwater that is lost to deep or inactive groundwater and to determine
the amount of active groundwater outflow. These two fluxes will in turn
affect the active groundwater storage.
Method of Determining Groundwater Fluxes
The quantity of direct infiltration plus percolation from the upper zone
which does not go to the lower zone (determined in subroutine LZONE) will be
inflow to either inactive or active groundwater. The distribution to active
171
-------
Module Section PWATER
X
0.5
1.0 I.B
LZRAT
Figure 4.2(1).3-5 Fraction of infiltration plus percolation
entering lower zone storage
and inactive groundwater is user designated by parameter DEEPFR. DEEPFR is
that fraction of the groundwater inflow which goes to inactive groundwater.
The remaining portion of the percolating water and all external lateral
inflow if any make up the total inflow to the active groundwater storage.
The outflow from active groundwater storage is based on a simplified model.
It assumes that the discharge of an aquifer is proportional to the product
of the cross-sectional area and the energy gradient of the flow. Further, a
representative cross-sectional area of flow is assumed to be related to the
groundwater storage level at the start of the interval. The energy gradient
is estimated as a basic gradient plus a variable gradient that depends on
past active groundwater accretion.
Thus, the groundwater outflow is estimated by:
KVARY*GWVS)*AGWS
AGWO = KGW*(1.0
where:
AGWO
KGW
KVARY
(24)
GWVS
AGWS
active groundwater outflow in in./interval
groundwater outflow recession parameter, per interval
parameter which can make active groundwater storage to outflow
relation nonlinear in per inches
index to groundwater slope in inches
active groundwater storage at the start of the interval in inches
GWVS is increased each interval by the inflow to active groundwater but is
also decreased by 3 percent once a day. It is a measure of antecedent
active groundwater inflow. KVARY is introduced to allow variable groundwater
recession rates. When KVARY is nonzero, a semilog plot of discharge vs. time
is nonlinear. This parameter adds flexibility in groundwater outflow
simulation which is useful in simulating many watersheds.
172
-------
Module Section PWATER
The parameter KGW is calculated by the Run Interpreter using the relationship:
KGW s 1.0 - (AGWRC)**(DELT60/2U.O) (25)
where:
AGWRC = daily recession constant of groundwater flow,
if KVARY or GWVS =0.0
That is, the ratio of current groundwater discharge
to groundwater discharge 24-hr earlier
DELT60 = hr/interval
.2(1).3.7 Simulate Evapotranspiration
(subroutine EVAPT)
Purpose
The purpose of EVAPT and its subordinate subroutines is to simulate
evaporation and evapotranspiration fluxes from all zones of the pervious
land segment. Since in most hydrologic regimes the volume of water that
leaves a watershed as evapotranspiration exceeds the total volume of
streamflow, this is an important aspect of the water budget.
Method of Determining Actual Evapotranspiration
There are two separate issues involved in estimating evapotranspiration
(ET). First, potential ET must be estimated. ET potential or demand is
supplied as an input times series, typically using U.S. Weather Bureau Class
A pan records plus an adjustment factor. The data are further adjusted for
cover and temperature in the parent subroutine PWATER. Second, actual ET
must be calculated, usually as a function of moisture storages and the
potential. The actual ET is estimated by trying to meet the demand from
five sources in the order described below. The sum of the ET from these five
sources is the total actual evapotranspiration from the land segment.
Subroutine ETBASE
The first source from which ET can be taken is the active groundwater
outflow or baseflow. This simulates effects such as ET from riparian
vegetation in which groundwater is withdrawn as it enters the stream. The
user may specify by the parameter BASETP the fraction, if any, of the
potential ET that can be sought from the baseflow. That portion can only be
fulfilled if outflow exists. Any remaining potential not met by actual
baseflow evaporation will try next to be satisfied in subroutine EVICEP.
173
-------
Module Section PWATER
Subroutine EVICEP
Remaining potential ET then exerts its demand on the water in interception
storage. Unlike baseflow, there is no parameter regulating the rate of ET
from interception storage. The demand will draw upon all of the interception
storage unless the demand is less than the storage. When the demand is
greater than the storage, the remaining demand will try to be satisfied in
subroutine ETUZON.
Subroutine ETUZON
There are no special ET parameters for the upper zone, but rather ET is
based on the moisture in storage in relation to its nominal capacity. Actual
evapotranspiration will occur from the upper zone storage at the remaining
potential demand if the ratio of UZS/UZSN, upper zone storage to nominal
capacity, is greater than 2.0. Otherwise the remaining potential ET demand
on the upper zone storage is reduced; the adjusted value depends on
UZS/UZSN. Subroutine ETAGW will attempt to satisfy any remaining demand.
Subroutine ETAGW
Like ET from baseflow, actual evapotranspiration from active groundwater is
regulated by a parameter. The parameter AGWETP is the fraction of the
remaining potential ET that can be sought from the active groundwater
storage. That portion of the ET demand can be met only if there is enough
active groundwater storage to satisfy it. Any remaining potential will try
to be met in subroutine ETLZON.
Subroutine ETLZON
The lower zone is the last storage from which ET is drawn.
Evapotranspiration from the lower zone is more involved than that from the
other storages. ET from the lower zone depends upon vegetation
transpiration. Evapotranspiration opportunity will vary with the vegetation
type, the depth of rooting, density of the vegetation cover, and the stage
of plant growth along with the moisture characteristics of the soil zone.
These influences on the ET opportunity are lumped into the LZETP parameter.
Unlike the other ET parameters LZETP can be input on a monthly basis to
account for temporal changes in the above characteristics.
If the LZETP parameter is at its maximum value of one, representing near
complete areal coverage of deep rooted vegetation, then the potential ET for
the lower zone is equal to the demand that remains. However, this is
normally not the case. Usually vegetation type and/or rooting depths will
vary over the land segment. To simulate this, a linear probability density
function for ET opportunity is assumed (Figure 4.2(1).3-6). This approach
174
-------
Module Section PWATER
is similar to that used to handle areal variations in
infiltration/pereolation capacity.
The variable RPARM, the index to maximum ET opportunity, is estimated by:
(26)
RPARM = 0.25/CI.O - LZETP)*(LZS/LZSN)*DELT60/24.0
where:
RPARM = maximum ET opportunity in in./interval
LZETP = lower zone ET parameter
LZS = current lower zone storage in inches
LZSN = lower zone nominal storage parameter in inches
DELT60 = hr/interval
• V
"jE REMPET
|1 Eya potranspiration 11
9O
IOO
Percent oi brea with
Opportunity Equal 1o or less 1h
-------
Module Section SEDMNT
4.2(1).4 Simulate Production and Removal of Sediment (Section SEDMNT of
Module PERLND)
Purpose
Module section SEDMNT simulates the production and removal of sediment from
a pervious land segment. Sediment can be considered to be inorganic,
organic, or both; the definition is up to the user.
Sediment from the land surface is one of the most common pollutants of
waters from urban, agricultural, and forested lands. It can muddy waters,
cover fish eggs, and limit the capacity of reservoirs. Nutritious and toxic
chemicals can be carried by it.
Approach
The equations used to produce and remove sediment are based on the ARM and
NPS Models (Donigian and Crawford 1976 a,b). The algorithms representing
land surface erosion in these models were derived from a sediment model
developed by Moshe Negev (Negev 1967) and influenced by Meyer and Wischmeier
(1969) and Onstad and Foster (1975). The supporting management practice
factor which has been added to the soil detachment by rainfall equation was
based on the "P" factor in the Universal Soil Loss Equation (Wischmeier and
Smith 1965). It was introduced in order to better evaluate agricultural
conservation practices. The equation which represents the scouring of the
matrix soil, which is not included in ARM or NPS, was derived from Negev's
method for simulating gully erosion.
Figure 4.2(1).4-1 shows the detachment, attachment, and removal involved in
the erosion processes on the pervious land surface, while Figure 4.2(1).4-2
schematically represents the fluxes and storages used to simulate these
processes. Two of the sediment fluxes, SLSED and NSVI, are added directly
to the detached sediment storage variable DETS in the parent subroutine
SEDMNT while the other fluxes are computed in subordinate subroutines. SLSED
represents external lateral input from an upslope land segment. It is a
time series which the user may optionally specify. NVSI is a parameter that
represents any net external additions or removals of sediment caused by
human activities or wind.
Removal of sediment by water is simulated as washoff of detached sediment in
storage (WSSD) and scour of matrix soil (SCRSD). The washoff process
involves two parts: the detachment/attachment of sediment from/to the soil
matrix and the transport of this sediment. Detachment (DET) occurs by
rainfall. Attachment occurs only on days without rainfall; the rate of
attachment is specified by parameter AFFIX. Transport of detached sediment
is by overland flow. The scouring of the matrix soil includes both pick up
and transport by overland flow combined into one process.
176
-------
5LSEO
lateral
input of
sediment
to sur-
face
X
/net vertical
sediment
input
v
DETS
detached
sediment
storage
/ft^T\
I sediment )
^attachment/
^/
WSSD
washoffof
detached
sediment
*y water
X
/ detachment
\ of soil bv
\ rainfall
i -
a
total
removal
of sol\ 4
sediment
from sur-
face by
water
\
soil matrix
(assumed to
have
unlimited
storage)
'9CR5D'
scourpf
matri>i
soil by
water
Fiqure4.ZO).4-Z. Flow diaqram for SEDMNT section
of PERLND Application Module.
178
-------
Module Section SEDMNT
Module section SEDMNT has two options for simulating washoff of detached
sediment and scour of soil. One uses subroutine SOSED1 which is identical to
the method used in the ARM and the NFS Models. However, some equations used
in this method are dimensionally nonhomogeneous, and it has only been used
with 15- and 5-min intervals. The results obtained are probably highly
dependent on the simulation time step. The other option uses subroutine
SOSED2 which is dimensionally homogeneous and is, theoretically, less
dependent on the time step. However, it has not been tested.
t.2(1).4.1 Detach Soil By Rainfall (subroutine DETACH)
Purpose
The purpose of DETACH is to simulate the splash detachment of the soil
matrix by falling rain.
Method of Detaching Soil by Rainfall
Kinetic energy from rain falling on the soil detaches particles which are
then available to be transported by overland flow. The equation that
simulates detachment is:
DET = DELT60*(1.0 - CR)*SMPF*KRER*(RAIN/DELT60)**JRER (1)
where:
DET = sediment detached from the soil matrix by rainfall in
tons/acre per interval
DELT60 = number of hr/interval
CR = fraction of the land covered by snow and other cover
SMPF = supporting management practice factor
KRER = detachment coefficient dependent on soil properties
RAIN = rainfall in in. /interval
JRER = detachment exponent dependent on soil properties
The variable CR is the sum of the fraction of the area covered by the
snowpack (SNOCOV), if any, and the fraction that is covered by anything else
but snow (COVER). SNOCOV is computed by section SNOW. COVER is a parameter
which for pervious areas will typically be the fraction of the area covered
by vegetation and mulch. It can be input on a monthly basis.
4.2(1).4.2 Remove by Surface Flow Using Method 1 (subroutine SOSED1)
Purpose
Subroutines SOSED1 and SOSED2 perform the same task *f >* ^^f scouring
methods. They simulate the washoff of the detached sediment and the scouring
of the soil matrix.
179
-------
Module Section SEDMNT
Method
When simulating the washoff of detached sediment, the transport capacity of
the overland flow is estimated and compared to the amount of detached
sediment available. The transport capacity is calculated by the equation:
STCAP = DELT60*KSER*((SURS + SURO)/DELT60)**JSER (2)
where:
STCAP = capacity for removing detached sediment in
tons/acre per interval
DELT60 = hr/interval
KSER = coefficient for transport of detached sediment
SURS = surface water storage in inches
SURO = surface outflow of water in in./interval
JSER = exponent for transport of detached sediment
When STCAP is greater than the amount of detached sediment in storage,
washoff is calculated by:
WSSD = DETS*SURO/(SURS + SURO) (3)
If the storage is sufficient to fulfill the transport capacity, then
the following relationship is used:
WSSD = STCAP*SURO/(SURS + SURO) (4)
where:
WSSD =. washoff of detached sediment in tons/acre per interval
DETS = detached sediment storage in tons/acre
WSSD is then subtracted from DETS.
Transport and detachment of soil particles from the soil matrix is simulated
with the following equation:
SCRSD = SURO/(SURS + SURO)*DELT60*KGER*((SURS + SURO)/DELT60)**JGER (5)
where:
SCRSD = scour of matrix soil in tons/acre per interval
KGER = coefficient for scour of the matrix soil
JGER = exponent for scour of the matrix soil
The sum of the two fluxes, WSSD and SCRSD, represents the total sediment
outflow from the land segment.
The same algorithms are used for the simulation whether or not areal blocks
are used. When blocks are used, block specific values for WSSD and SCRSD
are calculated from the surface water fluxes and storages which correspond
to each block.
Subroutine SOSED1 differs from SOSED2 in that it uses the dimensionally
nonhomogeneous term (SURS + SURO)/DELT60 in the above equations, while
SOSED2 uses the homogeneous term SURO/DELT60.
130
-------
Module Section SEDMNT
U.2(1).4.3 Remove by Surface Flow Using Method 2 (subroutine SOSED2)
Purpose
The purpose of this subroutine is the same as SOSED1. They only differ in
method.
Method of Determining Removal
This method of determining sediment removal has not been tested. Unlike
subroutine SOSED1, it makes use of the dimensionally homogeneous term
SURO/DELT60 instead of (SURO+SURS)/DELT60.
The capacity of the overland flow to transport detached sediment is
determined in this subroutine by:
STCAP = DELT60*KSER*(SURO/DELT60)**JSER (6)
When STCAP is more than the amount of detached sediment in storage, the flow
washes off all of the detached sediment storage (DETS). However, when STCAP
is less than the amount of detached sediment in storage, the situation is
transport limiting, so WSSD is equal to STCAP.
Direct detachment and transport of the soil matrix by scouring (eg.
gullying) is simulated with the equation:
SCRSD = DELT60*KGER*(SURO/DELT60)**JGER (7)
Definitions of the above terms can be found in subroutine SOSED2. The
coefficents and exponents will have different values than in subroutine
SOSED1 because they modify different variables.
M.2(1).4.4 Simulate Re-attachment of Detached Sediment (subroutine ATTACH)
Purpose
Subroutine ATTACH simulates the re-attachment of detached sediment (DETS) on
the surface (soil compaction).
Method
Attachment to the soil matrix is simulated by merely reducing DETS. Since
the soil matrix is considered to be unlimited, no addition to the soil
matrix is necessary when this occurs. DETS is diminished at the start of
each day that follows a day with no precipitation by multiplying it by (1.0
- AFFIX), where AFFIX is a parameter. This represents a first order rate of
reduction of the detached soil storage.
181
-------
Module Section PSTEMP
4.2(1).5 Estimate Soil Temperatures (Section PSTEMP of Module PERLND)
Purpose
PSTEMP simulates soil temperatures for the surface, upper, and
lower/groundwater layers of a land segment for use in module section PWTGAS
and the agri-chemical sections. Good estimates of soil temperatures are
particularly important for simulating first order transformations in the
agri-chemical sections.
Method
The two methods used for estimating soil temperatures are based on the
regression equation approach in the ARM Model (Donigian, et al. 1977) and
the smoothing factor approach used in HSP QUALITY (Hydrocomp 1977) to
simulate the temperatures of subsurface flows.
Simulation of soil temperatures is done by layers which correspond to those
specified in the agri-chemical sections. The surface layer is the portion
of the land segment that affects overland flow water quality
characteristics. The subsurface layers are upper, lower, and groundwater.
The upper layer affects interflow quality characteristics while the lower is
a transition zone to groundwater. The temperature of the groundwater layer
affects groundwater quality transformations and outflow characteristics.
Lower layer and groundwater temperatures are considered approximately equal;
a single value is estimated for both layers.
Surface layer soil temperatures are estimated by the following regression
equation:
SLTMP = ASLT + BSLT*AIRTC (1)
where:
SLTMP = surface layer temperature in degrees C
ASLT = Y-intercept
BSLT = slope
AIRTC = air temperature in degrees C
Temperatures of the other layers are simulated by one of two methods. If
TSOPFG is set equal to one in the User's Control Input, the upper layer soil
temperature is estimated by a regression equation as a function of air
temperature (similar to equation above), and the lower layer/groundwater
layer temperature is specified by a parameter which can vary monthly. This
method is similar to that used in the ARM Model except that ARM relates the
upper layer temperature to the computed soil surface temperature instead of
directly to air temperature.
If TSOPFG is set equal to zero, both the upper layer and the lower layer/
groundwater layer temperatures are computed by using a mean departure from
132
-------
Module Section PWTGAS
air temperature plus a smoothing factor. The same basic equation is used
with separate state variables and parameters for each layer:
IMP = IMPS + SMO*(AIRTCS + TDIF - IMPS) (2)
where:
IMP = layer temperature at the end of the current interval in
degrees C
SMO = smoothing factor (parameter)
AIRTCS = air temperature at the start of the current interval, Deg C
TDIF = parameter which specifies the difference between the mean air
temperature and the mean temperature of the soil layer, Deg C
TMPS = layer temperature at the start of the current interval in
degrees C
The values of the parameters for any of the layer computations can be
linearly interpolated from monthly input values to obtain daily variations
throughout the year. If this variation is not desired, the user may supply
yearly values.
Estimate Water Temperature and Dissolved Gas Concentrations
(Section PWTGAS of Module PERLND)
Purpose
PWTGAS estimates the water temperature and concentrations of dissolved
oxygen and carbon dioxide in surface, interflow, and groundwater outflows
from a pervious land segment.
Method
The temperature of each outflow is considered to be the same as the soil
temperature of the layer from which the flow originates, except that water
temperature can not be less than freezing. Soil temperatures must either be
computed in module section PSTEMP or supplied directly as an input time
series. The temperature of the surface outflow is equal to the surface
layer soil temperature, the temperature of interflow to the upper layer soil
temperature, and the temperature of the active groundwater outflow equals
the lower layer and groundwater layer soil temperature.
The dissolved oxygen and carbon dioxide concentrations of the overland flow
are assumed to be at saturation and are calculated as direct functions of
water temperature. PWTGAS uses the following empirical nonlinear equation
to relate dissolved oxygen at saturation to water temperature (Committee on
Sanitary Engineering Research 1960):
SODOX = (14.652 + SOTMP*(-0.11022 +
SOTMP*(0.007991 - 0.000077774»SOTMP)))*ELEVGC (1)
133
-------
Module Section PQUAL
where:
SODOX = concentration of dissolved oxygen in surface outflow in tng/1
SOTMP = surface outflow temperature in degrees C
ELEVGC = correction factor for elevation above sea level
(ELEVGC is calculated by the Run Interpreter dependent upon
mean elevation of the segment)
The empirical equation for dissolved carbon dioxide concentration of the
overland flow (Harnard and Davis 1943) is:
SOC02 = (10**(2385.73/ABSTMP - 14.0184 + 0.0152642*ABSTMP))
*0.000316*ELEVGC*12000.0 (2)
where:
SOC02 = concentration of dissolved carbon dioxide in
surface outflow in mg C/l
ABSTMP = absolute temperature of surface outflow in degrees K
The concentrations of dissolved oxygen and carbon dioxide in the interflow
and the active groundwater flow cannot be assumed to be at saturation.
Values for these concentrations are provided by the user. He may specify a
constant value or 12 monthly values for the concentration of each of the
gases in interflow and groundwater. If monthly values are provided, daily
variation in values will automatically be obtained by linear interpolation
between the monthly values.
4.2(1).7 Simulate Quality Constituents Using Simple Relationships with
Sediment and Water Yield (Section PQUAL of Module PERLND)
Purpose
The PQUAL module section simulates water quality constituents or pollutants
in the outflows from a pervious land segment using simple relationships with
water and/or sediment yield. Any constituent can be simulated by this
module section. The user supplies the name, units and parameter values
appropriate to e'ach of the constituents that he wishes to simulate. However,
more detailed methods of simulating sediment, heat, dissolved oxygen,
dissolved carbon dioxide, nitrogen, phosphorus, soluble tracers, and
pesticide removal from a pervious land segment are available, in other
module Sections.
Approach
The basic algorithms used to simulate quality constituents are a synthesis
of those used in the NPS Model (Donigian and Crawford 1976b) and HSP QUALITY
(Hydrocomp 1977). However, some options and combinations are unique to HSPF.
134
-------
Module Section PQUAL
Figure 4.2(1).7-1 shows schematically the fluxes and storages represented in
module section PQUAL. The occurrence of a water quality constituent in both
surface and subsurface outflow can be simulated. The behavior of a
constituent in surface outflow is considered more complex and dynamic than
the behavior in subsurface flow. A constituent on the surface can be
affected greatly by adhesion to the soil and by temperature, light, wind,
and direct human influences. Section PQUAL is able to represent these
processes in a general fashion. It allows quantities in the surface outflow
to be simulated by two methods. One approach is to simulate the constituent
by association with sediment removal. The other approach is to simulate it
using basic accumulation and depletion rates together with depletion by
washoff; that is, constituent outflow from the surface is a function of the
water flow and the constituent in storage. A combination of the two methods
may be used in which the individual outfluxes are added to obtain the total
surface outflow. These approaches will be discussed further in the
descriptions of the corresponding subroutines. Concentrations of quality
constituents in the subsurface flows of interflow and active groundwater are
supplied by the user. The concentration may be linearly interpolated to
obtain daily values from input monthly values.
The user has the useful option of simulating the constituents by any
combination of these surface and subsurface outflow pathways. The outflux
from the combination of the pathways simulated will be the total outflow
from the land segment. In addition, the user is able to select the units to
be associated with the fluxes. These options give the user considerable
flexibility. For example, he may wish to simulate coliforms in units of
organisms/acre by association with sediment in the surface runoff and using
a concentration in the groundwater which varies seasonally. Or he may want
to simulate total dissolved salts in pounds per acre by direct association
with overland flow and a constant concentration in interflow and groundwater
flow.
PQUAL allows the user to simulate up to 10 quality constituents at a time.
Each of the 10 constituents may be defined as one or a combination of
the following types: QUALSD, QUALOF, QUALIF, and/or QUALGW. If a
constituent is considered to be associated with sediment, it is called a
QUALSD. The corresponding terms for constituents associated with overland
flow, interflow, and groundwater flow are QUALOF, QUALIF, and QUALGW,
respectively. However, no more than seven of any one of the constituent
types (QUALSD, QUALOF, QUALIF, or QUALGW) may be simulated in one operation.
The program uses a set of flag pointers to keep track of these associations.
For example, QSDFP(3) = 0 means that the third constituent is not associated
with sediment, whereas QSDFP(6) = 4 means that the sixth constituent is the
fourth sediment associated constituent (QUALSD). Similar flag pointer
arrays are used to indicate whether or not a quality constituent is a
QUALOF, QUALIF, or QUALGW.
185
-------
/ removal >
[by cleaning,
Idecav/ 4 wind/
CO
accumulation
SOQO
by
owerland
flou),
QUAL. means quality constttuervh
QQM
associated with
detached
•sediment
aswciaied
detached
sed.
WA5UQ5
SOQUAL
of QUM"
wsociated
urith -soil
•hrt-al
outflouj
ot QUM
associaied
ediment
\OQUAtL
of-QUAL
with
SOQS
XX
out-flow
ofQUM
with
active
water
erf QUAL
Bqure 4.ZCD.7-I Flow diaqram for PQUAL section of PERLND Application
Module.
-------
Module Section PQUAL
4.2(1).7.1 Remove by Association with Sediment (subroutine QUALSD)
Purpose
QUALSD simulates the removal of a quality constituent from a pervious land
surface by association with the sediment removal determined in module
section SEDMNT.
Method
This approach assumes that the particular quality constituent removed from
the land surface is in proportion to the sediment removal. The relation is
specified with user-input "potency factors." Potency factors indicate the
constituent strength relative to the sediment removed from the surface.
Various quality constituents such as iron, lead, and strongly adsorbed
toxicants are actually attached to the sediment being removed from the land
surface. Some other pollutants such as ammonia, organics, pathogens, and BOD
may not be extensively adsorbed, but can be considered highly correlated to
sediment yield.
For each quality constituent associated with sediment, the user supplies
separate potency factors for association with washed off and scoured
sediment (WSSD and SCRSD). Typically, the washoff potency factor would be
larger than the scour potency factor because washed off sediment is usually
finer than the scoured material and thus has a higher adsorption capacity.
Organic nitrogen would be a common example of such a constituent. The user
is also able to supply monthly potency factors for constituents that vary
somewhat consistently during the year. For instance, constituents that are
associated with spring and fall fertilization may require such monthly input
values.
Removal of the sediment associated constituent by detached sediment washoff
is simulated by:
WASHQS = WSSD*POTFW (1)
where:
WASHQS = flux of quality constituent associated with
detached sediment washoff in quantity/acre per interval
WSSD = washoff of detached sediment in tons/acre per interval
POTFW = washoff potency factor in quantity/ton
Removal of constituents by scouring of the soil matrix is similar:
SCRQS = SCRSD*POTFS (2)
where:
SCRQS = flux of quality constituent associated with scouring
of the matrix soil in quantity/acre per interval
187
-------
Module Section PQUAL
SCRSD = scour of matrix soil in tons/acre per interval
POTFS = scour potency factor in quantity/ton
WASHQS and SCRQS are combined to give the total sediment associated flux of
the constituent from the land segment, SOQS.
The unit "quantity" refers to mass units (pounds or tons in the English
system) or some other quantity, such as number of organisms for coliforms.
The unit is user specified.
1.2(1).7.2 Accumulate and Remove by a Constant Unit Rate and by Overland
Flow (subroutine QUALOF)
Purpose
QUALOF simulates the accumulation of a quality constituent on the pervious
land surface and its removal by a constant unit rate and by overland flow.
Method
This subroutine differs from the others in module section PQUAL in that the
storage of the quality constituent on the land surface is simulated. The
constituent can be accumulated and removed by processes which are
independent of storm events such as cleaning, decay, and wind erosion and
deposition, or it can be washed off by overland flow. The accumulation and
removal rates can have monthly values to account for seasonal fluctuations.
A pollution indicator such as fecal coliform from range land is an example
of a constituent with accumulation and removal rates which may need to vary
throughout the year. The concentration of the coliform in the surface runoff
may fluctuate with the seasonal grazing density, and the weather.
When there is surface outflow and some quality constituent is in storage,
washoff is simulated using the commonly used relationship:
SOQO = SQO*(1.0 - EXP(-SURO»WSFAC)) (3)
where:
SOQO = washoff of the quality constituent from the land
surface in quantity/acre per interval
SQO = storage of available quality constituent on the surface
in quantity/acre
SURO = surface outflow of water in in./interval
WSFAC = susceptibility of the quality constituent to washoff
in units of 1/in.
EXP = Fortran exponential function
188
-------
Module Section PQUAL
The storage is updated once a day to account for accumulation and removal
which occurs independent of runoff by the equation:
SQO = ACQOP + SQOS*(1.0 - REMQOP) (4)
where:
ACQOP = accumulation rate of the constituent,
quantity/acre per day
SQOS = SQO at the start of the interval
REMQOP = unit removal rate of the stored constituent,
per day
The Run Interpreter computes REMQOP and WSFAC for this subroutine according
to:
REMQOP = ACQOP/SQOLIM (5)
where:
SQOLIM = asymptotic limit for SQO as time approaches infinity
(quantity/acre), if no washoff occurs
and
WSFAC = 2.30/WSQOP (6)
where:
WSQOP = rate of surface runoff which results in a 90 percent
washoff in one hour, in./hr
Since the unit removal rate of the stored constituent (REMQOP) is computed
from two other parameters, it does not have to be supplied by the user.
1.2(1).7.3 Simulate by Association with Interflow Outflow (subroutine
QUALIF)
Purpose
QUALIF is designed to permit the user to simulate the occurrence of a
constituent in interflow.
Method
The user simply specifies a concentration for each constituent which is a
QUALIF. An option permits him to supply 12 monthly values, to account for
seasonal fluctuations. In this case, the system interpolates a new value
each day.
4.2(1).7.1 Simulate by Association with Active Goundwater Outflow
(subroutine QUALGW)
Purpose
QUALGW is designed to permit the user to simulate the occurrence of a
constituent in ground water outflow.
Method
Identical to that for QUALIF
139
-------
Intro to Ag Chem Sections
Introduction to the Agri-chemical Sections
The introduction of agricultural chemicals into streams, lakes, and
groundwater from agricultural land may be detrimental. For example,
persistant fat soluble pesticides, such as DDT, have been known to
concentrate in the fatty tissue of animals causing toxic effects. Nitrogen
and phosphorus are essential plant nutrients which when introduced into
certain surface waters will increase productivity. This may or may not be
desirable depending upon management objectives. Significant productivity
results in algal blooms, but some increase in productivity will increase
fish production. Drinking water containing high nitrate concentrations may
cause methomoglobinemia in small children.
Pesticide, nitrogen, and phosphorus compounds are important to agricultural
production, but prediction of their removal from the field is necessary for
wise management of both land and water resources. HSPF can be used to
predict such outflows. The agri-chemical sections of the PERLND module of
HSPF simulate in detail nutrient and pesticide processes, both biological
and chemical, and the movement of any nonreactive tracer in a land segment.
These chemicals can also be simulated in module section PQUAL but in a
simplified manner. The dynamic and continuous processes that affect the
storages and outflow of pesticides and of nutrients from fertilized fields
should be simulated in detail to fully analyze agricultural runoff. If the
situation does not require full representation of these processes, or if
data are not available, the PQUAL subroutines could be used.
The basic algorithms in the agri-chemical sections of HSPF were originally
developed for use on agricultural lands, but can be used on other pervious
areas where pesticides and plant nutrients occur, for example, orchards,
nursery land, parks, golf courses, and forests. All pervious land contains
nitrogen and phosphorus in the soil; it is possible to use this module to
simulate the behavior of agricultural chemicals in any such area.
Comparison of HSPF and ARM
The methods used to simulate pesticide processes in the agri-chemical
sections were developed orginally for the Pesticide Transport and Runoff
(PTR) Model (Crawford and Donigian 1973), then expanded to include nutrients
in the Agricultural Runoff Management (ARM) Model (Donigian and Crawford
1976) and tested and modified in ARM Version II (Donigian, et al. 1977). In
HSPF the ARM Version II algorithms were recreated with some additional
options. (For more detail on the basic methods, refer to the above
reports.)
The differences between HSPF and ARM Model Version II should, however, be
discussed. The biggest difference is the availability of new options to
simulate soil nutrient and pesticide adsorption and desorption. Ammonium and
190
-------
Intro to Ag Chem Sections
phosphate adsorption/desorption in HSPF can be accomplished by using
Freundlich isotherms as well as by first order kinetics. Pesticides can be
adsorbed and desorbed by the two Freundlich methods used in the ARM Model or
by first order kinetics. In addition, the pesticide parameter values are
now input for each separate soil layer instead of inputting one parameter
set for all the layers. HSPF also allows the user to simulate more than one
pesticide in a run. (The ARM Model only simulates one per run.) In addition
to the percolation factors which can still be used to retard any solute
leaching from the upper layer and lower layer, a multiplication factor has
been introduced that can reduce leaching from the surface layer. Also, in
HSPF, nitrogen and phosphorus chemical and biochemical transformations can
each be simulated at different time steps to save computer time. Plant
uptake of ammonium is another new option in HSPF.
The topsoil layers can still be divided into areal blocks (see the PWATER
section). They have been used in the ARM Model to represent areal variation
in chemical concentrations over the land surface. They divide the hyrologic
responses and chemical storages in the surface and upper layers into
conceptual areal zones based on hydrologic variations. HSPF allows the user
to use the model with from one to five blocks (Section 4.2(1).3); whereas
the ARM Model uses a fixed five block setup. The capability to vary the
number of blocks allows experimentation. In some situations a land segment
may not need to be subdivided into more than one block, resulting in
substantial savings in computer costs.
Units
The fluxes and storages of chemicals modeled in these module sections are in
mass per area units. The user must supply his input in appropriate units;
kg/ha if he is using the Metric system, and Ib/ac for the English system.
Internally, most of the code does not differentiate between the unit
systems. Fluxes are determined by either proportionality constants,
fractions of chemicals in storage, or unitless concentrations. First order
kinetics makes use of proportionality constants for determining reaction
fluxes. Chemicals are transported based on the fractions of that in
storage. Freundlich adsorption/desorption is based on ppm concentrations.
Module Sections
There are five agri-chemical module sections. They are shown in the
structure chart of PERLND (No. 4.2(1)). Module section MSTLAY manipulates
water storages and fluxes calculated in module section PWATER. This section
must be run before the following sections can be run, since it supplies them
with data for simulating the storage and movement of solutes. Module section
PEST simulates pesticide behavior while NITR and PHOS simulate the plant
nutrients of nitrogen and phosphorus. Simulation of a nonreactive solute
(tracer) is accomplished in module section TRACER.
191
-------
Module Section MSTLAY
4.2(1).8 Estimate Moisture Content of Soil Layers and Fractional Fluxes
(Section MSTLAY of Module PERLND)
Purpose
This module section estimates the storages of moisture in the four soil
layers with which the agricultural chemical sections deal (Figure
4.2(1).8-1); and the fluxes of moisture between the storages. MSTLAY is
required because the moisture storages and fluxes computed by module section
PWATER can not be directly used to simulate solute transport through the
soil. For example, in PWATER, some moisture which infiltrates can reach the
ground water in a single time step (Figure 4.2(1).3-2). While this
phenomenon does not have any serious effect in simulating the hydrologic
response of a land segment, it does seriously affect the simulation of
solute transport.
Thus, MSTLAY takes the fluxes and storages computed in PWATER and adapts
them to fit the storage/flow path picture in Figure 4.2(1).8-1. The revised
storages, in inches of water, are also expressed in mass/area units (that
is, Ib/acre or kg/ha) for use in the adsorption/desorption calculations.
Method
Figure 4.2(1).8-1 schematically diagrams the moisture storages and fluxes
used in subroutine MSTLAY. Note that the fluxes are represented in terms of
both quantity (eg. IFWI, in inches/interval) and as a fraction of the
contributing storage (eg. FII, as a fraction of UMST/interval)
The reader should also refer to Figure 4.2(1).3-2 in module section PWATER
when studying this diagram and the following discussion.
For the agri-chemical sections the moisture storages (the variables in
Figure 4.2(1).8-1 ending in MST) are calculated by the general equation:
MST = WSTOR + WFLUX (1)
The variable WSTOR is the related storage calculated in module section
PWATER (Figure 4.2(1).3-2). For example, in the calculation of the lower
layer moisture storage (LMST), WSTOR is the lower zone storage (LZS). The
variable WFLUX generally corresponds to the flux of moisture through the
soil layer. For the computation of LMST, WFLUX is the sum of water
percolating from the lower zone to the inactive (IGWI) and active
groundwater (AGWI) as determined in section PWATER. Note that these
equations are dimensionally non-homogeneous, because storages (inches) and
fluxes (inches/interval) are added together. Thus, the results given are
likely to be highly dependent on the simulation time step. The ARM Model,
from which the equations come, uses a step of 5 minutes. Extreme caution
should be exercised if the agricultural chemical sections (including MSTLAY)
are run with any other time step. For more details on the calculation of
the layer moisture storages, the reader should consult the pseudo code.
192
-------
surface
lay***
SMST
surface
layer
storqqe
f
ttrcptatinq }
moisture /
upper
layer
SURO
surface
outflow
F5O
UMST
upper layer
principal
storaqe
/UDOWN\
[percolatinq >
y moisture /
VFWI
motciure
qotnqto
transitory
storaqe
Fll
I3MST
upper layer
transitory
(interflow)
storaqe
IFWO
interflow
outflow
no
FUP
lower
layer
l&WI
moisture
LMST
lower
layer
storaqe
layer
FLDP
FLP
contain* the identifier
for the solute f luxe*
which are expressed
as -fraction of contri-
buting storaqes
S AGWI \
/moisture percolattnq
\ to active qtound-
N. water j \~
AMST
active
qroundw/ater
storage
FAO
Rqure
AfiWO
active
qround-
water
fluxes
of +he -hrar«spor+ of mois+ure I solutes
05 e&timcrfed in the MSTLAY sec+ion of ihe PERLND
Appl ication Mcdule .
193
-------
Module Section MSTLAY
The upper layer has been subdivided into two storages, principal and
transitory. The transitory (interflow) storage is used to transport
chemicals from the upper layer to interflow outflow. The chemicals in it do
not undergo any reactions. However, reactions do occur in the principal
storage.
The fluxes shown in Figure M.2(1).8-1 are the same as those in Figure
4.2(1).3-2 with the exceptions of SDOWN and UDOWN. SDOWN encompasses all
the water that moves downward from the surface layer storage. It is the
combination of the water infiltrating from the surface detention storage
directly to the lower zone (INFIL), the inflow to the upper zone (UZI), and
the water flowing into interflow storage (IFWT). UDOWN is all the water
percolating through the upper layer. It is INFIL plus the percolation from
the upper zone storage to the lower zone storage (PERC).
Each fractional solute flux is the appropriate moisture flux divided by the
contributing storage. For example, the fraction of chemical in solution
that is transported overland from the surface layer storage (FSO) is the
surface moisture outflow (SURO) divided by the surface layer moisture
storage (SMST).
The above estimates are based on the assumption that the concentration of
the solute being transported is the same as that in storage. They also
assume uniform flow through the layers and continuous mixing of the solutes.
However, these assumptions may need to be revised or implemented differently
for some of the transport. Past testing has shown that the above method
leads to excessive leaching of solutes (Donigian, et al. 1977). Factors that
retard solute leaching were added in the ARM Model Version II to remedy this
problem. For the surface layer, the percolation factor (SLMPF) affects the
solute fraction percolating (FSP) by the relationship:
FSP = SLMPF*SDOWN/SMST (2)
The variables SDOWN and SMST are defined in Figure 4.2(1).8-1. FSP will
typically be between 0 and 1.
For the upper or lower layer percolating fraction (FUP, FLDP, or FLP), the
retardation factor only has an influence when the ratio of the respective
zonal storage to the nominal storage times the factor (ZS/(ZSN*LPF)) is less
than one. The relationship under this condition is:
F = (ZS/(ZSN*LPF))*(PFLUX/MST) (3)
where:
F = layer solute percolating fraction
ZS = zonal moisture storage, either UZS or LZS
ZSN = zonal nominal moisture storage, either UZSN or LZSN
LPF = factor which retards solute leaching for the layer,
either ULPF or LLPF
PFLUX = percolation flux, either UDOWN, IGWI, or AGWI
MST = layer moisture storage, either UMST or LMST
194
-------
Module Section MSTLAY
4.2(1).9 Simulate Pesticide Behavior in Detail (Section PEST of Module
PERLND)
Purpose
Because of the complexity of pesticide behavior on the land, simulation of
the processes frequently requires considerable detail. Pesticide
applications vary in amount and time during the year. Various pesticides
adsorb and degrade differently. Some, like paraquat, attach themselves
strongly to the soil thereby appearing in low concentrations in water but in
high concentrations on soil particles. Others, like atrazine, undergo
complex interactions with the soil and are found in higher concentrations in
the runoff water than on the eroded sediment.
Section PEST models pesticide behavior by simulating the processes of
degradation and adsorption as well as transport. The pesticides are
simulated in the soil and runoff in three forms: dissolved, adsorbed, and
crystallized. These phases in the soil affect the forms and amounts in the
runoff.
Method
Pesticides are simulated by using the time series generated by other PERLND
module sections to transport and influence the adsorption and degradation
processes. Pesticides move with water flow or by association with the
sediment. They also may be adsorbed to the soil in varing degree as a
function of the chemical characteristics of the toxicant and the exchange
capacity of the soil layer. Pesticide degradation occurs to varying degrees
depending upon the susceptability of the compound to volatilization and
breakdown by light, heat, microorganisms and chemical processes. The
subroutines in module section PEST consider these transport and reaction
processes.
All the subroutines described in this module section except NONSV and DEGRAS
are accessed by other agri-chemical module sections because many of the
basic transport and reaction processes are similar. The subroutines are
described here because they are physically located in this subroutine group.
Subroutine AGRGET is first to be called. This subroutine has no computing
function; it obtains any required time series from the INPAD that is not
already available.
Subroutine SDFRAC determines the fraction of the surface layer soil that has
eroded. The amount eroded is the total sediment removed by scour and
washoff as determined in module section SEDMNT. The mass of soil in the
surface layer is a parameter value which does not vary even when material is
removed. The chemical which is associated with the sediment is assumed to
be removed from the surface layer storage in the same proportion that the
layer has eroded. Chemical removal is simulated in subroutine SEDMOV. A
sediment associated chemical is one that may be attached to the eroding soil
or one which may move with the soil. With pesticides the adsorbed form will
be attached to the soil particle, while the crystalline form will move with
195
-------
Module Section PEST
the soil particle being eroded but will not be attached to it. Both forms
are taken from their respective surface layer storages in proportion to the
fraction of the surface soil layer removed by overland flow.
Subroutines TOPMOV and SUBMOV perform a function similar to SEDMOV except
they move the solutes. Chemicals in solution move to and from the storages
according to the fractions calculated in section MSTLAY. Figure 4.2(1).9-1
schematically illustrates the fluxes and storages used in these subroutines.
The fractions (variables beginning with the letter "F") of the storages are
used to compute the solute fluxes. The equations used to compute the solute
transport fluxes from the fractions and storages are given in the figure.
Subroutine TOPMOV performs the calculations of the fluxes and the resulting
changes in storage for the topsoil layers (surface and upper), while SUBMOV
performs them for the subsurface layers (lower and active groundwater).
Biological and chemical reactions are performed on the pesticides (and other
chemicals) in each layer storage. Chemicals in the upper layer principal
storage undergo reactions while those in the transitory (interflow) storage
do not. The upper layer transitory storage is a temporary storage of
chemicals on their way to interflow outflow. Subroutine PSTRXN is called to
perform reactions on the pesticide in each layer.
4.2(1).9.5 Perform Reactions on Pesticides (subroutine PSTRXN)
Purpose
This code simulates the degradation and adsorbtion/desorption of pesticides.
This subroutine is called for each of the four soil layers and each
pesticide.
Method of Reacting Pesticides
The user has the option of adsorbing/desorbing the pesticide by one of three
methods. The first method is by first order kinetics. This method assumes
that the pesticide adsorbs and desorbs at a rate based on the amount in soil
solution and on the amount on the soil particle. It makes use of a
proportionality constant and is independent of the concentration. The
second method is by use of the single value Freundlich isotherm. This
method makes use of a single adsorption/desorption curve for determining the
concentration on the soil and in solution. The third method is by use of
multiple curves based on a varying Freundlich K value. Further details of
these methods can be found in the discussion of the individual subroutines
that follows and in the ARM Model reports (Donigian and Crawford 1976;
Donigian, et al. 1977).
Degradation is performed once a day by subroutine DEGRAS for each of the
four layers that contain pesticide. The amount degraded is determined
simply by multiplying a decay rate parameter specified for each soil layer
by each of the three forms (adsorbed, solution, and crystalline) of
pesticide in storage. The degraded amounts are then subtracted from their
respective storages. This method of simulating degradation lumps complex
processes in a simple parameter.
196
-------
SSCM
dissolved
chemical in
surface
storaqe
i
a
USCM
dissolved chenx
in upper layer
principal
storaqe
USGM*FUP* UPCM
/percobtion
to lower
layer
storaqe i
LSCM
chemical in
lower layer
storaqe
it
FLDP-LOPCM
/pereoUitionio
inocii^ aroun
SOCKA»SSttA*F50
chemical
in
surface
outflow
15CM*FIO*
transfer
from prin-
cipal to
transitory
storaqe
ISCM
dissolved diem.
in upper layer
transitory
(interflow)
storage
chemical
in
interflow
outflow
LPCM*L5CM*FLP
percolation
to active
qround-
water
ASCM
ASCM^FAO
chemical
in active
qround-
water
in active
qroundwater
Fiqure 4.10 ").
-------
Module Section PEST
4.2(1).9.5.1 Adsorb/Desorb Using First Order Kinetics (subroutine FIRORD)
Purpose
The purpose of this subroutine is to calculate the adsorption and desorption
reaction fluxes of chemicals using temperature dependent first order
kinetics. These fluxes are calculated every simulation interval when the
subroutine is called by section PEST, but they are determined only at the
designated chemical reaction frequency when called by sections NITR and
PHOS.
Method
The calculation of adsorption and desorption reaction fluxes by first order
kinetics for soil layer temperatures less than 35 degrees C takes the form:
DES = CMAD*KDS*THKDS**(TMP-35.0) (1)
ADS = CMSU*KAD*THKAD**(TMP-35.0) (2)
where:
DES = current desorption flux of chemical in mass/area per interval
CMAD = storage of adsorbed chemical in mass/area
KDS = first order desorption rate parameter, per interval
THKDS = temperature correction parameter for desorption
TMP = soil layer temperature in degrees C
ADS = current adsorption flux of chemical in mass/area per interval
CMSU = storage of chemical in solution in mass/area
KAD = first order adsorption rate parameter, per interval
THKAD r temperature correction parameter for adsorption
THKDS and THKAD are typically about 1.06
All of the variables except the temperature coefficients may vary with the
layer of the soil being simulated. The soil temperatures are time series
which may be input (eg. using field data) or simulated in module section
PSTEMP. The temperature correction of the reaction rate parameter is based
on the Arrhenius equation. At temperatures of 35 degrees C or above no
correction is made. When the temperature is at 0 degrees C or below or the
soil layer is dry, no adsorption and desorption occurs.
The storage of the solution chemical is updated every simulation interval in
the calling subroutine, that is, in PSTRXN, NITRXN, or PHORXN, by adding DES
minus ADS. Likewise, the storage of the adsorbed chemical is updated there
also by adding ADS minus DES. An adjustment is made in the calling
subroutine, if any of the fluxes would cause a storage to go negative. When
this happens a warning message is produced and fluxes are adjusted so that
no storage goes megative. This usually occurs when large time steps are used
in conjuction with large KAD and KDS values.
198
-------
Module Section PEST
4.2(1).9.5.2 Adsorb/Desorb Using the Single Value Freundlich Method
(subroutine SV)
Purpose
Subroutine SV calculates the adsorption and desorption and the resulting new
storages of a chemical using the single value Freundlich method.
Method
The Freundlich isotherm methods, unlike first order kinetics, assume
instantaneous equilibrium. That is, no matter how much chemical is added to
a particular phase, equilibrium is assumed to be established between the
solution and adsorbed phase of the chemical. These methods also assume that
for any given amount of chemical in the soil, the equilibrium distribution
of the chemical between the soil solution and on the soil particle can be
found from an isotherm. Figure 4.2(1).9-2 illustrates such an isotherm.
Three phases of the chemical are actually possible; crystalline, adsorbed,
and solution. The crystalline form is assumed to occur only when the soil
layer is dry, or when there is more chemical in the layer than the combined
capacity to adsorb and hold in solution. When the soil is dry, all the
chemical is considered to be crystalline salt. When there is more total
chemical in the soil layer than the the soil adsoption sites can contain and
more than that saturated in solution, then the chemical content which
exceeds these capacities is considered to be crystalline salt. Module
section PEST considers crystalline phase storage, but in module sections
NITR and PHOS this is not so. Instead, any crystalline phosphate or
ammonium predicted by an isotherm is added to the adsorbed phase storage.
The adsorbed and solution phases of the chemical are determined in this
subroutine by the standard Freundlich equation as plotted by curve 1 in
Figure 4.2(1).9-2. When the amount of chemical is less than the capacity of
the soil particle lattice to permanently bind the chemical (XFIX), then all
the material is consider fixed. All the fixed chemical is contained in the
adsorbed phase of the layer storage. Otherwise, the Freundlich equation for
curve 1 is used to determine the partitioning of the chemical into the
adsorbed and solution phases:
X = KF1*C**(1/N1) + XFIX (3)
where:
X = chemical adsorbed on soil, in ppm of soil
KF1 = single value Freundlich K coefficient
C = equilibrium chemical concentration in solution,
in ppm of solution
N1 = single value Freundlich exponent
XFIX = chemical which is permanently fixed, in ppm of soil
The above equation is solved in subroutine ITER by an iteration technique.
The parameters used in the computation can differ for each layer of the
soil.
199
-------
XMAX
XJCT
X,ppm
XDIF
XFIX
1
C,ppm
I
CMAY
Rqure 4.
.<)-2 Freundlich Isotherm Calculations. X ift the
chemical adsorbed on the soil and C is the
chemical in solution.
200
-------
Module Section PEST
4.2(1);9.5.3 Msorb/Desorb Using the Non-single Value Freundlich Method
(subroutine NONSV)
Purpose
The purpose of this subroutine is to calculate the adsorption/desorption of
a chemical by the nonsingle value Freundlich method. The single value
Freundlich method was found to inadequately represent the division of some
pesticides between the soil particle and solution phases, so this method was
developed as an option in the ARM Model (Donigian and Crawford 1976). This
subroutine is only available for use by the PEST module section.
Method
The approach in this code uses the same algorithms and solution technique as
subroutine SV for determining curve 1 in Figure M.2(1).9-2. However, curve
1 is used solely for adsorption. That is, only when the concentration of the
adsorbed chemical is increasing. When desorption occurs a new curve (curve
2) is used:
X = KF2*C**(1/N2) + XFIX (4)
KF2 = (KF1/XDIF)**(N1/N2) * XDIF (5)
where:
KF2 = nonsingle value Freundlich coefficient
N2 = nonsingle value Freundlich exponent parameter
XDIF = XJCT - XFIX
XJCT = the adsorbed concentration where curve 1 joins curve 2
(ie. where desorption started)
as shown in Figure 4.2(1).9-2, in ppm of soil
The other variables are as defined for subroutine SV.
Once curve 2 is used, both desorption and adsorption follow it until the
adsorbed concentration is less than or equal to XFIX or until it reaches
XJCT. Then, adsorption will again take place following curve 1 until
desorption reoccurs, following a newly calculated curve 2. The solution of
the Freundlich equations for curves 1 and 2 utilizes the same iteration
technique introduced in subroutine SV (subroutine ITER).
201
-------
Module Section NITR
4.2(1).10 Simulate Nitrogen Behavior in Detail (section NITR of module
PERLND)
Purpose
NITR, like section PEST, simulates the behavior of chemicals in the soil
profile of a land segment. Section NITR handles the nitrogen species of
nitrate, ammonia, and organic nitrogen. This involves simulating nitrogen
transport and soil reactions. Nitrogen, like phosphorus, may be a limiting
nutrient in the eutrophication process in lakes and streams. Nitrates in
high concentrations may also pose a health hazard to infants.
Method of Simulating Nitrogen
Nitrogen species are transported by the same methods used for the pesticide
forms. The subroutines that are called to transport nitrogen are located in
and described with the PEST module section. Organic nitrogen and adsorbed
ammonium are removed from the surface layer storage by association with
sediment (subroutines SDFRAC and SEDMOV). Nitrate and ammonium in the soil
water are transported using the simulation subroutines TOPMOV and SUBMOV.
Nitrogen reactions are simulated separately for each of the soil layers. The
method is discussed below.
4.2(1).10.1 Perform Reactions on Nitrogen Forms (subroutine NITRXN)
Purpose
The purpose of NITRXN is to simulate soil nitrogen transformations. This
includes plant uptake of nitrate and ammonium, denitrification or reduction
of nitrate-nitrite, immobilization of nitrate-nitrite and ammonium,
mineralization of organic nitrogen, and the adsorption/desorption of
ammonium.
Method of Nitrogen Transformations
Nitrogen reactions can be divided between those that are chemical in nature
and those that are a combination of chemical and biological reactions. The
adsorption and desorption of ammonium is a chemical process. The user has
the option of simulating ammonium adsorption and desorption by first order
kinetics with subroutine FIRORD or by the Freundlich isotherm method with
subroutine SV. These subroutines are described in the PEST module section.
The user has the option of specifying how often the adsorption and
desorption rates are calculated. When adsorption/desorption is simulated by
the Freundlich method, the solution and adsorbed storages of ammonium are
determined instantaneously at the specified frequency of reaction. However,
when the first order method is used, the temperature corrected reaction
202
-------
Module Section NITR
fluxes (Figure 4.2(1).10-1) are recomputed intermittently, but the.storages
are updated every simulation interval.
The other reactions are a combination of biological and chemical
transformations. They are accomplished by first order kinetics only. The
optimum first order kinetic rate parameter is corrected for soil
temperatures below 35 degrees C by the generalized equation:
KK = K*TH**(TMP-35.0) (1)
where:
KK = temperature corrected first order transformation rate
in units of per simulation interval
K = optimum first order reaction rate parameter
TH = temperature coefficient for reaction rate correction
(typically about 1.06)
TMP = soil layer temperature in degrees C
When temperatures are greater than 35 degrees C, the rate is considered
optimum, that is, KK is set equal to K. When the temperature of the soil
layer is below 4 degrees C or the layer is dry, no biochemical
transformations occur.
Identifiers with a leading "K" (eg. KDNI) are the optimum rates; those for
corrected rates have both a leading and trailing "K" (eg. KDNIK).
The corrected reaction rate parameters are determined every biochemical
reaction interval and multiplied by the respective storages as shown in
Figure 4.2(1).10-1 to obtain the reaction fluxes. Plant uptake can vary
monthly and can be distributed between nitrate and ammonium by the
parameters N03UTF and NH4UTF. These parameters are intended to designate
the fraction of plant uptake from each-species of N; the sum of N03UTF and
NH4UTF should be 1.0.
The first order reaction rate fluxes that are shown in Figure 4.2(1).10-1
are coupled, that is, added to and subtracted from the storages
simultaneously. The coupling of the fluxes is efficient in use of computer
time but has a tendency to produce unrealistic negative storages when large
reaction intervals and large reaction rates are used jointly. A method has
been introduced which will modify the reaction fluxes so that they do not
produce negative storages. A warning message is issued when this
modification occurs.
4.2(1).11 Simulate Phosphorous Behavior in Detail (Section PHOS of Module
PERLND)
Purpose
Module section PHOS simulates the behavior of phosphorus in a pervious land
segment. This involves modeling the transport, plant uptake,
adsorption/desorption, immobilization, and mineralization of the various
forms of phosphorus. Because phosphorus is readily tied to soil and
203
-------
harvettinq
to
atmosphere
M03*KDNIK- DEMI
(derutrifi-\
cation J
PLTN
plant
nitroqen
H103UTF*W03*VCPLMK=UTNI
AMSU*KNIK
N03
nitrate
(plus
rite)
(pl
nit
I nitrification
/immobiU-
( zationot
\ nitrate
AMSU^KIMAMK= I MM AM
immobi-
lization
of
ammonium
ORGN
orqanic
nitroq«n
ammo-
nification
of
orqanic
nitroqen
I plant
uptake of
nrtrate
S plant \
/ uptake of I
\ ammonium /
AM5U
ammonium
in
solution
DESAM=AMAD* KDSAMK
or
sinqle value
Freundlich
meihod
(instaneous)
ammonium
desorphon
AMAD
ammonium
adsorbed
ammonium
adsorption
ADSAM=AMSU* KADAMK
or
»inqle value
Freundlich
method
Rqure 4.ZCO.IO-J Flow diaqram for nitroqen reactions
204
-------
Module Section PHOS
sediment, it is usually scarce in streams and lakes. In fact, in many cases
it is the limiting nutrient in the eutrophication process. Because of its
scarcity accurate simulation is particularly important.
Method of Simulating Phosphorus
The method used to transport and react phosphorus is the same as that used
for nitrogen in module section NITR. The subroutines used to transport
phosphorus are described in module section PEST. Organic phosphorus and
adsorbed phosphate are removed on or with sediment by calling subroutine
SEDMOV. Phosphate in solution is transported in the moving water using
subroutines TOPMOV and SUBMOV. Phosphorus reaction is simulated in the soil
by subroutine PHORX.N.
In subroutine PHORXN, phosphate is adsorbed and desorbed by either first
order kinetics or by the Freundlich method. The mechanics of these methods
are described in module section PEST. As with the simulation of ammonium
adsorption/desorption, the frequency of this chemical reaction for phosphate
can also be specified. Unlike ammonium, typically phosphate includes a
large portion which is not attached to the soil particle but is combined
with cations. This is because phosphate is much less soluble with the ions
found in soils than ammonium.
Other reactions performed by subroutine PHORXN include mineralization,
immobilization, and plant uptake. These are accomplished using temperature
dependent first order kinetics; the same method used for the nitrogen
reactions. The general description of this process is in module section
NITR. Figure U.2(1).11-1 shows the parameters and equations used to
calculate the reaction fluxes for phosphorus. Reactions are simulated for
each of the four soil layers using separate parameter sets for each layer.
As with nitrogen, the biochemical phosphate reaction fluxes of
mineralization, immobilization, and plant uptake can be determined at an
interval less frequent than the basic simulation interval.
I.2(1).12 Simulate Movement of a Tracer (Section TRACER of Module PERLND)
Purpose
The purpose of this code is to simulate the movement of any nonreactive
tracer (conservative) in a pervious land segment. Chloride, bromide, and
dyes are commonly used tracers which can be simulated by section TRACER.
Also, total dissolved salts could possibly be modeled by this section.
Typically, this code is applied to chloride to calibrate solute movement
through the soil profile. This involves adjustment of the percolation
retardation factors (see section MSTLAY) until good agreement with observed
chloride concentrations has been obtained. Once these factors have been
calibrated, they are used to simulate the transport of other solutes, such
as nitrate.
Method of Simulating Tracer Transport
Tracer simulation uses the agri-chemical solute transport subroutines TOPMOV
and SUBMOV which are described in section PEST. No reactions are modeled.
205
-------
PLTP
plant
phosphorus
UTP4*P4SU*W>U>K.
1MMP4*M5U*KIMPK.
phosphate
immobili-
zation
ORGP
organic
phosphorus
/ plant
uptake of
\ phosphorus
draxption
of
phosphate
P4-SU
phosphate
in
solution
Or
sinqle value
Freundlich
method
organic
phosphonu
mineral-
iiation
P4AD
phosphate
adsorbed
adsorption
of
phosphate
AD5P4= P45U* KADPiC
or
single value
Freundlich method
(instantaneous)
Rqure 4.i(l).ll-l Flow diagram for phosphorus reactions
206
-------
Module IMPLND
4.2(2) Simulate an Impervious Land Segment (Module IMPLND)
In an impervious land segment, little or no infiltration occurs. However,
land surface processes do occur as illustrated in Figure 4.2(2)-1. Snow may
accumulate and melt, and water may be stored or may evaporate. Various
water quality constituents accumulate and are removed. Water, solids, and
various pollutants flow from the segments by moving laterally to a downslope
segment or to a reach/reservoir.
Module IMPLND simulates these processes. The sections of IMPLND and their
functions are given in Structure Chart 4.2(2) (Part D). They are executed
from left to right. Many of them are similar to the corresponding sections
in the PERLND module. In fact, since sections SNOW and ATEMP perform
functions that can be applied to pervious or impervious segments, they are
shared by both modules. IWATER is analogous to PWATER in module PERLND;
SOLIDS is analogous to SEDMNT; IWTGAS is analogous to PWTGAS; and IQUAL is
analogous to PQUAL. However, the IMPLND sections are simpler since they
contain no infiltration function and consequently no subsurface flows.
IPTOT, IBAROT, and IPRINT service the IMPLND module similarly to the
corresponding code in PERLND.
4.2(2).3 Simulate the Water Budget for an Impervious Land Segment
(Section IWATER of Module IMPLND)
Purpose
Section IWATER simulates the retention, routing, and evaporation of water
from an impervious land segment.
Method
Section IWATER is similar to section PWATER of the PERLND module. However,
IWATER is simpler because there is no infiltration and consequently no
subsurface processes. IWATER is composed of the parent subroutine plus
three subordinate subroutines: RETN, IROUTE, and EVRETN. RETN is analogous
to ICEPT, IROUTE is analogous to PROUTE, and EVRETN is analogous to EVICEP
in module section PWATER. The time series requirements are
the same as for section PWATER.
Figure 4.2(2).3-1 schematically represents the fluxes and storages simulated
in module section IWATER. Moisture (SUPY) is supplied by precipitation, or
under snow conditions, it is supplied by the rain not falling on the
snowpack plus the water yielded by the snowpack. This moisture is available
for retention; subroutine RETN performs the retention functions. Lateral
surface inflow (SURLI) may also be retained if the user so specifies by
setting the flag parameter RTLIFG=1. Otherwise, retention inflow (RETI)
207
-------
spswsir
\ \\ Runoff, Solids, Water
\ sy \ Quality Constituents
Figure 4.2(2)-l Inpervious land segment processes
208
-------
Module Section IWATER
equals SUPY. Moisture exceeding the retention capacity overflows the storage
and is available for runoff.
The retention capacity, defined by the parameter RETSC, can be used to
designate any retention of moisture which does not reach the overland flow
plane. RETSC may be used to represent roof top catchments, asphalt wetting,
urban vegetation, improper drainage, or any other containment of
water that will never flow from the land segment. The user may supply the
retention capacity on a monthly basis to account for seasonal variations, or
may supply one value designating a fixed capacity.
Water held in retention storage is removed by evaporation (IMPEV). The
amount evaporated is determined in subroutine EVRETN. Potential
evaporation is an input time series.
Retention outflow (RETO) is combined with any lateral inflow when RTLIFG=0
producing the total inflow to the detention storage (SURD. Water remaining
in the detention storage plus any inflow is considered the moisture supply.
The moisture supply is routed from the land surface in subroutine IROUTE.
4.2(2).3.2 Determine How Much of the Moisture Supply Runs Off
(subroutine IROUTE)
Purpose
The purpose of subroutine IROUTE is to determine how much of the moisture
supply runs off the impervious surface in one simulation interval.
Method of Routing
A method similar to that used in module PERLND (Section 4.2(1).3.2.1.3) is
employed to route overland flow.
4.2(2).4 Simulate Accumulation and Removal of Solids (Section SOLIDS of
Module IMPLND)
Purpose
Module section SOLIDS simulates the accumulation and removal of solids by
runoff and other means from the impervious land segment. The solids outflow
may be used in section IQUAL to simulate quality constituents associated
with particulates.
210
-------
Module Section SOLIDS
Method
The equations used in this section are based on those in the NFS Model
(Donigian and Crawford 1976b). Figure 4.2(2).4-1 schmatically represents
the fluxes and storages simulated by section SOLIDS. Lateral input of solids
by water flow is a user designated option which is unique to HSPF. Washoff
of solids may be simulated by one of two ways. One subroutine is similar to
the method used in the NFS Model. However, this method is dimensionally
nonhomogeneous. That is, a flux and a storage are added making the answer
more interval dependent. This technique has only been used with 15-min
intervals. The other subroutine is dimensionally homogenous, since only a
flux term is used in the solution. However, this method has not been
tested!
The accumulation and removal of solids which occurs independently of runoff
(eg. by atmospheric fallout, street cleaning) is handled in subroutine
ACCUM.
1.2(2).4.1 Washoff Solids Using Method 1
(subroutine SOSLD1)
Purpose
Subroutines SOSLD1 and SOSLD2 perform the same task but by different
methods. They simulate the washoff of solids from an impervious land
segment.
Method
When simulating the washoff of solids, the transport capacity of the
overland flow is estimated and compared to the amount of solids available.
The transport capacity is calculated by the equation:
STCAP = DELT60*KEIM*((SURS + SURO)/DELT60)**JEIM (1)
where:
STCAP = capacity for removing solids in
tons/acre per interval
DELT60 = hours per interval
KEIM r coefficient for transport of solids
SUES = surface water storage in inches
SURO = surface outflow of water in in./interval
JEIM = exponent for transport of solids
When STCAP is greater than the amount of solids in storage, washoff is
calculated by:
SOSLD = SLDS*SURO/(SURS + SURO) (2)
If the storage is sufficient to fulfill the transport capacity, then the
following relationship is used:
211
-------
nentoval N^ |
an
wind,
SLSLD
input" of
solids
'ACZ5DP \
accumubcfion)
J
SLDS
F igure 4. z (z). 4 -J Plow diagram c?f -the Sou DS
eeciion of the XMH.ND application
module.
212
-------
Module Section SOLIDS
SOSLD = STCAP*SURO/(SURS + SURO) (3)
where:
SOSLD = washoff of solids in tons/acre per interval
SLDS = solids storage in tons/acre
SOSLD is then subtracted from SLDS.
Subroutine SOSLD1 differs from SOSLD2 in that it uses the dimensionally
nonhomogeneous term (SURS + SURO)/DELT60 in the above equations, while
SOSLD2 uses the homogeneous term SURO/DELT60.
4.2(2).4.2 Washoff Solids Using Method 2
(subroutine SOSLD2)
Purpose
The purpose of this subroutine is the same as SOSLD1. They only differ in
method.
Method of Determining Removal
This method of determining sediment removal has not been tested. Unlike
subroutine SOSLD1, it makes use of the dimensionally homogeneous term
SURO/DELT60 instead of (SURO+SURS)/DELT60 in the following equation.
STCAP = DELT60*KEIM*(SURO/DELT60)**JEIM (4)
When STCAP is more than the amount of solids in storage, the flow washes off
all of the solids storage (SLDS). However, when STCAP is less than the
amount of solids in storage, the situation is transport limiting, so SOSLD
is equal to STCAP.
4.2(2).4.3 Accumulate and Remove Solids Independently of Runoff
(subroutine ACCUM)
Purpose
Subroutine ACCUM simulates the accumulation and removal of solids
independently of runoff; for example, atmospheric fallout and street
cleaning.
Method
The storage is updated once a day, on those days when precipitation did not
occur during the previous day, using the equation:
213
-------
Module Section SOLIDS
SLDS = ACCSDP + SLDSS*(1.0 - REMSDP) (5)
where:
ACCSDP = accumulation rate of the solids storage,
Ibs/acre per day
SLDS = solids in storage at end of interval in tons/acre
SLDSS = solids in storage at start of interval in tons/acre
REMSDP = unit removal rate of solids in storage
(i.e. fraction removed per day)
ACCSDP and REMSDP may be input on a monthly basis to account for seasonal
variations.
Note that, if no runoff occurs, equation 5 will cause the solids storage to
asymptotically approach a limiting value. The limit, found by setting SLDS
and SLDSS to the same value (SLDSL), is:
SLDSL = ACCSDP/REMSDP (6)
4.2(2).5 Estimate Water Temperature and Dissolved Gas Concentrations
(Section IWTGAS of Module IMPLND)
Purpose
IWTGAS estimates the water temperature and concentrations of dissolved
oxygen and carbon dioxide in the outflow from the impervious land segment.
Method
Outflow temperature is estimated by the following regression equation:
SOTMP = AWTF + BWTF*AIRTC (1)
where:'
SOTMP = impervious surface runoff temperature in degrees C
AWTF = Y-intercept
BWTF = slope
AIRTC = air temperature in degrees C
The parameters AWTF and BWTF may be input on a monthly basis. When snowmelt
contributes to the outflow, SOTMP is set equal to 0.5.
The dissolved oxygen and carbon dioxide concentrations of the overland flow
are assumed to be at saturation and are calculated as direct functions of
water temperature. IWTGAS uses the following empirical nonlinear equation
to relate dissolved oxygen at saturation to water temperature (Committee on
Sanitary Engineering Research 1960):
SODOX = (14.652 + SOTMP*(-0.41022 +
SOTMP*(0.007991 - 0.000077774*SOTMP)))*ELEVGC (2)
where:
214
-------
Module Section IWTGAS
SODOX = concentration of dissolved oxygen in surface outflow in mg/1
SOTMP = surface outflow temperature in degrees C
ELEVGC = correction factor for elevation above sea level
(ELEVGC is calculated by the Run Interpreter dependent
upon mean elevation of the segment)
The empirical equation for dissolved carbon dioxide concentration of the
overland flow (Harnard and Davis 1943) is:
SOC02 = (10**(2385.73/ABSTMP - 14.0184 + 0.0152642*ABSTMP))
*0.000316*ELEVGC*12000.0 (3)
where:
SOC02 = concentration of dissolved carbon dioxide in
surface outflow in rag C/l
ABSTMP = absolute temperature of surface outflow in degrees K
4.2(2).6 Simulate Washoff of Quality Constituents Using Simple Relationships
with Solids and Water Yield (Section IQUAL of Module IMPLND)
Purpose
The IQUAL module section simulates water quality constituents or pollutants
in the outflows from an impervious land segment using simple relationships
with water yield and/or solids. Any constituent can be simulated by this
module section. The user supplies the name, units and parameter values
appropriate to each of the constituents that he wishes to simulate. Note
that more detailed methods of simulating solids, heat, dissolved oxygen, and
dissolved carbon dioxide removal from the impervious land segment are
available in other sections of IMPLND.
Approach
The basic algorithms used to simulate quality constituents are a synthesis
of those used in the NPS Model (Donigian and Crawford 1976b) and HSP QUALITY
(Hydrocomp 1977). However, some options and combinations are unique to HSPF.
Figure 4.2(2).6-1 shows schematically the fluxes and storages represented in
module section IQUAL. A quality constituent may be simulated by two methods.
One approach is to simulate the constituent by association with solids
removal. The other approach is to simulate it by using basic accumulation
and depletion rates together with depletion by washoff; that is, constituent
outflow from the surface is a function of the water flow and the constituent
in storage. A combination of the two methods may be used in which the
individual outfluxes are added to obtain the total surface outflow. These
approaches will be discussed further in the descriptions of the
corresponding subroutines.
215
-------
Module Section IQU.AL
IQUAL allows the user to simulate up to 10 quality constituents at a time.
If a constituent is considered to be associated with solids, it is called a
QUALSD. The corresponding term for constituents associated directly with
overland flow is QUALOF. Each of the 10 constituents may be defined as
either a QUALSD or a QUALOF or both. However, no more than seven of any one
of the constituent types (QUALSD or QUALOF) may be simulated in one
operation. The program uses a set of flag pointers to keep track of these
associations. For example, QSDFP(3)=0 means that the third constituent is
not associated with solids, whereas QSDFP(6)=U means that the sixth
constituent is the fourth solids associated constituent (QUALSD). Similar
flag pointer arrays are used to indicate whether or not a quality
constituent is a QUALOF.
A
removal/X.
by cleaninqA
decoy, vwind.y
\/
ujccumulation
sao
storage of QUAL
on surface iof
direct washoff
by cuertand
storage of
QUAL
associated
\Nith solid 5
direct
wa*hof f of
QUM by
overland
flow
wcrthoff
of QUAL
associated
with
soiids
SOQU/kL
total
wa*ho
-------
Module Section IQUAL
4.2(2).6.1 Remove by Association with Solids
(subroutine WASHSD)
Purpose
WASHSD simulates the removal of a quality constituent from the impervious
land surface by association with the solids removal determined in module
section SOLIDS.
Method
This approach assumes that the particular quality constituent removed from
the land surface is in proportion to the solids removal. The relation is
specified by user-input "potency factors." Potency factors indicate the
constituent strength relative to the solids removal from the surface. For
each quality constituent associated with solids, the user supplies separate
potency factors. The user is also able to supply monthly potency factors for
constituents that vary somewhat consistently throughout the year;
Removal of the solids associated constituent by solids washoff is simulated
by:
SOQS = SOSLD*POTFW (1)
where:
SOQS = flux of quality constituent associated with
solids washoff in quantity/acre per interval
SOSLD = washoff of detached solids in tons/acre per interval
POTFW = washoff potency factor in quantity/ton
The unit "quantity" refers to mass units (pounds or tons in the English
system) or some other quantity, such as number of organisms for coliforms.
The user specifies the units of "quantity".
4.2(2).6.2 Accumulate and Remove by a Constant Unit Rate and'by Overland Flow
(subroutine WASHOF)
Purpose
WASHOF simulates the accumulation of a quality constituent on the impervious
land surface and its removal by a constant unit rate and by overland flow.
Method
This subroutine differs from subroutine WASHSD in that the storage of the
quality constituent is simulated. The stored constituent can be accumulated
and removed by processes which are independent of storm events, such as
cleaning, decay, and wind deposition, and it is washed off by overland flow.
217
-------
Module Section IQUAL
The accumulation and removal rates can have monthly values to account for
seasonal fluctuations.
When there is surface outflow and some quality constituent is in storage
then washoff is simulated using the commonly used relationship:
SOQO = SQO*(1.0 - EXP(-SURO*WSFAO) (2)
where:
SOQO = washoff of the quality constituent from the land
surface in quantity/acre per interval
SQO = storage of the quality constituent on the surface
in quantity/acre
SURO = surface outflow of water in in./interval
WSFAC = susceptibility of the quality constituent to washoff
in units of 1/inch
EXP = Fortran exponential function
The storage is updated once a day to account for accumulation and removal
which occurs independent of runoff by the equation:
SQO = ACQOP + SQOS*(1.0 - REMQOP) (3)
where:
ACQOP = accumulation rate of the constituent,
quantity/acre per day
SQOS = SQO at the start of the interval
REMQOP = unit removal rate of the stored constituent,
per day
The Run Interpreter computes REMQOP and WSFAC for this subroutine according
to:
REMQOP = ACQOP/SQOLIM (4)
where:
SQOLIM = asymptotic limit for SQO as time approaches
infinity, (quantity/acre), if no washoff occurs
and
WSFAC = 2.30/WSQOP (5)
where:
WSQOP = rate of surface runoff which results in a 90 percent
washoff in one hour, in./hr
Since the unit removal rate (REMQOP) is computed from two other parameters,
it is not supplied directly by the user.
218
-------
Module RCHRES
4.2(3) Simulate a Free-flowing Reach or Mixed Reservoir
(Module RCHRES)
This module simulates the processes which occur in a single reach of open or
closed channel or a completely mixed lake. For convenience such a
processing unit is referred to as a RCHRES throughout this documentation.
In keeping with the assumption of complete mixing, the RCHRES consists of a
single zone situated between two nodes, which are the extremities of the
RCHRES.
Flow through a RCHRES is assumed to be unidirectional. The inflow and
outflow of materials through a RCHRES are illustrated in Figure 4.2(3)-1.
Water and other constituents which arrive from other RCHRES's and local
sources enter the RCHRES through a single gate (INFLO). Outflows may leave
the RCHRES through one of several gates (OFLO). A RCHRES can have up to 5
OFLO gates. Precipitation, evaporation, and other fluxes also influence the
processes which occur in the RCHRES but do not pass through the gates.
The ten major subdivisions of the RCHRES module and their functions are
presented in Structure Chart 4.2(3). RPTOT, RBAROT, and RPRINT perform the
storage and printout of results from the other module sections of RCHRES
(HYDR through RQUAL). Within a module section, simulation of physical
processes (longitudinal advection, sinking, benthal release) is always
performed before simulation of biochemical processes.
The user specifies which-module sections are active. If any "quality"
sections (CONS through RQUAL) are active, section ADCALC must also be
active; it computes certain quantities needed to simulate advection of the
quality constituents. Besides fulfilling this requirement, the user must
ensure that all the time series required by the active sections are
available, either as supplied input time series or as data computed by
another module section. For example, if RQUAL is active, the water
temperature must be supplied, either as an input time series or by
activating section HTRCH which will compute it.
4.2(3).1 Simulate Hydraulic Behavior (Section HYDR of Module RCHRES)
Purpose
The purpose of this code is to simulate the hydraulic processes occurring in
a reach or a mixed reservoir (RCHRES). The final goal of the process may be
to route floods, study reservoir behavior, or analyze constituents dissolved
in the water.
219
-------
gate OFLO 1
OVOL-water
point discharges
NJ
O
gate IN FLO
water -IVOL
votft/e- ICON
etc.
torage withinthe RCHRES
water-VOL
con5erv