EPA 910/8-79-106	United States	Region 10
Environmental Protection	1200 Sixth Avenue
Agency	Seattle, WA 98101
Surveillance & Analysis Division	
<&EFft The Application of a
Dynamic
Estuary Water Quality Model
To Grays Harbor, Washington

-------
EPA 910/8-79-106
April 1980
THE APPLICATION OF A DYNAMIC ESTUARY
WATER QUALITY MODEL
TO
GRAYS HARBOR, WASHINGTON
Prepared by
John R. Yearsley
8ruce R. Cleland
EPA--Region 10
1200 Sixth Avenue
Seattle Washington 98101

-------
THIS DOCUMENT IS AVAILABLE IN LIMITED QUANTITIES
THROUGH THE U.S. ENVIRONMENTAL PROTECTION AGENCY,
SURVEILLANCE AND ANALYSIS DIVISION, 1200 SIXTH
AVENUE, SEATTLE, WASHINGTON 98101

-------
A WORKING PAPER PRESENTS RESULTS OF
INVESTIGATIONS WHICH ARE, TO SOME EXTENT,
LIMITED OR INCOWLETE. THEREFORE,
CONCLUSIONS OR RECOMMENDATIONS 	
EXPRESSED OR IMPLIED -- MAY BE TENTATIVE.

-------
TABLE OF CONTENTS
CHAPTER	PAGE
FINDINGS AND CONCLUSIONS	 1
INTRODUCTION	 5
PREVIOUS INVESTIGATIONS	 9
ESTUARY HYDRODYNAMIC MODEL	 11
Theoretical Considerations	 11
Model Implementation	 15
ESTUARY ECOLOGIC MODEL	 37
RESULTS	 45
LIST OF REFERENCES	 59

-------
LIST OF FIGURES
FIGURE	PAGE
1.	Grays Harbor and freshwater tributaries	6
2.	Junction configuration for Dynamic Estuary Model	13
3.	Comparison fo tidal stages at Montesano	23
4.	Comparison of simulated and observed tidal
stages at Port Dock			24
5.	Comparison of simulated and observed salinities
in Grays Harbor	52
6.	Comparison of simulated and observed salinities
in Grays Harbor	53
7.	Comparison of simulated and observed temperatures
near Light 15	54
8.	Comparison of simulated and observed temperatures
near UPRR Bridge	55
9.	Comparison of simulated and observed dissolved
oxygen in Grays Harbor	*	56
10. Comparison of simulated and observed dissolved
oxygen in Grays Harbor	57

-------
LIST OF TABLES
TABLE	PAGE
1.	Difference between National Geodetic Vertical Datum
and mean lower low water in Grays Harbor			16
2.	Time/stage relationship between tide gauge reading at
Port Dock and other locations in Grays Harbor	21
3.	Freshwater discharges in Grays Harbor (July 25-29, 1977)	31
4.	Major waste discharges in Grays Harbor (July 25-29, 1977)	32
5.	Calculations of Manning's "n" to account for energy
losses due to river bends	34
6.	Water quality parameter Interactions in Dynamic Estuary Model..39

-------
FINDINGS AND CONCLUSIONS
Water quality data collected in Grays Harbor during the period
July 25-29, 1977, a period of low flow throughout the Pacific Northwest,
was used to evaluate a dynamic estuary water quality model. The model,
developed by Water Resources Engineers, Inc. (1975) simulates estuarine
hydrodynamics, based upon tidal boundary conditions and freshwater
inflow. The hydrodynamic simulation is used, in conjunction with
relations defining other physical, chemical and biological processes, to
describe the mass balance equations for as many as 22 water quality
parameters. In this study, due to resource limitations, the only
parameters simulated were temperature, salinity dissolved oxygen and
biochemical oxygen demand.
The evaluation of the hydrodynamic model was based upon tidal height
data collected at the Port Dock in Aberdeen, Washington by the U.S. Army
Corps of Engineers. No velocity data was collected during the period
July 25-29, 1977. It was, therefore, necessary to rely upon estimates of
velocities from field studies done by Beverage and Swecker (1969) and
hydraulic model studies done by the U.S. Army Corps of Engineers (Brogdon
(1972)) for evaluating velocities simulated by the model.
Estimates of slack tide times at various locations in Grays Harbor
(National Ocean Survey (1975)), were also used to analyze the adequacy of
the hydrodynamic simulations, as well as to adjust observed tidal times
and heights at Aberdeen, Washington to tidal times and heights at the
Point Chehalis, the entrance to Grays Harbor.
For the resulting hydrodynamic simulations the range of difference
between published data and simulated values was jjjO.5 for tidal stages
(feet), +15i6 for harbor volumes and +105S for cross-sectional areas.
-1-

-------
Simulated velocities, however, differed by as much as 30056 from maximun
observed flood and ebb tide velocities. Part of this problem was
attributed to the fact that the simulated velocities were cross-sectional
averages, while the velocities published by Beverage and Swecker (1969)
were measured at single points, generally where the velocity was at a
maximum. Simulated velocities were approximately 70% of those obtained
from the physical hydraulic model study performed by the U.S. Army Corps
of Engineers.
Two other areas in which difficulty was encountered in hydrodynamic
model simulations were associated with predicting the timing of slack
tides at Montesano achieving the proper flow in the South Channel of
Grays Harbor. The timing of slack tide at Montesano was improved by
accounting for energy dissipated in the river bends between Aberdeen and
Montesano. The flow problem in the South Channel was corrected by
decreasing the channel depths and increasing the rate of aiergy
dissipation due to friction.
The analysis of the ecologic model was somewhat more difficult.
Observed salinity profiles were used to evaluate the success of the model
with regard to advection and diffusion-1 ike processes. Considerable
difficulty was encountered in the temperature simulations. This was
attributed to several factors, including the fact that the model does not
provide for heat exchange between the water and the exposed tidal flats,
insufficient knowledge of the dynamic boundary conditions at the ocean
boundary and numerical dispersion.
For the dissolved oxygen model, the difference between observed and
simulated values was +0.8 mg/1, or less. The average difference was 0.3
mg/1. In general, the maximun error occurred upstream from Aberdeen.
-2-

-------
The dynamic estuary model requires additional study before it can be
used for river basin planning in the Grays Harbor-Chehalis River system.
Nunerical dispersion appears to be the major source of error, although
lack of knowledge of exchange characteristics at the entrance to Grays
Harbor also contribute to the difference between simulated and observed
salinity, temperature and dissolved oxygen.
-3-

-------
-4-

-------
INTRODUCTION
Conditions of extreme low river flows, which occurred in the Pacific
Northwest during 1977, provided an opportunity to evaluate water quality
models under critical conditions. Grays Harbor, an estuary on the coast
of Washington (Figure 1), experienced water quality problems during such
critical periods. The National Pollution Discharge Elimination System
(NPDES), written for the Industries which discharge to the estuary,
contain limitations based upon freshwater discharge to the estuary.
Because the estuary is water quality limited, it was the subject of a
comprehensive water quality survey (Yearsley (1979)) and subsequent
analysis of the data, using mathematical models of water duality. Due to
time and resource limitations, the models were limited to those which
simulate temperature, salinity, dissolved oxygen and biochemical oxygen
demand.
Two different modeling concepts were applied to the analysis of the
data. The one described in this report was based upon a dynamic, quasi-
two-dimensional hydraulic and water quality model which accounts for
diurnal variations of water transport and surface elevations, as well as
processes which drive biochemical reactions, such as solar radiation.
_c_

-------
Figure 1. Grays Harbor and fresh water tributaries.

