SMIG -- Boston Harbor Sewage Outfall Feature Article: September, 1999

Predicting the Physical Effects of Relocating Boston's Sewage Outfall

by Richard P. Signell¹, Harry L. Jenter², and Alan F. Blumberg³

¹ U.S. Geological Survey, Woods Hole, MA, U.S.A.
² U.S. Geological Survey, Reston, VA, U.S.A.
³ HydroQual, Inc., Mahwah, NJ, U.S.A.

Please direct correspondence to:
Richard P. Signell
USGS, Woods Hole Field Center
384 Woods Hole Road
Woods Hole, MA 02543-1598
Internet: rsignell@usgs.gov
Phone: (508) 457-2229
FAX: (508) 457-2309

Editor's note:
This paper was written for a special issue on visualization in marine sciences for the Journal of Estuarine, Coastal and Shelf Science. This special issue of ECSS will be published in the year 2000 and will include a CD containing the visualization data sets. This paper is presented with the permission of the publisher. SPECIAL NOTE: This presentation of the paper is a mirror site for the original which can be found at: http://smig.usgs.gov/SMIG/features_0999/bharbor_inline.html

Some of the "movies" presented in this paper are large files. For your convenience, links to both unzipped and zipped versions of two of the movies are provided. Three of the movies use Indeo 5 compression. To view those movies, you may need to download a driver from: http://developer.intel.com/ial/indeo/video/driver.htm

The movie files are in .avi format. A good player for .avi files (under Windows 95/98/NT) that allows you to step forward and backward one frame at a time may be downloaded from this site [216 Kbytes, zipped].

Citation:
Signell, R.P., Jenter, H.L., and Blumberg, A.F., 2000, "Predicting the Physical Effects of Relocating Boston's Sewage Outfall," J. Estuarine, Coastal and Shelf Science, vol. xx, p xx-xx.

Contents

Abstract

Introduction

Environmental Setting

Boston Harbor

Massachusetts Bay

Numerical Model
Model Results

Model-Data comparison

Tidal Currents

Subtidal Currents

Seasonal Stratification

Effluent Simulations

Winter Effluent Simulations

Summer Effluent Simulations

Effect of Outfall Relocation on Salinity and Velocity Fields

Discussion
Conclusions
Acknowledgements
References

Abstract

Boston is scheduled to cease discharge of sewage effluent in Boston Harbor in Spring 2000 and begin discharge at a site 14 km offshore in Massachusetts Bay in a water depth of about 30 m. The effects of this outfall relocation on effluent dilution, salinity, and circulation are predicted with a three-dimensional hydrodynamic model.

The simulations predict that the new bay outfall will greatly decrease effluent concentrations in Boston Harbor (relative to the harbor outfall) and will not significantly change mean effluent concentrations over most of Massachusetts Bay. With the harbor outfall, previous observations and these simulations show that effluent concentrations exceed 0.5% throughout the harbor, with a harbor wide average of 1-2%. With the bay outfall, effluent concentrations exceed 0.5% only within a few km of the new outfall, and harbor concentrations drop to 0.1-0.2%, a ten-fold reduction. During unstratified winter conditions, the local increase in effluent concentration at the bay outfall site is predicted to exist throughout the water column. During stratified summer conditions, however, effluent released at the seabed rises and is trapped beneath the pycnocline. The local increase in effluent concentration is limited to the lower layer, and as a result, surface layer effluent concentrations in the vicinity of the new outfall site are predicted to decrease (relative to the harbor outfall) during the summer.

Slight changes are predicted for the salinity and circulation fields. Removing the freshwater associated with the effluent discharge in Boston Harbor is predicted to increase the mean salinity of the harbor by 0.5 psu and decrease the mean salinity by 0.10-0.15 psu within 2-3 km of the outfall. Relative to the existing mean flow, the buoyant discharge at the new outfall is predicted to generate density-driven mean currents of 2-4 cm/s that spiral out in a clockwise motion at the surface during winter and at the pycnocline (15-20 m depth) during summer. Compensating counterclockwise currents are predicted to spiral in toward the source at the bottom. Because the scale of the residual current structure induced by the new discharge is comparable to or smaller than typical subtidal water parcel excursions, Lagrangian trajectories will not follow the Eulerian residual flow. Thus, mean currents measured from moorings within 5 km of the bay outfall site will be more useful for model comparison than to indicate net transport pathways.

Introduction

Boston is scheduled to cease sewage discharge from locations within Boston Harbor in Spring 2000, and begin sewage discharge through a 2 km long diffuser at a site 14 km offshore in Massachusetts Bay (Figure 1). Predicting the impact of outfall relocation on the regional marine environment requires understanding the water quality, contaminated sediment, biological metabolism, and a myriad of other processes. Simulations of the physical effects of relocation (the effects on effluent dilution, water currents, and salinity), however, give considerable insight into the nature of the predicted change, and form the basis for modeling of other processes. This paper describes simulations of the existing and future outfall locations using a fully three-dimensional, time-dependent hydrodynamic model that simulates currents, effluent concentration (treated as a conservative tracer), temperature, salinity and vertical and horizontal mixing coefficients.

fig. 1
Figure 1. Massachusetts and Cape Cod Bays, present sewage outfalls in Boston Harbor (triangles), location of new ocean outfall in Massachusetts Bay, the USGS long-term mooring, the National Weather Service buoy (NWS), and the location of the Stellwagen Bank National Marine Sanctuary. The 40- and 80-m depth contours are also shown.

The Bay Outfall location was selected primarily on the basis of 2D depth-averaged modeling studies (Walton et al, 1990). These studies determined that the site would result in dramatically lower effluent levels in Boston Harbor and would lead to greater dilution of effluent compared to other potential inshore sites. Because the receiving waters of Massachusetts Bay are strongly stratified during the summer months and the new outfall diffuser was designed to trap effluent beneath the seasonal pycnocline, there was a clear need for further studies that could test and refine the predictions of the 2D models.

The U.S. Geological Survey (USGS), as part of the National Coastal and Marine Geology Program, has implemented and refined a 3-D hydrodynamic model for Massachusetts Bay to be used in polluted sediment transport studies. Through a cooperative agreement with the Massachusetts Water Resources Authority (MWRA), the USGS has provided hydrodynamic model results that were used for input to a 22 component water quality model (HydroQual and Normandeau, 1995). Model assessment was performed primarily from January 1, 1990 to July 1, 1991, a period when intensive moored and shipboard information were obtained under funding from the Massachusetts Bays Program (Geyer et al, 1992). Simulations with existing and future outfall locations were run and compared over a three-year period (January 1, 1990 to January 1, 1993). This work predicts the influence of the new outfall location on effluent concentrations, salinity and circulation in Boston Harbor and Massachusetts Bay.

Environmental Setting

Boston Harbor

The existing outfalls are located in Boston Harbor, a shallow, small embayment located on the western side of Massachusetts Bay. The surface area of the harbor is 125 km² and the average water depth is 4.9 m. The average tidal range is 2.7 m, which is exchanged through two 15-m-deep passages that connect Boston Harbor to Massachusetts Bay. Currents over most of the harbor are dominated by tides. M₂ tidal current amplitudes of 80-100 cm/s are observed in the passages, and perigean spring tides are about 40% stronger. The vigorous tidal exchange between Boston Harbor and Massachusetts Bay results in flushing times of 3-10 days (Adams et al., 1990; Signell & Butman, 1992). The strength of the tides suffices to keep the passages vertically well mixed throughout the year. About 50% of the 41 m³/s annual mean freshwater input to the harbor is supplied by the sewage outfalls.

Massachusetts Bay

The new outfall is located 14 km offshore in Massachusetts Bay, at a water depth of about 30 m (Figure 1). Massachusetts Bay is a roughly 100 x 50 km semi-enclosed embayment on the western side of the Gulf of Maine. Stellwagen Bank bounds the bay on the east, rising to within 20 m of the sea surface (Movie 1). The average depth of the bay is 35 m, with the western and southern sections of the bay substantially shallower. The southern part of Massachusetts Bay is also referred to as Cape Cod Bay. The deepest water of 80-100 m is found just west of Stellwagen Bank in Stellwagen Basin. There are no major rivers entering directly into Massachusetts Bay. The largest local source of fresh water is the 41 m³/s contributed by the local rivers and sewage outfalls in Boston Harbor. River systems in the Gulf of Maine, within 200 km north of Massachusetts Bay, discharge about 1100 m³/s. Near-surface M₂ tidal current amplitudes are 5-15 cm/s over most of the Bay, but exceed 50 cm/s between Stellwagen Bank and Race Point and in the harbor entrances (Blumberg et al, 1993; Irish and Signell; 1992). M₂ tidal currents of about 10 cm/s at the Bay Outfall result in a tidal excursion of about 2 km.

Screenshot of Movie 1

Movie 1 [1.7 Mbytes, using Indeo 5 compression]
Movie 1. Fly-by visualization of color-shaded relief bathymetry of Massachusetts Bay. The vertical exaggeration is 200:1. Stellwagen Bank, which partially encloses Massachusetts Bay from the Gulf of Maine, is the most prominent feature.

Circulation in Massachusetts Bay is controlled by the large-scale circulation in the Gulf of Maine, local wind forcing, and intrusions of low salinity water associated with Gulf of Maine rivers. The observations of Geyer et al. (1992) reveal a mean annual flow in the surface waters of Massachusetts Bay (Figure 2) that is consistent with earlier drift bottle measurements of Bigelow (1927) and Bumpus & Lauzier (1965). The Maine coastal current flows south with typical speeds of 5-15 cm/s along the Maine and New Hampshire coasts (Normandeau Associates, 1975; Vermersch et al., 1979). When it reaches Cape Ann there is a branch point (Blumberg et al., 1993; Lynch et al., 1997) and much of the flow follows the topography southward past Stellwagen Bank and to the east of Cape Cod. A weaker branch of this current (2-5 cm/s) flows into Massachusetts Bay around Cape Ann, southward along the western shore of Massachusetts Bay, and flows out of the Bay at Race Point (Bumpus, 1973; Geyer et al., 1993). The mean circulation pattern is often altered by wind and runoff events and, as evident by the variability ellipses shown in Figure 2, the sub-tidal fluctuations are typically stronger than the mean, except at stations U2 and RP. The Bay Outfall (near station BB) is in a region of very weak mean flow west of the persistent southward current.

fig. 2
Figure 2. Observed mean flow (arrows) and the sub-tidal variability (ellipses centered on the tips of the mean flow arrows) for near-surface currents (4-8 m below surface) measured between December 1989 and September 1991. The arrows and ellipses have been scaled to correspond to the distance a particle moving with that current would travel in one day. The large gray arrows indicate a conceptual picture of the annual mean circulation in the system.

There is seasonal variation in stratification in Massachusetts Bay, with well-mixed conditions during winter and strong stratification during summer (Geyer et al., 1992). The pycnocline is generally found at 15-20 m depth during the summer, with typical density differences of 3-6 kg/m³. The stratification strongly influences both water properties and dynamics, greatly reducing the vertical exchange between surface and bottom waters and isolating the bottom layer from the direct influence of wind stress and river runoff.