-------
The other concept (Yearsley and Hess (1979)) was based upon the
assumption that the estuary could be considered as a steady-state,
one-dimensional system governed by processes of horizontal diffusion,
non-tidal advection and first-order biochemical reactions, such as that
developed by Crim and Lovelace (1973).
The intent of these model studies was to examine the processes which
govern the distribution of dissolved oxygen in Grays Harbor and to
provide the basis for comparing a simple steady-state model and a more
sophisticated dynamic model. The final objective of the study was to
improve our methods for developing and evaluating WOES permits in Grays
Harbor.
-7-

-------
-3-

-------
PREVIOUS INVESTIGATIONS
Previous studies conducted on Grays Harbor provided greater insight
about the hydraulic, chemical, and biological processes occurring in the
estuary as well as to supply additional data used in the verification of
the models. The earliest, and probably the most comprehensive study was
done by Eriksen and Townsend (1940). Beverage and Swecker (1969) have
provided the most comprehensive report detailing the results of studies
undertaken by cooperative agreement between the U.S. Geological Survey
(USGS) and the Washington State Pollution Control Commission. The scope
of the Beverage and Swecker (1969) report focused on describing the
general circulation patterns and qual 1 ty characteristics of the water
mass as well as the influence of bottom materials on water quality
conditions.
The Washington State Department of Ecology (1970, 1971) supervised
the Cooperative Grays Harbor Surveillance Program. Specific data
resulting from these studies were used by Lorenzen et al (1974) to
develop a mathematical water quality model for Grays Harbor. The results
of this study provided the impetus for the use of the dynamic model in
this investigation. The work of Lorenzen et al (1974) does give a fairly
complete synopsis of other water quality investigations which have been
conducted on Grays Harbor in recent years. In addition, the U.S. Army
Corps of Engineers (Brogdon (1972)), has constructed a physical model of
Grays Harbor used primarily to evaluate the effects of dredging in the
estuary and to provide a better understanding of the hydraulic
characteristics of the harbor.
-9-

-------
-10-

-------
ESTUARY HYDRODYNAMIC MODEL
Theoretical Consideration
The primary objective of this project focused on the application of
a dynamic water quality estuary model to Grays Harbor, Washington (Figure
1). The model selected for this study was developed by Water Resources
Engineers, Inc. (WRE) (1975)of Walnut Creek, California under contracts
with the Office of Water Resources Research and the U.S. Environmental
Protection Agency (EPA). The Dynamic Estuary Model (DEM) as formulated
by WRE, has two modules which are operated sequentially.
The first module, ECOHYD, the estuary hydrodynanic model, simulates
the movement of water in an estuary that is under tidal influence. The
program solves the equations of longwave propagation in a shallow water
system. This 1s accomplished by a link-node model with a set of
one-dimensional motion equations for the channels or links and continuity
equations for the junctions or nodes. ECOHYD computes flows, voltmes,
depths, and cross-sectional areas as well as surface overflow rates,
reaeration coefficients, and dispersion coefficients for the entire
estuarial network at specific time intervals over the simulation period.
The concepts upon which ECOHYD are based were first described by
Shubinski et.al. (1965). Complete docunentation of the computer
software, as well as considerable technical analysis of the hydraulic
model was done by Feigner and Harris (1970).
In order to simulate the movement of water in an estuary that is
under tidal influence, ECOHYD solves the equations of longwave
propagation in a shallow water system. By segnenting the estuary into

-------
junctions (nodes) which represent storage volimes and into channels
(links) which connect the junctions and thereby allow the movement of
water throughout the estuary, the one-dimensional equations of continuity
and mcmentun are solved to predict the hydrodynamic characteristics of
Grays Harbor. The WRE model is analogous to a network of buckets and
hoses where the junctions are discrete volunes of water assumed to be
uniformly mixed at any given instant and the channels provide the
computational mechanism for the transfer of properties from junction to
junction. Figure 2 depicts the portion of Grays Harbor that was studied,
segmented into the system of 79 nodes and 100 links used for the model
calculations. Although the model gives the initial impression of being
two-dimensional in nature, the hydraulic behavior of the estuary is
actually described by the one-dimensional form of the equations of
continuity and momentum for open channel flow.
-12-

-------
Figure 2. Junction configuration for application of the {Dynamic Estuary Model in Grays Harbor. |

-------
-]A.