The seasonal variations in stratification, wind stress, and river discharge change the nature of transport and dispersion processes in the bay over the course of the year. During the winter, generally northerly wind events drive basin-scale flows that enhance the counter clockwise circulation, southward flow along the western shore of Massachusetts Bay and northward flow against the wind in the deeper central regions of the bay (Butman, 1975; Brickley, 1994). In the spring, shallow (5-15 m) fresh plumes associated with river runoff in the Gulf of Maine intrude into the bay, commonly generating strong currents of 20-30 cm/s and flow structures with 10-30 km spatial scales (Butman, 1976; Lee, 1992). As the summer progresses, stratification builds, and the wind response of the surface layer is more efficient due to reduced friction. The surface currents at the new outfall location, for example, are strongest during summer, when the wind stress is weakest (Blumberg et al, 1993). The southwesterly winds that frequently occur during the summer result in upwelling along the western and northern coast of the bay, adding additional variability to the density field and density-driven circulation. During September and October, the mean currents along the western shore reverse and flow northward as the result of strong cooling that occurs near the coast (Geyer et al, 1992).

Numerical Model

The hydrodynamic model used for this study is a variant of the model described by Blumberg & Mellor (1987). It solves the primitive equations using finite differences on a curvilinear orthogonal grid in the horizontal plane and discrete sigma levels in the vertical dimension. The Mellor & Yamada (1982) level 2.5 turbulence scheme is used to represent vertical mixing, with modifications by Galperin (1988). This version of the model, called ECOM-si (Estuarine, Coastal and Ocean Model, semi-implicit), contains modifications by Vincenzo Casulli and Ralph Cheng that eliminate the need for separate internal and external mode time steps by solving the external mode implicitly (Blumberg, 1991). Time steps of 414 s were used during the winter and 207 s during the summer (due to stronger summer currents). The curvilinear grid is 68 x 68, resulting in a horizontal spacing ranging from 0.6 km in the vicinity of Boston Harbor to as large as 6 km along the open boundary (Figure 3).

fig. 3
Figure 3. Model grid of Massachusetts and Cape Cod Bays. The curvilinear orthogonal grid allows the mesh resolution to vary spatially, with grid spacing of 600 m in Boston Harbor, 1000 m in the vicinity of the new Bay Outfall site (S3), and up to 6 km at the northeastern open boundary. The grey contour indicates the 40 m isobath. The red transect labeled S1-S4 indicates the location of the vertical sections shown in Figures 7 and 8.

There are twelve vertical layers in the model, with three layers in the top 10% of the water column and the remaining nine layers evenly distributed over the remaining 90% of the water column (runs with increased vertical resolution did not qualitatively change results). The model domain extends well offshore of Massachusetts Bay to facilitate exchange with the Gulf of Maine, and extends to the north to include the Merrimack River, a major source of fresh water to Massachusetts Bay.

Elevation on the open boundaries is specified by a radiation condition with a restorative term that allows the boundary to relax back to specified elevation conditions (Blumberg and Kantha, 1985). In the model runs described here, a very small relaxation time is chosen along the eastern boundary, resulting in essentially a "clamped" condition. The elevation is specified along this boundary as a combination of M₂ tides obtained from a Gulf of Maine tidal model (Lynch & Naimie, 1993) and subtidal sea level fluctuations obtained from a western Gulf of Maine model (Signell et al, 1996). Along the southern open boundary off Cape Cod, a gravity wave radiation condition is applied to allow traveling wave energy to exit from the domain.

Salinity and temperature on the open boundary are specified by advecting interior values out on outflow, and relaxing over a period of several days toward specified values on inflow. The values on inflow are specified by the western Gulf of Maine model along the northern part of the eastern boundary, and by climatology along the rest of the open boundary.

Heat flux and wind stress are specified every 4 hours from meteorological data obtained at the National Weather Service Buoy near the site of the proposed outfall in western Massachusetts Bay (Figure 1). The heat flux components are calculated using techniques described by Weller et al. (1995). Observations are used for all variables except sea surface temperature, which is obtained from the model. Using the model sea surface temperature allows strong spatial gradients, such as those that occur in upwelling areas, to affect the heat flux in a realistic manner. It also provides a feedback mechanism to prevent very shallow water from becoming too warm. This formulation provides improved comparisons with time series temperature measurements, especially in Cape Cod Bay (Signell et al., 1996). Daily river discharges and effluent discharges are specified with data obtained from the USGS and the MWRA.

Sewage effluent is simulated as a conservative tracer, introduced at the grid cells that correspond to the appropriate outfall locations. For the Harbor Outfall simulations, effluent is released at two locations within the harbor corresponding to the main discharges. For the Bay Outfall simulations, the total amount of effluent discharge measured at the Deer and Nut Island plants is released at the offshore location of the new outfall. The effluent is released in the lowest grid cell, but rises due to positive buoyancy. The model grid cell into which the effluent is discharged is approximately 1 km x 1 km, while in reality the effluent will be discharged as a 2 km long line source (through a 55-riser-pipe diffuser system). The relatively coarse grid and hydrostatic formulation of the large scale model does not allow explicit simulation of the small-scale convective mixing processes that occur as less dense water ascends through the water column. Blumberg et al (1996), however, found that the trap height of the plume simulated in the model compared closely with predictions from small-scale buoyant plume models. Zhang et al (1999) attributed the success to three factors: (1) the total dilution is partly due to large scale density exchange flow that the far field model can resolve; (2) the strong pycnocline provides a natural ceiling for the plume; (3) there is beneficial feedback such that if the entrainment is overpredicted, the trap height will be underpredicted and less total dilution will occur. Obtaining the correct trap height is important, since horizontal mixing and transport processes vary greatly with depth during the stratified season.