-------
Model Implementation
Hydraulic Geometry
This task focused upon an analysis of the physical and spatial
dimensions of the prototype wate- body in an effort to discretize the
estuary into the channels and junctions depicted in Figure 2. Specific
details describing the orientation of a link-node system are documented
in Feigner and Harris (1970). As a starting point, a link-node network
based on the work of Lorenzen et al (1974) was formulated. In order to
yield a better comparison between the simulated and observed values of
the study discussed herein, the initial network was modified so that EPA
survey stations were situated in the center of model junctions. The ocean
boundary for the network was taken to be the harbor entrance at Pt.
Chehalis. The upstream boundary for the dynamic model on the Chehalis
River was at Porter, Washington, where the tidal influence on flow and
surface elevation is negligible.
The channel parameters required for the operation of the model are
length, width, cross-sectional area, frictional resistance coefficient or
Manning's "n", and channel bottom elevation referenced to the model
datijn. Channel lengths tend to be controlled by the computational
stability requirements discussed in Feigner and Harris (1970). The model
elevation datun selected was the Mean Lower Low Water (MLWl line at
Point Chehalis near Westport, Washington.
Channel lengths, widths, and depths were taken directly from
Navigation Chart 18502 (Grays Harbor) published by the National Ocean
Survey. Because depths on the chart were referenced to MLLW at that
particular location, it was necessary to convert those elevations to
-15-

-------
reference MLLW at Point Chehalis. A table of bench mark elevations of
various local MLLW lines referenced to the National Geodetic Vertical
Datun (NGVD) has been published by the National Ocean Survey (1976).
Table 1 summarizes the bench mark elevations extracted from the National
Ocean Survey report which were used in the Grays Harbor investigation. A
straight linear interpolation between the points was performed in order
to obtain other intermediate MLLW elevations. In addition, because all
tidal stages published by the National Ocean Survey and the Corps of
Engineers also reference local MLLW lines, a similar conversion was
required to allow comparison between simulated and observed tidal heights.
Table 1


The difference between the National Geodetic Vertical Datim
lower low water at selected locations in Grays Harbor.
and mean
N05 Index
Map Nunber
Locality
NGVD-MLLW
(feet)
33
Point Chehalis and Vicinity,
Grays Harbor
4.60
34
8ay City, Grays Harbor
5.06
36
Aberdeen, Grays Harbor
5.38
37
Montesano, Chehalis River, Grays
Harbor
2.76
One adjustment to the depths of channels between Moon Island and
Montesano merits special notation. The channel depths input to the model
represent the average channel depth across a section while the depths
appearing on the National Ocean Survey chart are measured at
mid-channel. The average channel depth was based upon the assunption
that the channel shape could be approximated by a harmonic curve:
-16-

-------
d = do
(i)
where:
x = distance from midchannel (ft)
d » channel depth at x (ft)
dQ = maximum channel depth (at midchannel) (ft)
b = channel width (ft)
Because the channel widths in the link-node network for Grays Harbor
were generally more than an order of magnitude greater than the channel
depths, it was assuned that the hydraulic radius could be reasonably
approximated by the mean channel depth.
Once the channels of the network were estahlished, the definition of
the junctions of the nunerical model was a fairly straightforward task.
For a two-dimensional body, drawing the perpendicular bisectors of the
channels outlined the junction boundaries. The properties associated
with the junctions are water surface elevation, depth, water surface
area, volune, inflows, outflows, and side slope (change in junction
surface area with change in water surface elevation). Junction surface
-17-

-------
areas and depths were scaled from the same navigation chart used to
define the channel geometry. Initial junction vol tones were calculated
internally in the computer program by multiplying the surface area of the
junction by the mean nodal depth.
From the navigation chart, it was possible to approximate the water
surface area at MLLW and at Mean Higher High Water (MHHW) for each
junction of the network. In order to establish the role of the tidal
flats in the numerical model, it was assumed that the change in surface
area as the water level rises and falls could be described by a linear
relationship. To obtain the side slope, the difference between the
surface area at MHHW and at MLLW was divided by the diurnal tidal range
published by the National Ocean Survey. During the model simulation,
this value was then multiplied by the change in water surface elevation
to determine the increase or decrease in the junction surface area for
the period of interest.
Freshwater Discharge
The main tributaries which affect the freshwater hydraulic and water
quality properties of the harbor are the Chehalis, Satsop, Wynoochee,
Wishkah, Hoquiam, and Himptulips Rivers. Streanrflow data for the
investigation was obtained from a provisional report provided by the USGS
(1977). Streamgage stations were located on the Chehalis River at
Porter, on the Satsop River, on the Wynoochee River, and on the
Humptulips River. Estimates were made by the USGS for freshwater flow on
the Chehalis River at its mouth at Hoquiam based on measurements taken at
Grand Mound. 8ecause the Hoquiam and Wishkah Rivers are not gaged, it
-18-

-------
was necessary to approximate these flows. It was assuned that the
difference between the estimated flow of the Chehalis Rive* at Hoaquiam
and the flow at the confluence of the Chehalis and Wynoochee Rivers
represented the freshwater contribution of the Hoquiam and Wishkah
Rivers. Based on tributary drainage areas provided by Richardson (1962),
it was assuned that the flow in the Hoquiam River was equal to the flow
of the Wishkah River due to the fact that the drainage areas were
approximately equal. A similar procedure of utilizing drainage areas was
applied to obtain flow estimates for the Elk and Johns Rivers based on
the recorded streamflow of the HunptuHps River. Table 2 summarizes the
flow data collected for the five-day survey period which was employed in
model simulation runs of the investigation discussed herein.
Initial Conditions
Before solving the equations of continuity and mcmentun, it was
necessary to supply the hydraulic program with a set of initial
conditions for flow rates in each channel and for water surface
elevations at each junction.
For this analysis, nodal water surface elevations were set at 9.5
feet and channel velocities were initialized to 0.0 feet per second
(fps). With these values, quasi-steady-state conditions were attained
within two tidal cycles.
-19-

-------
Boundary Conditions
Tidal boundary conditions supplied to the model exert a significant
influence on the nunerlcal simulation of the hydraulic characteristics of
the harbor. Because the hydrodynamic behavior of an estuary is dominated
by tidal action, any change that occurs in the tidal state at the ocean
boundary, either spatial or temporal, affects the entire estuarial
network. For this reason, the tidal conditions imposed at the seaward
boundary should be as closely representative as possible to those found
during the simulation period.
In ECOHYD, the four extreme stages of the tide (higher high, lower
low, lower high, higher low) are used in a harmonic least-squares
regression analysis to obtain the best fit described by a seven term
sinusoidal equation. Data from the National Ocean Survey Tide Tables for
July 25-29, 1977, at the entrance to Grays Harbor (Pt. Chehalis) were
used in the regression analysis.
The tidal input to the model was formulated from tide gauge readings
recorded on strip-charts at Port Dock obtained frcm the Corps of
Engineers. To estimate the four extreme stages of the tide at the ocean
boundary (Point Chehalis), it was necessary to apply a time/stage
difference equation to the Port Dock values. The relationship employed
was reconmended in the National Ocean Survey Tide Tables (1976) and is
shown in Table 2.
-20-

-------
Table 2.
Tidal time/stage relationship between tide gauge reading at
Port Dock and other locations in Grays Harbor, Washington
Location
Relative Time of	Water Surface Elevation
Stage of Occurrence (in minutes) Relative to M-LW (in feet)
Tide (PDT=Time at Port Dock) (P0H=Height at Port Dock)
Point Chehalis High water
Point Chehalis Low water
Morrtesano High water
Montesano Low water
POT - 32
POT - 43
POT + 81
PUT + 108
PO -1.1
PDH - 0.1
0.80 * PDH
0.53 * PDH
The final boundary conditions required for the operation of ECOHYD
were the daily inflow-outflow conditions at each junction. These
conditions include net flows from import, export, and waste effluents as
well as evaporation and wind stress. In the Grays Harbor investigation,
water losses due to evaporation were assuned to be of negligible
consequence. This assumption was based on National Weather Service
estimates of annual evaporation for the Washington coast. Using the
surface areas of the harbor, the potential change in volune due to
evaporation was calculated to be substantially less than volune
deviations caused by the uncertainty associated with average junction
depths taken from the nautical charts. Daily river flows and waste
discharges were input to the hydraulic model at the junctions reported in
Tables 3 and 4.
-21-

-------
Simulated Water Surface Evaluations
For initial runs, the model output was compared to the results of
studies conducted by 8everage and Swecker (1969). Data published in
their report included estimates of the longitudinal variation of the
cumulative volune of the harbor at Mean High Water (MHW), at Mean Low
Water (MLW), and at Mean Tide Level (MTL), estimates of the longitudinal
variation of the harbor cross-sectional area at MLW and at MTL,
measurements of mean velocities at selected sites for different
freshwater flows, and predicted tidal stages for the locations where
velocity measurements were performed. U.S. Coast and Geodetic Survey
tide tables (1966) as well as US6S flow records (1966) were referenced to
supplement the data published by Beverage and Swecker. Data reported by
Eriksen and Townsend (1940) also provided further information about the
physical characteristics of Grays Harbor which was used for model
calibration. Simulated tide heights are compared to predicted heights at
Montesano and observed heights at Aberdeen (Port Dock) in Figures 3 and
4, respectively.
Simulated Channel Velocities
For the initial hydraulic simulation runs which employed the inputs
outlined in the preceding section, the difference between published data
and simulated values was jO.5 feet for tidal stages +15% for harbor
volumes, and +10% for cross-sectional areas was noted. Considering the
roughness of the observed values as well as the simplifying assunptions
required for the nimerical computations, the hydraulic behavior of the
model for those parameters seemed adequate. However, the correspondence
between measured and simulated velocities left much to be desired.
-22-

-------
Figure 3
?~rTidais--
For Olfforent va|U08 Manning's "n"
* »	4"	L>	.
36	48			
00	72
time in HOURS
Manning's "n" Iron, Table 4.3
Manning's of 0 040
T«te Table Piethctions from Table 4-2

-------
Figure 4
At Port Dock1 S,niU,ated and Observed r*
(July 25-29, 1977)	Sta00s
12.6
Calibrated Model
0 Observed Values

-------
Maximim velocities reported in Beverage and Swecker (1969) ranged from
3.0 fps on flood tides to 4.5 fps on ebbtides while the maximum simulated
velocities varied from 1 fps on floodtides to 1.5 fps on ebbtides. One
potential cause of differences stenned from the fact that the velocity
measurements presented 1n Beverage and Swecker (1969) were made at the
edge of the navigation channel while velocities calculated in the model
were cross-sectional averages.
The physical model of Grays Harbor, built by the U.S. Army Corps of
Engineers Waterways Experiment Station, Vicksburg, Mississippi, provided
another source for evaluating the ECOHYD simulations. In 1972, the Corps
of Engineers completed construction and initial base testing of the
model. The model contained stations at which measurements of water
surface elevation, channel velocity, salinity, and dye concentration
(representing potential pollutant discharges) were made. Most
cross-sectional average velocities for the channels of interest ranged
fran 1.5 to 2.0 fps both on floodtides as well as on ebbtides for a
freshwater flow of 1,270 cfs (cubic feet per second) on the physical
model. (The freshwate- discharge during the period of investigation was
about 1,200 cfs).
The observations pertaining to the channel widths did lead to one
modification. The version of the WRE model that was used for Grays
Harbor provided for a junction side slope. As was observed, exposure of
the tidal mud flats reduced the network channel widths. To account for
this, channel widths were increased or decreased proportionally to the
average chaige in surface area of the upstream and dwnstream junctions.
This modification greatly improved the simulated velocities relative to
the physical model.
-25-

-------
Advection and Circulation Problems
To gain insight about the circulation pattern in Grays Harbor and,
in particular, to provide a comparison between flow rates in the North
and South Channels, photographs made by the Corps of Engineers of their
physical model were studied. Analysis of the pictures indicated that the
velocities in both channels were roughly equal. However, initial
simulations with the numerical model implied that the wave was propagated
up the South Channel at a much faster rate than through the North
Channel. Channel depths and Manning's coefficients in the South Channel
were adjusted such that the wave was propagated through both channels at
roughly the sane rates.
Model Sensitivity to Manning's Coefficient
According to the National Ocean Survey (Table 3), the tidal wave
normally requires one hour and twenty minutes to reach Montesano from
Aberdeen under the average annual flow conditions 4,000 cfs). However,
the hydraulic model indicated that high slack at Montesano appeared less
than thirty minutes after its occurrence at Aberdeen. Because flow
conditions during the survey period were lower than the annual mean, it
would be expected that the wave would be propagated somewhat faster than
the rate predicted by the National Ocean Survey. Unfortunately, no tide
gauge was installed at Montesano that could have established the time of
occurrence and the water surface elevation associated with slack tide.
Furthermore, the model simulated a diurnal tidal range of 11 to 12 feet
at Montesano while the National Ocean Survey predictions indicated that
the actual range was around eight feet.
-26-