Model Results

The model was run for the time period corresponding to October 1, 1989-January 1, 1993. Two simulations were performed: the "Harbor Outfall" run, with effluent loading at the outfalls in Boston Harbor and the "Bay Outfall" run, with the same loading at the new offshore location. The flushing time of the bay is several months during winter (Signell et al, 1996), so the first three months of results were considered a spin-up period, and only the results of the three-year period January 1, 1990 to January 1, 1993 were used for analysis.

Model-Data Comparison

Model assessment was performed primarily from January 1, 1990 to July 1, 1991, a period when intensive moored and shipboard information were obtained under funding from the Massachusetts Bays Program (Geyer et al, 1992). Hourly data from the Harbor Outfall run were saved at each of the mooring sites and were analyzed using the same methods as the field data.

Tidal Currents

The simulated tidal currents are close to those observed at the mooring sites, with less than 10% errors in magnitude, and errors of a few degrees in orientation and phase (Figure 4). The tidal currents in the vicinity of the outfall site and throughout much of Massachusetts Bay is 5-15 cm/s, resulting in tidal excursions of 1-3 km extent. These oscillations cause little net advection or dispersion of material in the water column. In Boston Harbor, however, stronger tidal currents interacting with the complex harbor geometry create ebb/flood asymmetries that are an important component of the harbor flushing process (Signell & Butman, 1992). Although the 600-800 m grid spacing in Boston Harbor is not sufficient to reproduce all the strongly nonlinear features observed by Signell & Butman (1992), the simulated flushing time of the harbor was determined to be 4-8 days, comparable to values obtained by observations and higher resolution harbor modeling (Adams et al, 1990; Signell & Butman, 1992).

Figure 4. Comparison of observed and modeled surface barotropic tidal current ellipses.

Subtidal Currents

Once effluent released in the harbor flushes into Massachusetts Bay, or when effluent is released at the Bay Outfall site, subsequent dispersion and transport is determined principally by subtidal currents. To characterize the nature of the subtidal currents, we found it useful to group the data into seasons based on dynamical regimes: Winter (November 1 to March 1) represents a wind-dominated, well-mixed regime; Spring (March 1 to June 1) represents a runoff-influenced, transitional stratification regime; Summer (June 1 to September 1) represent an upwelling, strongly stratified regime; and Fall (Sepember 1 to November 1) represents a transitional overturning regime. Here we discuss the results from unstratified winter conditions (November-February) and stratified summer conditions (June-August) as useful end members, as spring and fall represent transitional stratification periods.

Qualitatively, the model agrees with the observed subtidal circulation during both stratified and unstratified conditions (Figures 5a and 5b). The counterclockwise mean flow through most of the bay, the weak mean flow at the Bay Outfall site, and the structure of the observed current variability in western Massachusetts Bay are all basic features of the system that are represented by the model. Quantitatively, however, the agreement between data and model varies substantially with season and region within the bay. During winter, relatively high correlation (0.6-0.8) between data and model time series of subtidal currents is found in western Massaschusetts Bay, with significantly lower correlation (0.4 or less) at stations U7 in central Cape Cod Bay and at the stations U2, U6 and RP near Stellwagen Bank. During summer, correlation between data and model is lower, never exceeding 0.6. The correlation between data and model currents tends to track the correlation between wind stress and observed current, indicating that the model can represent the timing of locally generated wind events. Due to the semi-enclosed, shallow nature of western Massachusetts Bay and the dominance of wind-driven events during winter, the model-data correlation during winter tends to be high. In the summer and in other regions of the bay, the response is more complex, as remote wind events and density driven events associated with surface trapped river plumes and upwelling play a larger role in the subtidal current response. Small errors in the temporal and spatial phasing of these events destroy the direct correlation between data and model (Signell et al, 1993; Signell et al, 1996). This indicates that the model should not be used to predict the event scale behavior of effluent plumes, except perhaps during strong winter wind events.

Figure 5a. Model-data comparison of mean and subtidal current variability in Massachusetts Bay the Bay Outfall site. (a) Winter upper and lower layer statistics. The conventions for the mean and variability ellipses are the same as in Figure 2. The complex correlation between data and model and between wind stress and data are also indicated at each site. For example, "0.82/0.73" indicates that the correlation between data and model was .82, while the correlation between wind stress and observed current was 0.73.

Figure 5b. Model-data comparison of mean and subtidal current variability in Massachusetts Bay the Bay Outfall site. (b) Summer upper and lower layer statistics. The conventions for the mean and variability ellipses are the same as in Figure 2. The complex correlation between data and model and between wind stress and data are also indicated at each site. For example, "0.82/0.73" indicates that the correlation between data and model was .82, while the correlation between wind stress and observed current was 0.73.

Reproducing the day to day variations of ocean "weather" at specific locations may be required to warrant the use of the model to address some issues, but for other issues it may only be necessary to provide a reasonable representation of the "climate" of the currents over a few weeks or months. If the "climate" of current variability is reasonable, the model should provide a useful prediction of the average plume behavior, or the mean effluent concentration field. In western Massachusetts Bay, the errors in the mean flow are typically less than 3 cm/s and the errors in current fluctuation intensity are less than 40%. More substantial errors in the current climate exist at the other locations. In particular, both the mean and variability in the summer surface flow at stations U2, SC and U7 are underrepresented by the model, suggesting that the counterclockwise flow through the bay was more active than simulated during this season. The error in the mean flow at RP is attributable to a model-data discrepancy in the location of a tide-induced residual eddy in this region.

Seasonal Stratification

The model qualitatively represented the seasonal cycle of stratification in the Bay (Figure 6), an important factor in determining the trapping height of the new outfall plume, as well as its subsequent transport and dispersal. The data and model both show stratification due to salinity starting in April, minimum surface salinities in June, maximum temperatures in August and September and destratification at the end of October. The modeled salinities are about 0.5 psu too fresh in both the upper and lower layers in April, and the modeled temperatures are a couple of degrees too cool in the summer. The modeled and observed density difference between the upper and lower levels (the strength of the stratification) is very similar except during August, when the modeled stratification is about 0.5 kg/m³ too weak, and during September, when the modeled stratification is about 0.5 kg/m³ too strong. Since the level of stratification during August and September suffices to trap the Bay Outfall plume in both the model and in the observations, this was not a significant issue for the effluent dilution modeling simulations.

Figure 6. Model-data comparison of temperature and salinity stratification at the Bay Outfall site: (a) upper and lower layer temperature; (b) upper and lower layer salinity and (c) density difference between upper and lower layers.

Effluent Simulations

Despite some discrepancies between observed and modeled currents and stratification, the agreement was judged sufficient by the authors and the Model Evaluation Group (see acknowledgements) to qualitatively predict the effects of the bay outfall. We focus on describing the predicted changes in effluent concentration, salinity and circulation that will occur when the outfall moves from Boston Harbor to western Massachusetts Bay. We present results in terms of percent effluent, such that a concentration of 1.0% effluent represents a 100 fold dilution of effluent with seawater. Based on long term observations of nutrient distributions in western Massaschusetts Bay (Kelly, 1991), an effluent level of 0.5% (200 fold dilution) represents a threshold above which nutrients introduced with the effluent should be clearly identifiable above background levels. This level is exceeded throughout the harbor with discharge from the harbor outfall, and thus also seems a conservative level to respresent conditions similar to those of the harbor prior to outfall relocation. To facilitate discussion of the comparitive change in effluent concentrations, we refer to the region where the effluent concentration exceeds 0.5% as "the region of significantly elevated effluent levels" or simply as "the plume," although we realize lower or higher concentrations of effluent may be significant depending on the issue of interest. As with the model-data comparisons, we focus on winter and summer simulations as useful end members representing unstratified and stratified conditions.

Winter Effluent Simulations

During the winter, the buoyant effluent discharged at either Harbor or Bay outfall location rises to the sea surface, then is subsequently vertically mixed and transported by tide and wind-driven currents. A comparison of the tidally averaged Harbor Outfall plume and the Bay Outfall plume (whose boundary we define as the 0.5% effluent level) over a two month winter period (January 17-March 17, 1991) illustrates both the rapid fluctuations of the plumes in response to wind-driven currents and the benefits of moving the outfall into deeper, less confined waters (Movie 2). The movie shows particularly interesting sequences February 6-13 and March 10-15 where periods of weak winds allow the plumes to build, and then northerly winds sweep the effluent southward along the coast. This is a mechanism that can carry significantly elevated effluent levels toward Cape Cod. It is also clear throughout the simulation that the areal extent of the region affected by significantly high effluent levels with the Bay Outfall is smaller than with the Harbor Outfall. The region where levels exceed 0.5% always includes all of Boston Harbor and often the coastal region south of Boston for the Harbor Outfall run, while this region is significantly reduced in size and mostly remains offshore in the vicinity of the outfall for the Bay Outfall run.

Screenshot of Movie 2

Movie 2 [2.6 Mbytes, using Indeo 5 compression]
Movie 2. Winter simulation of effluent "plumes" from the Harbor and Bay Outfalls from January 17-March 17, 1991. The plume boundaries are defined as the 0.5% effluent isosurfaces (three-dimensional contours) determined from tidally-averaged results. The region landward of the purple isosurface is the plume from the Harbor Outfall, while the region enclosed by the blue isosorface is the plume from the Bay Outfall. The yellow "pole" shown at the Bay Outfall is not the outfall pipe, but merely a visualization device to give an additional depth cue. The yellow arrow is the wind stress vector with magnitudes reaching 2.0 dyn/cm² during the February 12-13 event and 4.5 dyn/cm² during the March 14-15 event. The vertical exaggeration is 300:1. Isobaths are shown in 20 m intervals.

To more fully assess the predicted change in effluent concentration caused by relocation of the outfall from the harbor to the bay, we compute the mean concentration field over all the winter months of the three-year simulation the Harbor and Bay Outfall scenarios. Comparing horizontal sections of mean winter effluent concentration, we see that when effluent is discharged from the Harbor Outfalls, surface concentrations are relatively high in Boston Harbor, decrease rapidly with distance from the Harbor entrance, and are low throughout most of Massachusetts Bay (Movie 3, upper left panel). The mean effluent concentration has a maximum of about 2% (1 part effluent to 50 parts seawater) in Boston Harbor, and decreases rapidly outside the harbor mouth to a level of about 0.15% at the location of the Bay Outfall. The effluent is further diluted and transported toward the southeast with the mean flow, resulting in effluent concentrations of about 0.125% in Cape Cod Bay (1 part effluent to 800 parts seawater). For the Harbor Outfall simulation, the 0.5% level is exceeded throughout Boston Harbor and extends 10 km to the southeast along the shore. A vertical section from Boston Harbor to Cape Cod shows that the effluent is mostly well mixed from top to bottom, with some slight stratification in western Massachusetts Bay (Figure 7a).