-------
In an attempt to isolate the problem, the hydraulic model was run
with a freshwater flow of 4,000 cfs, i.e., the mean annual flow rate upon
which the National Ocean Survey predictions were based. The time of
travel of the simulated tidal wave increased by only ten minutes, while
the diurnal range was reduced by only 0.5 foot. Consequently, this
analysis indicated that a difficulty existed with the numerical
representation of the Chehalis River between Aberdeen and Montesano.
Nunerous approaches were taken to rectify the problem. In the original
model, the upstream boundary was located at the State Highway 107 bridge
near Montesano. Because the water surface elevation at that point is not
fixed, it was felt that extending the boundary beyond the tidal fall line
would reduce the problem. Hence, the boundary used in ECOHYD was
relocated, 14 miles upstream frcm Montesano, where fluctuations in the
water surface elevation are caused by changes in freshwater flow rather
than by tidal forces. This modification resulted in only a small
improvanent of the model behavior.
While Lorenzen et al (1974) made no reference to model predictions
of the tidal stages at Montesano, the Corps of Engineers reported
difficulties with the physical model at that location which were similar
to those encountered in the study discussed herein. Because their model
did not appear to account for energy losses encountered in the bends,
attention was next focused on improving the model description of this
phenomena.
No known studies have been conducted on Grays Harbor which provide
measuranents of Manning's roughness coefficients. Consequently, the
values were initially estimated from the work of Henderson (1966) for
natural river channels. Downstream from Moon Island, the values selected

-------
for Manning's "n" appeared to be of gnall consequence due to the strong
dominance of tidal effects at the ocean boundary. However, model
calculations were fairly sensitive to resistance coefficients between
Montesano and Aberdeen. The modeling problem is complicated by the
presence of ninerous bends and sloughs in this reach. It is well known
that energy losses are created by changes in the direction of flew as
well as by lateral inflows and outflows. However, information which
could aid in quantifying these losses is scarce, especially in natural
channels where precious little experimentation has been conducted. When
attempting to calculate the losses, a multitude of variables must be
taken into consideration including channel width, depth of flow, angle of
deflection, strength of secondary currents, and suspended sediments to
nane a few.
Henderson (1966) has recommended that one method to include the loss
coefficient, C^, due to bends when utilizing segnented channels for
hydraulic computations is to raise the value of Manning's "n"
sufficiently to account for CL< One set of experiments has shown that
losses due to bends may be estimated by:
Cl = ik
(31
-23-

-------
where:
= energy-loss coefficient for channel bends
(dimension!ess).
b 3 average channel width (ft)
r„ = radius of channel curvature (ft)
c
Otter sets of data have indicated that decreases as the ratio
of channel depth to width decreases. As an extremely crude attempt to
include head losses due to bends in ECOHYD, was estimated by the
following. The empirical relationship based on the findings outlined
above:
Cl - zA.
r.
(4)
Following an estimate of from Equation (4), Manning's "n" was
increased by the amount that would generate an equivalent increase in the
friction term of the mamentun equation. In addition, channel lengths in
those regions were segnented such that the head loss was applied solely
through the duration of the bend. Table 5 sunmarizes the calculations
utilized to increase Manning's "n" for the bends in the Chehalis River as
-29-

-------
depicted by the representation of Grays Harbor in ECOHYD. While this
technique 1s highly empirical, its implementation moderately improved the
correspondence between simulated and predicted (from tide tables) water
surface elevations at Montesano.
Figure 3 depicts the improvement that was attained by increasing the
Manning's coefficients to the values finally used as input to the Estuary
Hydrod>namic Model and listed in Table 5. The dashed line represents the
hydraulic model simulation in which Manning's coefficient in all channels
between Aberdeen and Montesano was equal to 0.040.
-30-

-------
Table 3
Freshwater Discharges in Grays Harbor, Washington (July 25-29, 1977)
Tributary
Recei ving
Node
7-25 Flow
7-26 Flow
(cf si
7-27 Flow
(cf s)
7-28 Flow
(cf si
7-29 Flow
(cf s)
Chehal is R 1ver at
Porter +
68
416.0
401.0
3^5.0
395.0
398.0
Satsop River +
67
378.0
374.0
374 .0
374 .0
369.0
Wynoochee River +
64
211.0
208.0
214.0
208.0
204 .0
Wlshkah RIver*
79
77.5
78.5
78.5
81.5
74.5
Hoqulam River*
72
38.8
39.3
39.3
40.8
37.3
E. Fork Hoqulam
R iver*
74
38.7
39.2
39.2
40.7
37.2
Chehal 1s R 1ver at
Hoquiam**
—
1160.0
1140.0
1140.0
1140.0
1120.0
Huiiptullps River +
29
164.0
160.0
156.0
151.0
148.0
Elk River*
20
12.2
11.9
11.6
11.2
11.0
Johns River*
19
21.0
20.4
19.9
19.3
18.9
+Reported in USGS provisional data
~Estimated by tributary drainage area approximation
**Est1mated by USGS

-------
Table 4
Major Waste Dischargers 1n Grays Harbor, Washington (July 25-29, 1977)*
Industri al/
Municipal
Source
Date
Recel vi ng
Node
Flow
(cfs)
Temp
(c)
BOD
(mg/1)
DO
(mg/1)
P04-P
(mg/1)
pH
NH3-N
(mg/1)
N03-N
(mg/1)
ITT
Rayonler (1)
7-25
11
44.3**
32.9
37.6**
4.0
0.030
7.0**
0.030
0.025
11
7-26
11
45.5
32.9
60.8
2.6
0.910
6.5
0.030
0.025
II
7-27
11
46.5
32.9
59.9
3.0
1.600
6.2
0. 580
0.018
II
7-28
11
46.4
32.9
74.7
3.0
1.500
7.1
0.172
0.004
II
7-29
11
45.6
32.9
2.0
3.0
1.500
7.1
0.172
0.004

Weyerhaeus er
#1 (2)
7-25
24
31.9**
23.3
101.7** 4.0
1.100
6.7**
0.030
0.025
II
7-26
24
34.3
23.3
125.6
4.0
1.400
7.0
0.030
0.025
II
7-27
24
31.3
25.7
113.7
4.0
1.100
6.8
0.580
0.018
II
7-28
24
35.8
©
•
C\i
143.5
6.2
1.600
6.3
0.172
0.004
II
7-29
24
33.3
28.9
128.6
6.2
1.500
6.6
0.172
0.004

-------
Table 4 (cont'd)
Major Waste Dischargers in Grays Harbor, Washington (July 25-29, 1977)*
Industrial/
Municipal
Source
Date
Receiving
Node
Flow
(cfs)
Temp
(c)
BOD
fmg/1)
DO
(mg/1)
P04-P
fmg/1)
pH
NH3-N
fmg/1)
N03-N
fmg/1)
Weyerhaeuser
#2 (2)
7-25
46
1.24**
37.0
4.5
7.6
0.050
8.4
0.030
0.025
<1
7-26
46
1.24
37.0
20.9
7.6
0.050
8.4
0.030
0.025
II
7-27
46
1.24
27.2
26.9
5.2
0.045
7.8
0.580
0.018
II
7-28
46
1.24
46.1
9.0
7.2
0.092
8.5
0.172
0.004
II
7-29
46
1.24
21.9
10.5
7.2
0.092
8.2
0.172
0.004