Screenshot of Movie 3

Movie 3, unzipped [711 Kbytes] or zipped [279 Kbytes]
Movie 3. Horizontal sections of mean effluent concentration contrasting the Harbor and Bay Outfalls during Winter (top panel) and Summer (bottom panel). The black contour indicates the 0.5% effluent level. The white contours indicate the 40 m isobath.

fig. 7
Figure 7. Vertical section of mean effluent concentration during winter (November-February) along the transect indicated in Figure 3. Boston Harbor is on the left, Cape Cod Bay is on the right. The red dashed contour indicates the 0.5% effluent level.

When effluent is discharged from the Bay Outfall, the buoyant effluent rises to the surface through the nearly uniform density water column and is then dispersed, vertically mixed, and eventually transported toward the southeast (Movie 3, upper right panel; Figure 7b). Away from western Massachusetts Bay, the levels are similar to the Harbor Outfall simulation. The effluent concentration exceeds 0.5% only within 2-3 km from the outfall, while the effluent concentration is less than 0.25% throughout Boston Harbor. Concentrations in Cape Cod Bay are about 0.125%, similar to the Harbor Outfall simulation.

Summer Effluent Simulations

During the summer, the buoyant effluent discharged in the harbor is vertically mixed in the entrance of the harbor by the strong tidal currents, then is subsequently mixed and transported into the surface waters of Massachusetts Bay. Effluent discharged in Massachusetts Bay, however, becomes trapped below the seasonal pycnocline. A comparison of the tidally averaged Harbor and Bay Outfall plumes over a two month summer period (June 10-August 10, 1990) illustrates the time varying behavior due to currents above and beneath the pycnocline (Movie 4). Particularly noticable are upwelling events caused by southwesterly winds around June 20, June 28, July 21 and August 9. These events bring relatively high effluent levels (>0.5%) to Boston's north shore. It is also clear that, as in winter, the areal extent of the region affected by significantly high effluent levels with the Bay Outfall is smaller than with the Harbor Outfall. The region where levels exceed 0.5% includes all of Boston Harbor for the Harbor Outfall run, while this region is significantly reduced in size and mostly remains offshore and below the pycnocline for the Bay Outfall run (only reaching shore during upwelling events).

Screenshot of Movie 4

Movie 4 [2.5 Mbytes, using Indeo 5 compression]
Movie 4. Summer simulation of effluent "plumes" from the Harbor and Bay Outfalls from June 10-August 10, 1990. The boundary of the regions exceeding 0.5% effluent is indicated by a purple isosurface for the Harbor Outfall and blue isosurface for the Bay Outfall.

Examining the mean effluent concentrations over the summer months, we find that for effluent released from the Harbor Outfalls, concentrations are high in the harbor, decrease rapidly offshore and are low throughout most of Massachusetts Bay (Movie 3, lower left panel). The maximum concentration in the Harbor is about 1.6%, slightly less than the 2% levels found during winter. When the effluent leaves the harbor during the summer, however, it remains in the surface layer, and therefore spreads more effectively across isobaths toward the east. As a result, effluent levels in excess of 0.25% are found in the surface waters directly above the location of the Bay Outfall (Figure 8a).

fig. 8
Figure 8. Vertical section of mean effluent concentration during summer (June-August) along the transect indicated in Figure 3. Boston Harbor is on the left, Cape Cod Bay is on the right. The red dashed contour indicates the 0.5% effluent level.

When effluent is discharged from the Bay Outfall during summer, the buoyant effluent rises to the base of the pycnocline (about 15 m from the surface) and then spreads laterally in the lower layer (Movie 3, lower right panel). Effluent levels at the surface are small, slightly exceeding 0.125% in Boston Harbor. The peak effluent level is about 1.4% at the outfall site, greater than the maximum winter level because stratification confines the effluent to the lower layer. Concentrations in the lower layer decrease rapidly away from the outfall site, with levels exceeding 0.5% again found only within 2-3 km of the outfall. In Cape Cod Bay, the Bay Outfall results in slightly higher effluent levels below the pycnocline, as the Harbor Outfall effluent is mostly trapped in the upper layer (Figure 8b). These effluent levels are still relatively low (<0.125%), however, comparable to winter levels in the lower layer.

Summarizing the winter and summer effluent results, the simulations indicate that the bay outfall location will dramatically reduce effluent concentrations in Boston Harbor and have a small impact on effluent concentrations in most of Massachusetts Bay. During each season, there is a greater than 10-fold reduction in the harbor averaged effluent concentration (Table 1).

Table 1: Harbor-Averaged Effluent Concentrations

Table 1. The average concentration of effluent in Boston Harbor, as a function of season, for the Harbor Outfall and Bay Outfall runs. The reduction factor shows that the Bay Outfall results in more than a 10-fold reduction in the harbor average effluent level throughout the year.

Effect of Outfall Relocation on Salinity and Velocity Fields