Aberdeen
STP (3)
7-25
37
3.5**
17.3
82.0
1.3
2.800
6.7
19.000
0.006
II
7-26
37
3.5
17.3
82.0
1.1
2.800
6.7
19.000
0.006
u
7-27
37
3.7
17.3
131.0
1.0
2.800
6.7
19.000
0.006
ii
7-28
37
3.7
16.9
97.0
0.9
3.800
6.8
14.000
0.004
7-29	37	3.7	16.8 97.0	1.0 2.500 6.8 13.000 0.010

-------
Table 4 (cont'd)
Major Waste Dischargers in Grays Harbor, Washington (July 25-29, 1977)*
Industri al/
Municipal
So ur ce
Date
Recel vi ng
Node
Flow
(cfs)
Temp
(cf
BOD
(mg/1)
DO
(rag/1)
P04-P
(rag/1)
PH
NH3 -N
(mg/1)
N03-N
(mg/1)
Hoqui am
SIP (4)
7-25
9
2.6**
18.5
53.0
16.6**
3.600
8.8**
1.730
0.041
II
7-26
9
1.9
18.5
53.0
9.0
3.600
8.8
1.730
0.041
It
7-27
9
1.9
18.5
53.0
12.7
3.600
8.8
1.730
0.041
II
7-28
9
4.6
18.5
53.0
9.5
3.600
8.8
1.730
0.041
IJ
7-29
9
4.6
18.5
53.0
18.0
3.600
8.8
1.730
0.041
~All values obtained from EPA survey unless otherwise noted.
**Average value for July 16-25.
(1)	Flow, BOD, and pH reported by ITT Rayonier to WSDOE (Washington State Department of
Ecology)
(2)	Flew, BOD, and pH reported by Weyerhaeuser to WSDOE.
(3)	Flow, DO, and pH reported by Aberdeen to WSDOE.
(4)	Flow, DO, and pH reported by Hoquiam to WSDOE.

-------
Table 5.
Calculations of Manning's "
n" to account for
energy i
osses due
to river
bends







Radius of
Channel
Manning's "n"
Channel
Downstrean
Upstream
Curvature, rc
depth, d
Used in Model
Number
Node*
Node *
(ft)
(ft)
0.040**+2drc
62
40
41
4,800
24.2
0.050
63
41
42
2,000
24.5
0.065
66
44
45
10,200
25.5
0.045
67
45
46
10,400
25.9
0.045
69
47
48
2,100
26.1
0.065
70
48
49
2,100
26.0
0.065
71
49
50
5,100
25.5
0.050
72
50
51
1,200
24.6
0.080
73
51
52
3,200
23.7
0.055
74
52
53
3,000
22.6
0.055
75
53
54
4,300
21.4
0.050
76
54
55
7,500
20.1
0.045
77
55
56
7,500
19.5
0.045
78
56
57
7,500
19.4
0.045
79
57
58
7,500
19.3
0,045
80
58
59
7,500
19.0
0.045
81
59
60
7,500
18.5
0.045
82
60
61
7,500
18.1
0.045
83
61
62
7,500
17.6
0.045
84
62
63
7,500
17.2
0.045
85
63
64
7,500
16.5
0.045
* Reference Figure 2
** Manning's "n" of H«iderson (1966) for natural river channels
Model Sensitivity to the Freshwater Flow in the Hoquiam and Wishkah Rivers
As stated earlier, the freshwater flow in the Chehalis River at its
mouth at Hoquiam is only approximated by the LEGS . The major sources
of freshwater flow below the farthest downstream gaging points of the
Chehalis River Basin are the Hoquiam and Wishkah Rivers.
-35-

-------
-36-

-------
ESTUARY ECOLOGIC MODEL
The second module of the Dynamic Estuary Model, ECOSIM, the estuary
ecological model, simulates water quality conditions in an estuarine
environment based upon one-dimensional mass balance equations which
include advection, diffusion, and chemical and biological reactions.
ECOSIM utilizes the output from ECOHYD to solve the mass balance
equations for these parameters at each junction. The model is dynamic
and can simulate processes whose periods are equal to, or greater than,
one day. The basic concepts were first described by Chen and Orlob
(1972). A more recent description, including the governing equations for
all reactions, is given by WRE (1975).
A pollutant discharged into a body of water 1s transported primarily
by one of two mechanisms, advection and diffusion. In addition, the
concentration of a constituent may also be influenced by importation and
exportation, by biological and/or chemical decay, and by mass transfer at
the air-water interface. For an estuarial situation such as Grays
Harbor, advection processes are normally an order of magnitude larger
than the other mechanisms, mainly due to tidal fluctuations. The current
trend in the mathematical modeling of water quality in estuaries is to
attempt to Improve the description of the hydraulic parameters
(velocities, flows, volumes) which play a key role in advective
transport. For this reason, ECOHYD, the hydrodynanic model, is run prior
to water quality simulation. Flows and volunes in the estuary, as
simulated by ECOHYD, are averaged over a time period corresponding to the
time period used for the ecologic simulations. The flows and volunes
-37-

-------
computed in this way are then included in the mass balance equations for
the various water quality sub-models in ECOSIM, the estuary ecologic
model. Because the literature is rich with descriptions of the different
water quality submodels, discussion of the theoretical considerations
associated with these mechanisms will be limited. Baca, et al (1973)
provided an excellent suimary of the mathematical rate equations used to
describe the reactions which influence the various water quality
parameters.
ECOSIM as used in this investigation contained finite-difference
equations for 14 different water quality constituents. The parameters,
and the interactions between parameters are shown 1n Table 6. However,
data availablillty and time limitations were such that it was feasible to
simulate temperature, salinity, dissolved oxygen and biochemical oxygen
demand, only.
-38-

-------
Table 6
Water Quality Parameter Interactions Utilized in WRE
Estuary Ecologic Model
OUTPUT
Water Temp.
Detritus
BOD
(Carbonaceous)
Oxygen
Phosphate
Ammonia
Nitrite
Nitrate
| Algae (2 types)
Zooplankton
BentUic
Animals
Sediment
Secchi Disk
Water
Temperature
X
X
X
X
X
X
X
X
X
X
X
X
X
Detritus

X

X
X
X





X
X
SOD
(carbonaceous)


X
X









Oxygen



X






X


Phosphate




X



X




Aumonia



X

X
X

X




Nitrite



X


X
X





Nitrate







X
X




Algae
(two types)



X
X
X

X
X
X

X
X
Zooplankton

X

X
X
X


X
X


X
Benthie
Animals



X
X
X




X
X

Sediment



X
X
X




X
X

Secchi Disk








X




-39-

-------
Meteorological Inputs
In addition to the hydraulic data generated by the ECOHYD, another
set of data required for the ecologic model was meteorological conditions
for the simulation period to be used in the temperature and phytoplankton
submodels. Values for atmospheric pressure, cloud cover, wind speed, dry
bulb temperature, and wet bulb temperature were needed to compute the
energy budget in the form:
where:
«t = °sn * "at - % - % + «c	f5)
on
= net energy flux across the air-water interface
Q = net short wave radiation
sn
Qat = net atmospheric (long-wave) radiati
Q = water surface radiation
fr
Qg = energy loss by evaporation
Qc = convective energy flux
10

-------
Each term of the energy budget is calculated by utilizing the
algebraic equation derived from heat transfer theory and/or empirical
relationships. Specific details of the equations which comprise the
temperature submodel may be found in Tennessee Valley Authority (1968)
and WRE (1968). The equations used for the ECOSIM are suranarized in WRE
(1975). Input weather data was obtained from the monthly surmary of local
cl imatological data published by the Environmental Data Service of NQAA
assessing the values at Astoria, Oregon; Olynpia, Washington, and
Quillayute, Washington as well as field notes collected during the
survey. In general, it was found that the average of weather
measurements from the Astoria and Quillayute stations were reasonable
Indicators for the weather in Grays Harbor during the period, July 25-29,
1977.
Initial Conditions
The next input required by the ecologic model was a set of Initial
conditions. Unlike the hydraulic model, initial water quality conditions
exerted a significant impact on the simulated results. In order to
minimize the effect of initial conditions, the model was run for a
ten-day period prior to actual simulation. Average waste discharge and
water quality values reported by the industries before the actual survey
(Table 4) were utilized as inputs.
Initial conditions for salinity, temperature, and dissolved oxygen
resulted fran a harmonic regression analysis using Hydrolab data obtained
during the survey period. As stated in the discussion of the hydraulic
model, water surface elevations in a system under tidal influence
generally can be described by a sinusoidal relationship. Because most
-41-

-------
water quality parameters are greatly influenced by advective transport
and because tidal action usually controls the movement of the water in an
estuary, it was felt that the fluctuations of the tidally dominated
constituents could also be described by a sinusoidal relationship. To
determine an equation describing tidal fluctuations of measured water
quality parameters, a harmonic regression analysis was used. In
addition, this analysis was used to smooth data points in order to
generate longitudinal profiles of the various water quality parameters at
"slack times" as well as to establish initial conditions for the dynamic
model. Details of this analysis can be found in Cleland f 1978).
Initial conditions for ultimate biological oxygen demand fBOD) were
obtained from laboratory analyses of water samples collected during the
EPA survey Yearsley (1979).
Point Source Loadings
Wastewater discharges from the paper products industry exert the
greatest influence on water quality in Grays Harbor. Waste loads frcm
municipal sewage treatment plants are small in comparison to discharges
from the Weyerhaeuser and the ITT-Rayonier pulp and paper mills. Table 4
outlines the major waste sources in Grays Harbor and the average values
of the quality of discharges for the five-day survey period. Grab
samples were normally taken at each location twice daily by EPA. Values
reported to the Washington State Department of Ecology (1977) by
ITT-Rayonier, by Weyerhaeuser, and by the Hoquiam and Aberdeen sewage
treatment plants were also utilized. In addition, EPA placed a
continuous sampler at each major discharge location once during the week
for a twenty-four hour period to obtain a composite average.
- '1 ^ .

-------
Deoxygenation Rate
From the regression analysis on the BOD lab samples the,
deoxygenation rate, K.j, was found to be 0.175 day (base e).
Reaeration Coefficient
A multitude of equations have been proposed to estimate the
reaeration coefficient, Kg. Beverage and Swecker (1969) mentioned the
Langbein-Ourun (1967) equation in their report, but gave no supporting
basis for selecting that particular formula. Many factors must be
considered when choosing a reaeration equation, but most researchers
agree that velocity and depth of flow play a dominant role in the
process. Grays Harbor presents a unique situation because upstream from
Aberdeen, river flow governs the hydraulic characteristics.
Consequently, many of the reaeration formulas proposed for natural
streans could be considered applicable. For this investigation, the
major reaeration equations utilized in most water quality models were
identified from Bennett and Rathbun (1972). Values for various flow
conditions found in Grays Harbor at different locations were then
substituted Into these equations and resulting reaeration coefficients
were compared. The one ultimately used in the model was the
Langbein-Ourun (1967) equation, though both the Langbein-Ourun (1967) and
the O'Connor-Dobbins (1958) formulations produced nearly Identical values
for K2. The formula appears as:
-43-

-------
K2 - T.GyZ lull
1 .5^
where:
= reaeration rate (base e)
= velocity in channel i (ft/sec)
= hydraulic radius or depth in channel i (ft)
Reaeration coefficients were calculated for each junction at
ecologic time intervals in the hydraulic model in order to account for
*
the velocity fluctuations that occur over a tidal cycle. was
computed for each channel in the network at ecologic time steps, steps
based on Equation (6). The original WRE model then derived the junction
reaeration rate by taking the weighted average of all channel K^'s
connected to the junction based on channel surface areas.
This computation was modified for the Grays Harbor study to base the
weighted average on channel flow rates, rather than channel surface
areas, in order to reflect the volune of water subjected to the different
reaeration rates.
-dd-

-------
RESULTS
Simulated concentrations of temperature, dissolved oxygen, salinity,
and BOD were compared with data obtained from a five-day intensive survey
conducted by the U.S. Environmental Protection Agency during the period
of July 25-29, 1977 (Yearsley (1979)). The survey was designed to obtain
data for evaluating both the dynamic estuary model and a steady-state
water quality model. Water quality measurements were made both in the
receiving waters, as well as in the major waste discharges and freshwater
sources tributary to Grays Harbor.
Salinity
It is generally easiest to begin the analysis of model simulations
with the prototype by using a conservative constituent. Conservative
constituents are, by definition, subject to change by mixing processes,
only. Because of this, conservative constituents can be used to evaluate
how well the model simulates processes such as diffusion and advection.
In the case of Grays Harbor, salinity was the most convenient and
accessible conservative constituent to use for this purpose.
For the salinity simulations we used a time step of 45 seconds in
the hydraulic computations. For the time step in the ecologic model we
initially used 60 minutes. However, this resulted in substantial
inaccuracies in the salinity simulation. Accuracy was improved by
reducing the time step in the ecologic model to 15 minutes.
-45-

-------
One other problem confronted in the salinity simulations was closely
related to the difficulties experienced with the hydraulic model. The
simulated tidal wave appeared to be propogated up the South Channel at a
greater rate than was expected based on salinity observations and on
velocity predictions from the Corps of Engineers physical model. After
analyzing simulated flow rates, it appeared that salt was being adverted
too rapidly through the South Channel and that a decrease in channel
velocities in that region would remedy the problem. As discussed in the
Estuary Hydrodynamic Model section, South Channel bottom elevations were
decreased and Manning coefficients were increased. As a result, greater
correspondence was also noted between simulated and observed salinities
in the South Channel.
Temperature
After obtaining reasonable results from the simulation of salinity,
attention was focused on temperature, since it affects growth and decay
of a nunber of water quality parameters. It is generally assuned that
estuarine temperatures are governed by advective and diffusive processes,
by heat exchange at the air-water interface (Equation (5)), and by the
temperature of the fresh water inflow, as well as conditions at the ocean
boundary.
The tidal mud flats play an important role in the movement and
mixing of esturial waters. In addition, these regions could influence
the heating of waters in Grays Harbor. Relatively little work has been
done on the effects of mud flats on the estuarial distribution of
temperature. It is well known that estuarial mud flats can absorb
-£e-