The Bay Outfall location represents not only a relocation of sewage effluent, but a relocation of half the freshwater input to the Harbor as well. It is of interest, therefore, to determine the impacts of outfall relocation not only on the effluent concentration fields, but on the salinity and velocity fields as well. By subtracting the Harbor Outfall salinity and velocity fields from the Bay Outfall salinity and velocity fields, we can clearly see the resulting impact of the relocation. During the winter, the Bay Outfall is a buoyant source, and the fresh water reduces the surface salinity at the outfall site by about 0.15 psu (Figure 9b). Although this is a modest change in salinity, this light rising water domes up on the surface and drives a 2-4 cm/s flow that spirals out in a clockwise fashion away from the Bay Outfall site. This influence is chiefly within about 5 km of the outfall. In the lower half of the water column the circulation induced by the outfall is reversed, with flow spiraling in a counterclockwise fashion toward the outfall site (not shown). The local circulation induced by the Bay Outfall is evident in the near surface mean flow, generating 6 cm/s of horizontal shear over less than 5 km. The Bay Outfall enhances eastward mean flow just north of the site, but enhances southward mean flow just south of the site (Figure 9a). In Boston Harbor, the outfall relocation has a much more substantial salinity effect. Salinity increases by an average of 0.5 psu, but there is little impact on Harbor circulation, since the circulation (and exchange) of the Harbor is dominated by nonlinear tidal processes (Signell & Butman, 1992).

Figure 9. (a) Mean winter salinity and velocity field at 2 m depth in western Massachusetts Bay for the Bay Outfall simulation. The concentric circles indicate the location of the Bay Outfall. (b) The difference between the Bay Outfall run and the Harbor Outfall run for the mean winter salinity and velocity field at 2 m depth. With the Bay Outfall, the salinity is reduced 0.14 psu at the Bay Outfall site, and increases from 0.3-0.7 psu in Boston Harbor.

During the summer, effluent from the Bay Outfall induces a circulation at the pycnocline that is similar to the induced near-surface winter circulation (Figure 10b). The influence of this outfall-induced circulation is even more noticeable in the mean flow at the pycnocline (Figure 10a), since the lower layer flow is less strongly driven by surface forcing.

Figure 10. (a) Mean summer salinity and velocity field at 16 m depth in western Massachusetts Bay for the Bay Outfall simulation. The concentric circles indicate the location of the Bay Outfall. (b) The difference between the Bay Outfall run and the Harbor Outfall run for the mean winter salinity and velocity field at 16 m depth. With the Bay Outfall, the salinity is reduced 0.14 psu at the Bay Outfall site.

Discussion

The model simulations predict that moving the outfall will greatly decrease effluent concentrations in Boston Harbor and have little effect on the existing low effluent concentrations that prevail throughout most of Massachusetts Bay. The exception is the region within a few km of the outfall, where locally increased effluent concentrations will occur. When the effluent is discharged from the Harbor Outfalls, the region of high concentration (>0.5%) covers approximately 150 km² (including all of Boston Harbor) and reaches maximum values of about 2% in the Harbor. This is consistent with the findings of Adams et al. (1990), who on the basis of effluent tracer studies determined that the harbor contains approximately 2% effluent. When the effluent is discharged from the Bay Outfall location, the region of relatively high concentration covers less than 10 km². Over Stellwagen Bank and in Cape Cod Bay, the effluent concentrations remain small, less than 0.125%. The reason for this improvement is simply that Boston Harbor is a shallow (4.9 m) confined system, and although the flushing time scale is 3-10 days, effluent levels are relatively high. The Massachusetts Bay outfall site is located in much deeper water (30 m), further from the coast, and allows more effective dilution due to subtidal current fluctuations, despite the weak mean flow.

During the summer stratified months, there is an additional benefit to the Bay Outfall location if the goal is to reduce the impact of anthropogenic nutrient loading. The effluent from the present Harbor Outfalls exits Boston Harbor in the lighter surface layer, discharging nutrients into the photic zone in western Massachusetts Bay. Effluent from the Bay Outfall released at approximately 30 m depth displays dramatically different behavior: it rises to the base of the pycnocline (about 15 m from the surface) and is trapped, spreading out below the pycnocline. This was a design consideration of the new outfall: the waters below the pycnocline have relatively high nutrient levels and are largely aphotic, thus it was determined that the nutrients in the discharged effluent would have less impact on the ecosystem if the plume received enough initial dilution to be trapped in the lower layer. This anticipated benefit is supported by the water quality modeling studies by HydroQual & Normandeau (1995) who used the hydrodynamic results presented herein as input. As the effluent from the new outfall is trapped below the pycnocline, surface layer values are very small. The level of nutrients in the surface layer in western Massachusetts Bay actually decreases with the new outfall location, even directly over the new outfall site.

The new discharge location not only influences the effluent concentrations, but has an effect on the salinity and circulation fields since it represents relocation of 21 m³/s of fresh water, approximately half the quantity that currently enters Boston Harbor. During winter months, the salinity in Boston Harbor is predicted to increase by about 0.5 psu when the outfall is moved offshore, and the near-surface salinity decreases by about 0.15 psu within a few km of the Bay Outfall location. That the salinity effect on Boston Harbor is larger than at the Bay Outfall site is consistent with the effluent dilution results: there is more dilution of fresh water and effluent at the Bay Outfall site. Increasing the salinity of Boston Harbor reduces the offshore gradient in salinity, thereby reducing the small but identifiable estuarine circulation between Boston Harbor and western Massachusetts Bay. Tidal flushing and wind-driven currents, however, dominate Boston Harbor exchange (Signell & Butman, 1992). Thus the reduction in estuarine exchange has little significance. Stronger circulation changes are seen in the vicinity of the new discharge, where the buoyant source creates a local mounding of fresh water. The mean surface currents spiral outward in a clockwise fashion, while the entrained water at the bottom spirals inward in a counterclockwise fashion, showing the influence of the Earth's rotation. The scale of this mean flow is comparable to the typical daily excursions of water parcels driven by wind and tide-driven currents, thus the mean flow pattern is not expected to indicate water parcel pathways (Zimmerman, 1979; Signell & Butman, 1992). Nonetheless, the small scale predicted mean flow modification due to introduction of freshwater at the new outfall would confound efforts to use Eulerian current measurements at the site to determine the general direction of effluent transport.