-------
significant amounts of heat which could be transmitted back to the water
during floodtide. For this reason, a crude submodel was constructed to
evaluate the potential impact of exposed mud flats on water temperature.
The mud flats energy balance was incorporated Into ECOSIM by
assuning that the net solar radiation on an exposed area would be stored
in the surface layer of the mud until the next high tide. The total
©lergy stored during the period the mud flats were exposed was then
transferred to the water col urn once the mud flats were covered with
water.
While the model did include the incident solar radiation, the time
of exposure of the mud flats, the amount of surface area exposed, and the
poteitial maximun amount of energy that could be absorbed, it did neglect
factors such as heat losses and transfer rates. However, the primary
goal of this analysis was to evaluate the poteitial effect of the heating
of tidal mud flats on the water temperature of Grays Harbor in order to
determine whether or not a more detailed thermodynamic analysis would be
a worthwhile endeavor. Using an absorption ratio for mud to water of
1.5, we found that the maximun increase in simulated water temperature at
any of the stations was approximately 0.6°C.
The last factor examined in the temperature submodel analysis was
the impact of the ocean boundary on the temperature of harbor waters.
ECOSIM software was not designed to allow time-varying boundary
conditions. On flood tides the incoming water contains a mixture of
ocean and estuary water. We have assuned that the mixture of ocean and
-47-

-------
estuary water changes as the flood tide advances. We used the following
empirical relationship to describe the mixture of ocean and estuarine
water:
Ti.
r
T0
i.o +
U	(7)
where:
T. = temperature of water advected into downstream node of
estuary
Tq = water temperature of ocean boundary ffixed)
= water temperature of downstream node of estuary
t = time from start of flood tide (hr)
6.25 = approximate duration of flood tide fhr)
r = recycle ratio
The recycle ratio, r, represents an empirical coefficient. For the
Grays Harbor study, we used a recycle ratio of 2. The use of Equation
(12) is justified only by the fact that the resulting simulations
correspond well with observed temperature values. However, it does serve
to demonstrate the need for establishing the time-varying water quality
characteristics at the ocean boundary. In the analysis of the thermal
model, the temperature at the ocean boundary was set to 13.5° based on
average temperatures for that particular time of the year reported by
Steffansson and Richards (1964).
-AJ3-

-------
Dissolved Oxygen
Dissolved oxygen 1s one of the most important water quality
parameters. Not only is it a key element 1n the activities of most
aquatic biota, but it is also a useful measure of the capacity of a water
body to assimilate waste materials. The oxygen submodel can beccme quite
complex because of the biological, chemical, and physical variables which
influence the concentration. Processes such as reaeration and
photosynthesis can add to the pool while the oxidation of pollutants,
vegetable decay, and plant respiration tend to deplete dissolved oxygen.
The factors which contribute to changes in dissolved oxygen levels other
than advection and diffusion which are included in the WRE Estuary
Ecologic Model are: 1) temperature variations, 2) BOD decay, 3) ammonia
oxidation, 4) nitrite oxidation, 5) reaeration, 6) photosynthesis, 7)
biota respiration, and 8) detritus oxidation. Initially the dissolved
oxygen budget for Grays Harbor included only advection, diffusion, BOD
decay, and reaeration. Considering only these processes, the dissolved
oxygen sag occurred between Elliot Slough and Port Dock and attained a
value of 6.6 mg/1. Survey measurements indicated minimun D.O. levels in
this region were around 5.5 mg/1. Other investigations (Yearsley and
Houck (1973), Lorenzen et al (1974)) of Grays Harbor have indicated that
sediment oxygen demand (SOD) can represent a significant sink for
dissolved oxygen resources. Generally, benthic demand originates from
the settling of suspended wastes which produces benthal deposits and
sludge banks. In developing the estuary ecologic model, WRE (1975) has
made an attempt to describe the processes contributing to SOD. However,
the submodel requires a significant amount of data (composition of
organic wastes, typical settling velocities, contribution frcm aquatic
biota, etc.).
-A9-

-------
Rather than assemble the vast quantities of data required for the
benthic submodel, an attempt was made to establish a relationship between
SOD rates and the results of other sediment studies, as done by Yearsley
and Houck (1973). Beverage and Swecker (1969) have reported data which
details the longitudinal distribution of the carbon content of bottom
materials. Their studies indicate that there was a direct relationship
between the SOD and the carbon content of the bottom materials.
Kreizenbeck (1973) determined SOD at three locations in Grays Harbor. It
was noted that both curves were similar in shape and both contained a
modest amount of uncertainty associated with the data. This similarity
was used to estimate the benthic oxygen demand in Grays Harbor.
Implementation of this approach did lead to a much more reasonable
agreement between simulated and observed values for dissolved oxygen.
One major pitfall of this technique for estimating SOD is the fact
that the scheme fails to account for the dredging activity that has taken
place in recent years. Dredging could potentially remove the oxygen
demand from the problem area by transferring the sediments to another
location. In this case, the dredge-spoil areas would suffer a
degradation in water quality resulting from a depletion in dissolved
oxygen. Essentially, the problem is not eliminated, but simply
relocated. However, it is felt that dredging activity has not
significantly changed the calculated SOD rates. Reasons for this stem
from the fact that most of the oxygen conscription by bottom materials
occurs at the sediment water interface regardless of the thickness of the
or gan i c 1 ayer.
-en.

-------
When dredging occurs, the benthic layer is disturbed and becomes
mixed with harbor waters. As a result, bottom material not previously in
contact with the water will also begin to consume oxygen. Due to the
increased denand, dissolved oxygen concentrations decrease far below
normal as has been observed during dredging operations. Upon completion
of dredging, a certain amount of SOD is exported with the dredge spoils.
However, the remainder which is suspended in the harbor waters settles
out covering approximately the same area as prior to dredging depending
on tidal velocity, tidal flow rate, and direction at the time. However,
as with the recycle ratio in the temperature submodel, this reasoning
remains solely hypothetical without further supporting evidence.
Discussion
Comparison of simulated and observed salinity, tanperature and
dissolved oxygen are shown in Figures 5-10. The maximun error associated
with the salinity simulations was 6.0 p.p.t., and occurred in that region
of the estuary with the steepest longitudinal salinity gradient (Figures
5 and 6). Simulated temperatures in this same part of the estuary
(Figure 8) differed from depth-averaged observed temperatures by 0.5 C to
1.5 C. Simulated dissolved oxygen concentrations (Figures 9 and 101
differed from observed values by less than 0.3 mg/1 near the ocean
boundary, but were as high as 0.8 mg/1 in the inner estuary. The average
difference between observed and simulated dissolved oxygen concentrations
for the entire estuary at both higfi and low tides was 0.3 mg/1.
-51-

-------
Figure 5
Comparison of Simulated and Observed Salinities
In Gray$ Harbor at Second Low
(July 25-29. 1977)
40
DISTANCE DOWNSTREAM FROM MONTESANO IN NAUTICAL MILES

-------
Figure 6
Comparison of Simulated and Observed Salinities
In Grays Harbor at First High
(July 25-29, 1977)
40
DISTANCE DOWNSTREAM FROM MONTESANO IN NAUTICAL MILES

-------
Figure 7
Comparison of Simulated and Observed Temperatures
Near Light 15 in Grays Harbor
(July 25-29. 1977)
24.0
	 Calibrated Model
O Observed Values
20.0
16.0
12.0
12
24
36
46
60
72
84
96
108
120
TIME IN HOURS

-------
Figure 8
Comparison of Simulated and Observed Temperatures
Near UPRR Bridge in Grays Harbor
(July 25-29. 1977)
20.0
19.0
y 18.0
17.0
16.0
15.0
14.0
13 0
Calibrated Modal
O Observed Values
12
24
36
48
60
72
84
96
108
120
TIME IN HOURS

-------
Figure 9
Comparison of Simulated and Observed Dissolved Oxygen
In Grays Harbor at First High
(July 25-29, 1977)
	Calibrated Mod si
O Observed 6-Day Average Values
io r-
X
X
X
X
5	10	15	20
DISTANCE DOWNSTREAM FROM MONTESANO IN NAUTICAL MILES
26

-------
Figure 10
Comparison of Simulated and Observed Dissolved Oxygen
In Grays Harbor at Second Low
(July 25-29, 1977)
	 Calibrated Model
O Observed 5-Day Average Values
10
UJ
&
Z 6
6	10	15	20	25
DISTANCE DOWNSTREAM FROM MONTESANO IN NAUTICAL MILES

-------
While the simulated and observed concentrations of dissolved oxygen
shew good agreement, there remain some questions regarding the
application of the dynamic estuary model to Grays Harbor. Of particular
concern is the discrepancy between simulated and observed valves of
salinity. The magnitude of these differences suggest that there are
still major questions associated with nunerical dispersion, tidal
exchange at the entrance to Grays Harbor and tidal transport within Grays
Harbor.
Another major question which has not been addressed in this study is
that of the impact of photosynthesis and respiration upon dissolved
oxygen. Measurements of primary productivity by the University of
Washington Department of Oceanography (reported in Yearsley (1979))
indicate that the average net production of dissolved oxygen varied from
2	2
0.37 gram 0£/meter /day to 8.13 grams ^/meter /day. If these
values are incorporated into the dissolved oxygen budget, it is likely
that the resulting simulations will predict dissolved oxygen levels
greater than the observed. Yearsley and Hess (1979) suggested that the
primary productivity data may not be representative of _in situ conditions
due to the elevated temperatures at which the samples were incubated.
Before the dynamic estuary model is used for planning purposes in
Grays Harbor these questions must be resolved. Field studies of tidal
transport and oxygen production would be useful, in this respect. Of
equal importance are tests of the nunerical scheme at various freshwater
discharge rates and tidal ranges.
-58-

-------
List of References
Baca, R. G., W. W. Waddel, C. R. Cole, A. Brandstetter, and
0. B. Cearlock, Explore I: A River Basin Water Quality Model, report to
the U.S. Envirormental Protection Agency, Battelle-Northwest,
Richland, Washington, August 1973.
Bennett, J. P. and R. E. Rathbun, Reaeration in Open Channel Flow. U. S.
Department of the Interior, Geological Survey Professional Paper 73 7,
Washington, D.C., 1972, 75 p.
Beverage, J. P. and M. N. Swecker, Estuarine Studies in Upper Grays
Harbor. Washington. U. S. Department of the Interior, Geological Survey
Water Supply Paper 1873-B, Washington, O.C., 196 9, 90 p.
Brogdon, N.0., Jr., Grays Harbor Estuary. Washington; Verification and
Base Tests; Hydraul ic Model Investigation. U.S. Army Engi neer Waterways
Experiment Station, Technical Report H-72-2, Report 1-6, Vicksburg, MS,
Jipril 1972.
Chen, C. W., and G.T. Orlab, Concepts and Utilities of Ecologic Model.
Journal of the Sanitary Engineering Division, ASCE, Volune 96, SA5
October, 1970, pp 1085-1097.
C lei and, B. R., Implementation and Calibration of a Dynamic Estuarial
Water Quality Model on Grays HarborT Masters Thesis, University of
Washington, 1978.
Crim, R. L. and N. L. Lovelace, AUTO-OUAL Model! ing System. Annapolis
Field Office, Region 3, U.S. Environmental Protection Agency, Technical
Report No. 54, March 1973, 301 p.
Eriksen, A. and L. D. Townsend, The Occurrence and Cause of Pollution in
Grays Harbor, Washington State Pollution Control Commission, Pollution
Service, Bulletin 2, 1940, 100 p.
Feigner, K. D. and H. S. Harris, Docimentation Report. FWQA Dynamic
Estuary Model. U.S. Department of the Interior, Federal Water Ouality
Acknini strati on, July 1970, 248 p.
Henderson, F. M., Open Channel Flow, MacMillan Company, Na* York, 1966,
522 p.
Kreizenbeck, R., Grays Harbor Sediment Oxygen Demand, memo to G. L.
O'Neal, Surveillance & Analysis Division, Region 10, U.S. Environmental
Protection Agency, August 1973, 4 p.
Langbein, W. B. And W. H. Qurtm, The Aeration Capacity of Streams, U.S.
Department of the Interior, Geological Survey Cireular 542, Washington,
D.C., 1967, 6 p.
-59-

-------
Lorenzen, M. W., W. W. Waddel, and P. A. Johanson, Development of a
Mathematical Water Quality Model for Grays Harbor and the Chehalis River.
Was hi ngton, report to the U. S. Environmental Protection Agency
Battel! e-Northwest, Richland, Washington, October 1974.
National Ocean Survey, Tidal Bench Marks - Washington, U.S. Departmeit of
Commerce, National Oceanic and Atmospheric Administration, February 1976,
3 p.
National Ocean Survey, Tide Tables 1977 West Coast of North and South
America, U.S. Department of Ccmme*ce, National Oceanic and Atmospheric
Adninistrati on, 1976, 222 pp.
O'Connor, D. J., and W.E. Dobbins, The Mechanism of Reaeration in Natural
Streams, ASCE Trans., Vol. 123, pp. 641-666, 1958
Richardson, D., Drainage-area Data for Western Washington, U.S.
Geological Survey open-file report, 1962, 244 p.
Shubinski, R. P., J. C. McCarty, and M. R. Lindorf, Computer Simulation
of Estuarial Networks, Journal of the Hydraulics Division, ASCE, Volune
91, HY5, September 1965, pp. 33-49
Stefausson, U., and F.A. Richards, Distribution of Dissolved Oxygen,
Density and Nutrients off the Washington and Oregon Coasts, Deep Sea
Research Volune 11, pp. 355-380, 1964.
Tennessee Valley Authority, Heat and Mass Transfer Between a Surface and
the Atmosphere, Division of Water Control Planning, Engineering
Laboratory, Norris, Tennessee, June 1968, 93 p.
U.S. Coast aid Geodetic Survey, Tide Tables: High and Low Water
Predictions, West Coast of North and South America. 1966
Washington State Department of Ecology, Interim Report, Cooperative Grays
Harbor Surveillance Program, 1970.
Washington State Department of Ecology, Interim Report No. 2. Cooperative
Grays Harbor Surveillance Progran. 1971.
Washington State Department of Ecology, unpublished effluent data, 1977.
Water Resources Engineers, Inc., Predictions of Thermal Energy
Distribution in Streams and Reservoirs, prepared for the Department of
Fish and Game, State of California, August 1968, 90 p.
Water Resources Engineers, Inc., Ecologic Modeling of Puqet Sound and
Adjacent Waters, prepared for the U.S. Environmental Protection Agency,
September 1975, 119 p.
Years ley, J. R., A Field Study of the Dissolved Oxygen Resources of Grays
Harbor During Sumner Low Flow, U.S. Environmental Protection Agency,
Region 10, Working Paper EPA 910/B-79-104.
_^n_

-------
Parsley, J.R. and W.C. Hess, Appl i cat ion of a Steady State Dissolved
Oxygen Model to Low Flow Conditions in Gravs Harhnr W^hin^nn, n.<
Environmental Protection Agency, Region 10, Working Paper EPA
910 $-79-105, 1579.
Years ley, J. R. and 0. Houck, A Cluster Analysis of Dissolved Oxygen in
Grays Harbor. Surveillance and Analysis Division, Region 10, U.S.
Environmental Protection Agency, December 1973.
-fil-

-------