Having predicted changes in the effluent concentration, salinity and velocity fields due to the new outfall, how much confidence can we have in the predictions? What are the error bars? These are difficult questions to address quantitatively. Certainly there is room for improvement in the match between modeled and observed hydrodynamics. Boundary conditions are a major source of error, and improving them through larger scale simulations, additional measurements, or data assimilation would lead to more accurate circulation in the bay. Bogden et al (1996), for example, showed that the model skill could be improved 50% during certain events by utilizing current meter data within the Bay to improve the boundary conditions. Also, comparison of the transfer function between low-frequency wind stress and 5 m currents in western Massachusetts Bay indicates that the modeled current response at 5 m is 30% stronger than observed, a condition that could possibly be improved by better wind stress determination or by better representation of near-surface mixing processes.

While a closer match to the observed hydrodynamics could plausibly be made with a better model, it is unclear how much these improvements would improve the simulated effluent concentrations. An attempt to simply relate the observed level of subtidal current variability at the Bay Outfall site with simulated plume size, for example, was not successful. The difficulty is that the processes that determine the dispersal of effluent are complex. Movie 5 shows the effluent and velocity fields from the Bay Outfall simulation at 5 m during 45 days of winter. The simulated winter flow field contains many features with spatial scales of order 10 km. It appears that the dispersion of the plume is strongly affected by horizontal shear dispersion operating at this scale. Determining the relationship between the hydrodynamics and the large-scale dispersion processes in this type of complex system remains an important research topic. If this relationship was determined, we might be able to estimate the error in our effluent concentrations predictions based on the model-data mismatch in hydrodynamics. We might also find, however, that our moored array was not sufficient to measure the dispersive characteristics.

Screenshot of Movie 5

Movie 5, unzipped [4.4 Mbytes] or zipped [1.9 Mbytes]
Movie 5. Animation of effluent discharge and simulated currents at 5 m depth during a 45 day period in winter. Observed velocities at 5 m depth from available data are also shown, although the currents at the Bay Outfall site are predicted to change somewhat with the new discharge, as discussed below. The red contour indicates the 0.5% effluent level.

Although a quantitative assessment of error is beyond our current level of understanding, there is reason to have qualitative confidence in the predicted results. As mentioned previously, the model does characterize the basic nature of circulation and stratification in the system and yields predictions of effluent from the harbor outfall that are consistent with observations. A quantitative assessment of the model performance will be conducted after the bay outfall comes on line, utilizing new model runs and data collected by special plume tracking studies.

Conclusions

Analysis of the simulations described herein support the conclusion that the relocation of the outfall to Massachusetts Bay decreases man's overall impact on the marine environment. Comparative effluent dilution simulations over a three year period indicate that relocation of the outfall from Boston Harbor to Massachusetts Bay greatly reduces effluent levels in Boston Harbor, produces mean effluent levels clearly identifiable above background only within 2-3 km of the Bay Outfall, and has little impact elsewhere in the Bay, where effluent concentrations remain low. The greatest salinity changes are predicted in Boston Harbor, yet the impact on circulation and flushing of the Harbor is small. The fresh water input in Massachusetts Bay associated with the Bay Outfall is predicted to induce a complex small-scale mean circulation. The magnitude of this flow is 2-4 cm/s, significantly affecting the mean flow in the region within 5 km of the outfall. Although this circulation will have little impact on the net transport of material, it will prevent mean flows observed near the site from being used to determine the large-scale mean circulation. In recognition of this, the USGS has added an additional current monitoring site 15 km to the south to supplement measurements made at the bay outfall site.

While challenges remain for obtaining more realistic simulations, the model in its present state has been an extremely useful tool for both scientists and managers. Visualization of the model results, particularly in the form of effluent concentration movies, has been an effective means for communicating results to managers and the interested public. Future work will include exploring the relationship between hydrodynamics and plume dispersion and testing the model predictions once the bay outfall comes on line.

Acknowledgements

The Massachusetts Water Resources Authority and the U.S. Geological Survey cooperatively funded this work. Special thanks to the MWRA Model Evaluation Group for their stimulating discussion and suggestions: Robert Beardsley, E. Eric Adams, Donald Harleman, Anne Giblin, John R. Kelly, John E. O'Reilly and John F. Paul.

References

Adams, E. E., Kossik, R. F., Cosler, D. J., Mac Farlane, J. K. & Gschwend, P. M., 1990, Calibration of a transport model using halocarbons, In Estuarine and Coastal Modeling I (Malcolm Spaulding, editor), ASCE, New York. 380-389.

Blumberg, A., 1991, A Primer for ECOM-si, Technical Report, HydroQual, Inc., 66 pp.

Blumberg, A., Ji, Z., & Ziegler, C. K., 1996, Modeling outfall plume behavior using a far field model, J. Hydraulic Engineering, 112, 610-616.

Blumberg, A., Signell, R.P., & Jenter, H. 1993, Modeling transport processes in the coastal ccean, J. Environmental Engineering, 1, 31-52.

Blumberg, A. F. & Kantha L. H., 1985, Open boundary condition for circulation models, J. Hydraulic Engineering, 111, 255-273.

Blumberg, A. F. & G. L. Mellor, 1987, A description of a three-dimensional coastal model, pp. 1-16, Coastal and Estuarine Sciences, Vol. 4. AGU, Washington, DC.