James J. Fitzpatrick
Dominic M. Di Toro
HydroQual, Inc.
1 Lethbridge Plaza
Mahwah, NJ 07430


The pollution of the Great Lakes from municipal and industrial sources was perceived to be a serious problem as early as the 1950's. This perception was even greater for Lake Erie, which during the 1960's and 1970's was known as the "Dead Sea of North America." Historically, a large part of Lake Erie's water quality problems appeared to be related to eutrophication. During the summer months, windrows of Cladophora glomerata covered large portions of the lake's shoreline. The surface waters of the western basin and parts of the central basin of Lake Erie were populated with Aphanizomenon flosaquae, giving the impression that someone had poured green paint on the water surface. The subsequent settling and decay of this and other algal biomass resulted in widespread hypoxia and anoxia in the bottom waters of the lake, further stressing commercial and sport fisheries, which were already suffering from the stresses of overfishing. While part of the eutrophication problems of the Great Lakes was due to increased phosphorus loadings associated with population growth, the post-World War II development of phosphorus-based detergents was also identified as a significant contributor to increased phosphorus loadings to the lakes. In its 1969 report to the United States and Canada, the International Joint Commission (IJC) concluded that phosphorus enrichment had accelerated the eutrophication of Lake Ontario and had resulted in a condition of advanced eutrophication in the western basin of Lake Erie. The IJC offered a number of recommendations to deal with the problem of eutrophication, including:

Subsequent to the 1969 IJC, the United States and Canada signed the Great Lakes Water Quality Agreement, which charted a course of action for the two governments in dealing with Great Lakes water quality problems. That Agreement, together with the 1972 Amendments of the Federal Water Pollution Control Act, which called for research and technical development with respect to Great Lakes water quality, including an analysis of present and projected future water quality under varying conditions of wastewater treatment and waste disposal, provided the impetus for water quality modeling of the Great Lakes.

The purpose of this paper is to provide a history of the development and application of mathematical models of the eutrophication processes in Lake Erie and to comment on the potential use of eutrophication models in an ecosystem approach to address Lake Erie water quality problems. However, we would be remiss if we did not include in this history the research and development efforts conducted in the other lakes of the Great Lakes system.


Initial Efforts Leading to LAKE 1

One of the first eutrophication modeling studies of Lake Erie was reported by Di Toro et al. (1973). This work, which focused on Western Lake Erie, grew out work performed by Hydroscience (1973) for the Great Lake Basin Commission. The original Hydroscience study provided an assessment of the feasibility of applying a limnological systems analysis to the water resources problems of the Great Lakes. During the study, Hydroscience developed and calibrated a demonstration water quality model for western basin of Lake Erie. The demonstration model included chlorides, coliform bacteria, eutrophication and a food chain model of cadmium.

The demonstration eutrophication model used seven vertically integrated segments and a steady-state circulation pattern (Figure 1a) to represent the physical and transport features of the western basin of Lake Erie. The estimates of the circulation pattern in western Lake Erie were based on hydrodynamic model computations developed by Gedney (1971) and on observed current information reported by Herdendorf (1967). The kinetic framework (Figure 1b) employed for the eutrophication model incorporated eight state-variables (phytoplankton chlorophyll, herbivorous zooplankton carbon, carnivorous zooplankton carbon, organic and inorganic phosphorus, organic nitrogen, ammonia nitrogen and nitrate nitrogen) and followed the phytoplankton-zooplankton-nutrient model first structured by Di Toro et al. (1971) for the Sacramento-San Joaquin Delta. Using water quality data collected by the Canadian Centre for Inland Waters (CCIW) and the USEPA, the eutrophication model was calibrated for the year 1970. Initial calibration results were encouraging. The magnitudes and shapes of the calculated curves for phytoplankton biomass, as indicated by chlorophyll-a (Figure 2a), inorganic phosphorus (Figure 2b), ammonia nitrogen (Figure 3a) and nitrate nitrogen (Figure 3b) were in reasonable agreement with the observed data, although some systematic deviations were observed. The calibrated model was also used to hindcast lake water quality for the year 1930. Estimates of chlorophyll-a, used to validate the 1930 hindcast calculations, were based on observations of algal cells counts collected by the U.S. Department of Interior, Fish and Wildlife Service for the period 1928-30 (Wright et al. 1955) and a relationship between chlorophyll-a and total algal cell counts (Figure 4a) based on EPA Western Lake Erie surveillance data for the period 1967-1968. Model computations of chlorophyll-a for the 1930 hindcast compared reasonably to the observed data (Figure 4b) and provided an encouraging test of the Lake Erie western basin eutrophication model.

Figure 1a Steady state transport for seven compartment western Lake Erie Model

Figure 1b Kinetic pathways of the endogenous variables

Figure 2a Chlorophyll verification comparison of model results and observed data

Figure 2b Inorganic phosphorus verification comparison of model results and observed data

Figure 3a Ammonia nitrogen verification comparison of model results and observed data

Figure 3b Nitrate nitrogen verification comparison of model results and observed data

Figure 4a Relationship between total algal cell counts and total chlorophyll measurements in western Lake Erie

Figure 4b Chlorophyll hindcast to 1930 comparison of model results and observed data

It was also during this time that the generalized water quality modeling code known as WASP (Water Analysis Simulation Program) came into being. The earlier Sacramento-San Joaquin Delta model developed by Di Toro et al. (1971) used an IBM 1130 version of the Continuous System Modeling Program (CSMP), a continuous system simulation language (cf. Wang and Chang 1976). While CSMP was very effective for "bread-boarding" models, it was not a practical tool for running multi-segment water quality models; hence, leading to the development of WASP. The WASP computer code was eventually documented (Di Toro et al. 1981) and transferred to the USEPA Large Lake Research Station, Grosse Ile, Michigan. Subsequent to this the USEPA Center for Exposure Assessment Modeling (CEAM), Athens, Georgia provided support for the model, eventually updating the eutrophication kinetics based on a eutrophication modeling study of the Potomac River Estuary (Thomann and Fitzpatrick 1982) and adding a riverine-based hydrodynamic model, DYNHYD, to the modeling package (Ambrose et al. 1988). The WASP computer code was also used as a platform for the development of WASTOX, which models the fate and transport of toxic substances (Connolly and Winfield 1984). The WASP modeling package is still supported by the USEPA (Ambrose et al. 1993) and is still being used in many water quality studies within the United States.

Lake Ontario was the next of the Great Lakes wherein eutrophication modeling was applied. This work was performed under an EPA grant to Manhattan College and resulted in a two-volume EPA Ecological Research Series Report (Thomann et al. 1975, Thomann et al. 1976). The major feature of the Lake Ontario eutrophication model (which used the eutrophication kinetic structure known as LAKE 1) that differed from the Lake Erie model was the incorporation of vertical segmentation to represent the epilimnion and hypolimnion of the lake. A relatively simple two layer version of Lake Ontario, which assumed horizontal homogeneity, was calibrated to a four year (1966-1969) CCIW data set and a satisfactory calibration was achieved. Thomann et al. (1975) also developed a preliminary calibration of a 67 segment three-dimensional representation of the lake. The 67 segments were distributed over five vertical layers and included a "ring" of segments, extending 10 km out from the shoreline, to represent the near-shore environment. The model used vertical casts of seasonal water column temperature to calibrate the vertical mixing coefficients between adjoining vertical segments. While initial calibration results were encouraging, the authors reported that, "The size of the LAKE 3 model makes it difficult to fully comprehend the output [from the model]..." and subsequent application of the LAKE 3 model was limited. Thomann et al. (1976) also performed some long-term simulations and noted that it would require a number of years for lake phytoplankton biomass and primary productivity to respond to reductions in phosphorus inputs to the lake.


The next phase in the development of Lake Erie eutrophication models (Di Toro and Connolly 1980) resulted in expansions in both the number of segments used to represent the spatial geometry of the lake (Figure 5) and the number of state-variables used in the kinetic framework (Figure 6). The revised model segmentation divided the lake into five segments, three to represent the epilimnion of the western, central and eastern basins of the lake and two segments to represent the hypolimnion of the central and eastern basins of the lake.

Figure 5 Lake Erie Model Segmentation of Western, Central and Eastern Basins: Water Segments 1-6, Epilimnion (top) and Hypolimnion (middle); Sediment Segments 7-10, (bottom).

Figure 6a Lake Erie State Variable Interactions. Representation of nitrogen (top> and silica (bottom) nutrient cycles.

Figure 6b Lake Erie State Variable Interactions. Representation of Dissolved Oxygen (top) and Phosphorus (bottom) Nutrient Cycles.

The number of state-variables was expanded from 8 to 15 and included: diatom and non-diatom chlorophyll-a, herbivorous and carnivorous zooplankton, detrital organic nitrogen, ammonia nitrogen, nitrate nitrogen, unavailable phosphorus and soluble reactive phosphorus, unavailable silica, soluble reactive silica, detrital organic carbon, dissolved inorganic carbon, alkalinity and dissolved oxygen. The kinetic framework of the model was expanded in part to better represent the observed seasonal patterns in phytoplankton growth (diatoms and non-diatoms) and the observation that silica appeared to be the nutrient limiting spring phytoplankton growth. The other reason for expanding the kinetic framework of the model by adding total inorganic carbon and alkalinity state-variables was to reduce the degrees of freedom in model calibration. Since a major concern of the authors was to ensure that the reactions involving oxygen were correctly incorporated and since carbon dioxide (inorganic carbon), as well as oxygen, is produced or consumed as a consequence of primary production, algal respiration and oxidation of organic carbon, by including dissolved inorganic carbon and alkalinity, it was possible to check model computations against observations of pH. Furthermore, since alkalinity production and destruction could be calculated for each of the relevant processes in the kinetic framework, using appropriate stoichiometric ratios, adding dissolved inorganic carbon and alkalinity to the calculation actually decreased the degrees of freedom, since there was no increase in the number of constants used during model calibration and since more data were available for model-data comparisons.

An additional feature of the LAKE 1A kinetic framework was the inclusion of a preliminary model of sediment oxygen demand (SOD) and sediment nutrient flux. The SOD/nutrient flux model was incorporated because it was realized that interactions between lake waters and sediments could have a profound effect on the concentrations of oxygen and nutrients in a comparatively shallow lake such as Lake Erie. Analysis of observed nutrient fluxes and measurements of SOD indicated that areal fluxes from the sediments of Lake Erie were substantial sources to the water column based on a volumetric basis. Additionally, it was observed that the occurrence of hypoxia and anoxia dramatically increased certain nutrient flux rates to the water column. The preliminary SOD/nutrient flux model utilized a single sediment layer and included a one-dimensional mass transport equation for the concentrations of dissolved substances (ammonia, CO2, and the oxygen-demanding end-products of organic carbon decomposition, i.e., sediment oxygen demand) and first-order decay for detrital organic nitrogen and organic carbon. At the time of their study, the authors believed that the complete analysis of phosphorus and silica fluxes would require a rather elaborate computation of solution-precipitate chemistry, which was beyond the scope of the study. Instead they used an empirical approach which relied on observed interstitial water concentrations of phosphate and silica. During periods of bottom water aerobic conditions, it was assumed that phosphorus and silica fluxes did not occur. During periods of bottom water anaerobic conditions, the diffusive exchange was set to the same value used for ammonia, CO2, and oxygen. Phosphorus and silica fluxes were then computed using this diffusion coefficient and gradients in concentration between bottom overlying water and observed interstitial concentrations of phosphorus and silica. This preliminary SOD/nutrient flux model subsequently provided the basis for a state-of-the-science model of sediment nutrient composition and nutrient flux (Di Toro and Fitzpatrick 1993), currently employed in a modern eutrophication model of the Chesapeake Bay system (Cerco and Cole 1993).

The LAKE 1A model was applied to Lake Erie (Di Toro and Connolly 1980) and, in general, the model reproduced most of the seasonal and spatial features of the observed data (Figure 7), Including:

Figure 7a Lake Erie Calibration Results, 1970. Eplimnion Western Basin (lefthand side); Central Basin (center); Eastern Basin (righthand side); chlorophyll a, µg/L, (top); C14 shipboard primary production mg C/m3/hr (middle); total phosphorus mg/L (bottom). Symbols: mean + standard deviation.

Figure 7b Lake Erie Calibration Results, 1970. Eplimnion Western Basin (lefthand side); Central Basin (center); Eastern Basin (righthand side); orhophosphorus (top); nitrate nitrogen (middle); reactive silica (bottom).

At the same time that Di Toro and Connolly were developing and calibrating the LAKE 1A model to Lake Erie, other researchers and water quality managers were investigating alternate model formulations for evaluating the interrelationships between phosphorus inputs and water quality in Lake Erie and the other Great Lakes. Chapra (1977, 1980) applied a general total phosphorus budget mass balance equation (shown in Equation 1) to each of the Great Lakes,

where V is the segment volume, p is the mean annual total phosphorus concentration of the segment, t is time, W is the mass loading rate of total phosphorus to the segment, Qb is the advective flow from an adjacent segment, pb is the total phosphorus concentration of the adjacent cell, Q is the advective flow leaving the segment, E is a bulk diffusion coefficient, v is the apparent settling velocity, As is the segment surface area, and F is the rate of feedback of phosphorus from the sediments of the lake. The purpose of the phosphorus budget model was to transform loadings into in-lake concentrations of total phosphorus. These concentrations were, in turn, related to other trophic variables, such as chlorophyll-a, primary production, and secchi depth (Figure 8a) using statistical correlation analysis. Chapra (1980) then developed a hypolimnetic oxygen model for Lake Erie and using observed epilimnetic and hypolimnetic temperature and oxygen data (Figure 8b) found the following relationship for the lake:

where DOmin is the minimum concentration of hypolimnetic oxygen just before turnover and Pr is the mean annual surface production, which could be estimated from the relationships shown in Figure 8a. Chapra (1980) then utilized this model to develop system response matrices for each Great Lakes model segment.

Figure 8a Correlations between surface water quality variables

Figure 8b Temperature and oxygen concentration for central Lake Erie from 1967 through 1972 for the summer stratified period (June through September). Triangels designate surface data (1m depth) and x's designate bottom data (>=22m depth).

Vollenweider et al. (1980) reported on the formulation of a phosphorus loading plot model (Vollenweider 1966) and its application by Rast and Lee (1978) to the Great Lakes system. The phosphorus loading model was then used to predict changes in the average chlorophyll-a concentrations and secchi depths that would result in each of the Great Lakes from the implementation of proposed phosphorus loading objectives contained in the 1978 Great Lakes Water Quality Agreement. For Lake Erie, the Vollenweider model projected a 33 percent reduction in average chlorophyll concentrations (from 5 ug/L to 3.3 ug/L).

Bierman (1980) made a comparison between the Lake Erie eutrophication models developed by Di Toro (Di Toro and Connolly 1980), Chapra (1980) and Vollenweider (Vollenweider et al. 1980) and found that model computations of projected phosphorus concentrations in response to changes in phosphorus loading to the lake were in reasonably good agreement with one another. In general, agreement among the model results for projected chlorophyll-a concentrations was not as good as for phosphorus. The Di Toro model computed significantly higher chlorophyll-a levels in the western basin than either the Chapra or Vollenweider models. However, there was better agreement in the central and eastern basins. The most important water quality indicator for Lake Erie was dissolved oxygen concentration in the hypolimnion of the central basin. It was here that the greatest divergence between the models was observed. The Di Toro model projected dissolved oxygen concentrations that were 2 mg/L higher than those projected by the Vollenweider model in the loading range of 8,000 to 10,000 metric ton/yr. The Di Toro and Chapra models were in reasonable agreement for phosphorus loads less than 12,000 metric ton/yr. However, results for the two models progressively diverged as phosphorus loads increase above 12,000 metric ton/yr. In part, this may have been due to the fact that the Vollenweider and Chapra models assumed a 3.3 meter hypolimnion depth at the end of summer stratification, while the Di Toro model represented average concentrations below a depth of 17 meters in the central basin at the end of summer stratification. In addition Bierman concluded that computations from the Di Toro model were more representative of volumetric conditions, while the Vollenweider and Chapra models were more representative of areal conditions.

The LAKE 1A modeling framework was also applied to Lake Ontario and a modeling analysis of 10 years of Lake Ontario data was reported by Thomann and Segna (1980). Thomann and Segna (1980) also used statistical analysis to help quantify the calibration of the model to the observed data and further compared on a statistical basis the improvement in model calibration provided by the LAKE 1A kinetic formulation versus the LAKE 1 kinetics. In summary, Thomann and Segna (1980) found model verification scores, as judged by a difference of means test for a pooled set of six state-variables (chlorophyll-a, total phosphorus, dissolved orthophosphorus, ammonia, nitrate and dissolved silica), improved from 15 to 20 percent (Figure 9a) when using the LAKE 1A kinetics versus the LAKE 1 kinetics. In addition, use of the LAKE 1A kinetics reduced the median relative error for chlorophyll-a from 42 percent to 30 percent. In general, relative error scores were reduced significantly for most variables when LAKE1A kinetics were employed in the eutrophication analysis (Figure 9b).

Figure 9a Verification scores, difference of means test, pooled six-state variables, 1967-1976

Figure 9b Relative error scores, siz state variables and across all variables.

During this time, researchers from the USEPA Large Lakes Research Station began to implement similar versions of the Di Toro and Thomann eutrophication models. In particular, Bierman et al. (1980) developed and calibrated a multi-class internal nutrient pool phytoplankton model for Saginaw Bay in Lake Huron. The Bierman model differed from the Di Toro and Thomann models in that phytoplankton biomass was partitioned into five functional groups: diatoms, greens, N2-fixing blue-greens, non-N2-fixing blue-greens and "others." In addition, a more detailed kinetic representation of phosphorus and nitrogen dynamics, which included internal nutrient pools within the phytoplankton was used.


The use of mathematical models in the analysis of Great Lakes eutrophication served two purposes: to synthesize the available information concerning the major physical, chemical and biological features that influenced eutrophication into a comprehensive framework; and to apply that framework to make projections of lake response to various changes, primarily reductions in nutrient inputs to the lakes. The evaluation the model performance in predicting water quality response to changes in nutrient inputs can only be performed after these reductions occur. The opportunity for such a post audit analysis was provided in Lake Erie during the period 1970 to 1980. During this period, phosphorus inputs to the lake declined from levels of 20,000 to 25,500 metric tons/yr in the early 1970's to approximately 13,500 metric tons/yr in 1980. Di Toro et al. (1987) reported on the results of a 10-year simulation from 1970 to 1980 using measured lake loadings. An examination of the model's long time scale predictive capability indicated that the model was able to reproduce some of the observed features of improved water quality resulting from total phosphorus reductions in the 1970-1980 decade (Figures 10 through 13). In particular, the model was able to predict the observed decrease in anoxia area of the central basin (Figure 14). However, the results also illustrated that short-term calibrations (e.g., 1 year) failed to capture long-term behavior of certain variables (e.g., nitrate nitrogen, as shown in Figure 15) if a small but significant source or sink (e.g., sediment denitrification) was not well calibrated in short-term computations.

Figure 10 Comparison of model predicted and 1970 to 1980 observed cruise mean total phosphorus data.

Figure 11 Comparison of model predicted and 1970 to 1980 observed cruise mean soluble reactive phosphorus

Figure 12 Comparison of model predicted and 1970 to 1980 observed cruise mean chlorophyll a - western central, and eastern basins of Lake Erie

Figure 13 Comparison of model predicted and 1970 to 1980 observed yearly dissolved oxygen statistics - western, central, and eastern basins of Lake Erie

Figure 14 Comparison of model predicted and 1970 to 1980 observed yearly total percent area anoxia statistics - central basin of Lake Erie

Figure 15 Comparison of model predicted and 1970 to 1980 observed yearly mean nitrite plus nitrate nitrogen statistics - western, central, and eastern basins of Lake Erie.

Bierman and Dolan (1986a) reported on additional calibration efforts conducted in Saginaw Bay in Lake Huron. One of their conclusions was that wind-induced sediment resuspension was an important mechanism for re-introducing phosphorus into the water column. In the calibrated model, the resuspension mechanism was found to account for 36% and 68% of the computed spring and fall average total phosphorus concentrations, respectively. Bierman and Dolan (1986b) also conducted a post audit of the Saginaw Bay eutrophication model. They compared a priori model predictions to an extensive set of survey data collected in 1980 and found that while the response of the bay was consistent with trends of model predictions, it was not consistent with their absolute values in all cases. In particular, observations were consistent with model predictions that threshold odor violations in the municipal water supply would be eliminated with total phosphorus tributary loadings that occurred in 1980. Model predictions were also consistent with the observation that in-bay total phosphorus concentrations decreased in a much smaller proportion than the decrease in tributary loadings, and that phytoplankton biomass, especially blue-greens, decreased by larger proportions that total phosphorus concentrations. It was hypothesized that discrepancies between model predictions and observations for total phosphorus and phytoplankton were due to wind-resuspension and changes in phytoplankton community structure, respectively.


One of the last applications of a eutrophication model to Lake Erie that the authors are aware of was performed by Blumberg and Di Toro (1980). In this application, Blumberg and Di Toro applied a coupled hydrodynamic and water quality model to examine the response of dissolved oxygen concentrations to warming of the central basin of Lake Erie. A synthesis of results from the coupled model, forced by climate warming scenarios from three atmospheric general circulation models, suggested that there would be a substantial decline in oxygen concentrations in the central basin. Losses of 1 mg/L of dissolved oxygen in the upper layers and of 1-2 mg/L in the lower layers of the central basin could be expected, along with an increase in the area of the lake that is anoxic. The modeled decline was predicted to be due to warmer lake temperatures, which would increase the rates of bacterial activity in the hypolimnion waters and sediment, rather than due to thermocline location and volume of water below the thermocline. While the 1 mg/L decline in the epilimnion would not greatly affect fish life, the projected declines to 3 mg/L or less in the hypolimnion could pose a threat to adult "coldwater" fish life.


In many ways the current state-of-the-science in eutrophication modeling has not changed appreciably from the models of Di Toro, Thomann and Bierman. The basic linkages between nutrients, phytoplankton biomass, primary production, nutrient recycle and dissolved oxygen in today's eutrophication models are quite similar to the Great Lakes' models. Most modern eutrophication models still use the Monod theory for which algal growth rates depend upon external nutrient concentration rather than the formulations in which growth rate depends on the internal cellular nutrient concentration (e.g., Bierman 1980). However, one recent modification to the structure of algal growth dynamics is the inclusion of variable nutrient stoichiometry for phytoplankton biomass. In their model application to Chesapeake Bay, Cerco and Cole (1993) included a function which permitted variation in the carbon to phosphorus (C:P) ratio for phytoplankton biomass. This empirical function was developed from observed C:P ratios and soluble reactive phosphorus data in the upper bay. HydroQual (1995) recently implemented a more process-based function for variable carbon to nutrient stoichiometry, which also included variable carbon to chlorophyll stoichiometry, in a modeling study of the Massachusetts Bays system. This model formulation was based on a model developed by Laws and Chalup (1990). A unique feature of this model is that it accounts for variable carbon to chlorophyll ratios in phytoplankton due to both light and nutrient status.

Perhaps the two most significant changes that have been incorporated into modern eutrophication models are: (1) the addition of a coupled water column/sediment nutrient flux submodel and (2) the direct coupling of time-variable three-dimensional hydrodynamic and water quality models. The sediment nutrient flux model framework accounts for the deposition of organic matter from the water column to the sediment bed of the water body, its subsequent diagenesis or decomposition, and the flux of resulting end-products back to the overlying water column. The model includes sediment processes for temperature and oxygen-dependent nitrification-denitrification, sorption of dissolved inorganic phosphorus and dissolved inorganic silica to sediment solids and the sorption of dissolved inorganic phosphorus to iron hydroxides in the aerobic layer of the sediment. The model also considers the generation of sediment oxygen demand, hydrogen sulfide and methane from the reduction of organic matter. The model was developed and calibrated in Chesapeake Bay to an extensive multi-year data set (Di Toro and Fitzpatrick 1993). The model has been further verified against a wide-range of nutrient conditions using an extensive nutrient flux data set obtained from the University of Rhode Island MERL mesocosms. The model successfully reproduced the observed sediment nutrient composition and nutrient flux data using essentially the same parameter set as was used for Chesapeake Bay. The only parameters to be changed between the two calibrations were the temperature-correction coefficients associated with diagenesis and the aerobic/anaerobic partition coefficients for phosphorus sorption. The MERL data were collected more frequently than the Chesapeake Bay data and, therefore, included more "cold weather" measurements than did the Chesapeake Bay data set. These additional data required minor adjustments to these temperature coefficients in order to reproduce the "cold weather" data. It also appeared reasonable that the phosphorus partition coefficients might be different between the two systems given differences in the iron content of the sediments in Chesapeake Bay and the MERL mesocosms (Narragansett Bay).

Examples of the coupling of water quality models to high resolution time-variable three-dimensional hydrodynamic models include the Chesapeake Bay system (Cerco and Cole 1993), the Massachusetts Bays system (HydroQual 1995, Blumberg and Fitzpatrick 1999). The use of well calibrated hydrodynamic models to drive water quality model computations is that it removes a degree of freedom in model calibration, i.e., the external specification of advective and dispersive transport. In addition, the use of high resolution hydrodynamic/water quality models removes, to a large degree, the "numerical errors" or numerical dispersion associated with coarse-grid or box water quality models.


While to a large degree the issue of eutrophication has been addressed in Lake Erie and signs of improvement in water quality and ecosystem health have been observed (Bertram 1993, Makarewicz 1993, Schloesser et al. 1995, Krieger et al. 1996), there still remain a number of environmental concerns within the lake. These include:

  1. primary production and total algal biomass,
  2. blue-green algal biomass,
  3. walleye production,
  4. invasion by non-native species (e.g., zebra mussels),
  5. fish body burdens of bioaccumulative chemicals (e.g., PCBs),
  6. richness and evenness of fish community trophic levels between algae and top predator fish.

At present it would seem that modern state-of-the-science eutrophication models could play an important role in helping to address some of these water quality issues through an "ecosystem approach." Certainly eutrophication models address the first two issues directly. Consideration should be given to updating the existing Lake Erie eutrophication water quality models to include the sediment nutrient flux model and coupling to time-variable three-dimensional hydrodynamic models of the lake. Primary production of organic carbon is also essential to understanding and modeling of the fate and transport of toxic materials in the lake, since many of these materials sorb to dissolved and particulate matter. The coupling of the eutrophication and sediment nutrient flux model may also help to understand the trapping and bio-availability of toxic metals in sediments, since sulfide is a state-variable in the sediment nutrient flux submodel. Some zebra mussel modeling has been conducted by LimnoTech, Inc. and may provide a useful starting point for integration into a eutrophication/zebra mussel ecosystem modeling package. HydroQual, Inc. has recently completed work on the development of a suspension and deposit feeder model for the Chesapeake Bay system. This submodel, which is linked to the water column/nutrient flux model of the Bay, has also successfully been applied to Jamaica Bay, New York with minimum changes in model parameters.

We believe that coupling eutrophication models to fisheries models becomes more speculative, although some efforts have been initiated in Chesapeake Bay and the Great Lakes community. While coupled hydrodynamic/water quality models of eutrophication can provide information (e.g., water temperature, available food, dissolved oxygen, etc.) to fisheries models, the projections of walleye production and fish community structure and diversity are made more complicated by the interactions between various predator/prey fish species, the impacts of overfishing and perhaps the relatively long life-spans of the fish themselves. Perhaps the area requiring the greatest research and monitoring effort is the area of the fishery.


Ambrose, R.B. et al., 1988. WASP4, A hydrodynamic and water quality model--Model theory, user's manual and programmer's guide. USEPA, Athens, GA. EPA/600/3-87-039.

Ambrose, R.B., T.A. Wool, J.L. Martin, 1993. The Water Quality Analysis Simulation Program, WASP5. Part A: Model documentation. USEPA ERL, Athens, GA 30613.

Bertram, P.E., 1993. Total phosphorus and dissolved oxygen trends in the central basin of Lake Erie, 1970-1991. J. Great Lakes Res. 19:224-236.

Bierman, V.J., Jr., 1980. A comparison of models developed for phosphorus management in the Great Lakes. Phosphorus Management Strategies for Lakes. Ed. R.C. Loehr, C.S. Martin, W. Rast. Ann Arbor Science Publishers Inc., Ann Arbor, MI.

Bierman, V.J., Jr., D.M. Dolan, E.F. Stoermer, J.E. Gannon and V.E. Smith, 1980. The development and calibration of a spatially simplified, multi-class phytoplankton model for Saginaw Bay, Lake Huron. Contribution No. 33, Great Lakes Environmental Planning Study, Great Lakes Basin Commission, Ann Arbor, MI.

Bierman, V.J., Jr., and D.M. Dolan, 1986a. Modeling of phytoplankton in Saginaw Bay: I. Calibration Phase. J. Environ. Engr. ASCE. 112:400-414.

Bierman, V.J., Jr., and D.M. Dolan, 1986b. Modeling of phytoplankton in Saginaw Bay: II. Post-Audit Phase. J. Environ. Engr. ASCE. 112:400-414.

Blumberg, A.F. and D.M. Di Toro, 1990. Effects of climate warming on dissolved oxygen concentrations in Lake Erie. Trans. Amer. Fisheries Soc. 119:210-223.

Blumberg, A.F. and J.J. Fitzpatrick, 1999. Invited lecture: Strategy for coupling hydrodynamic and water quality models for addressing ling time scale environmental impacts. In Environmental Hydraulics. Ed. J.H.W. Lee, A.W. Jayawardena, and Z.Y. Wang. Proceedings of the second international symposium on environmental hydraulics. Hong Kong, China. A.A. Balkema, Rotterdam, Netherlands.

Cerco, C.F., and T. Cole, 1993. Application of the three-dimensional eutrophication model CE-QUAL-ICM to Chesapeake Bay. U.S. Army Engineer Waterways Experiment Station, Vicksburg, MS.

Chapra, S.C., 1977. Total phosphorus model for the Great Lakes. J.Eng.Div., Am.Soc. Civil Engr. 103:147-161.

Chapra, S.C., 1980. Application of the phosphorus loading concept to the Great Lakes. In Phosphorus Management Strategies for Lakes. Ed. R.C. Loehr, C.S. Martin, W. Rast. AnnArbor Science Publishers Inc., Ann Arbor, MI.

Connolly, J.P. and R.P. Winfield, 1984. A user's guide for WASTOX, a framework for modeling the fate of toxic chemicals in aquatic environments. Part 1. Exposure concentration. USEPA Gulf Breeze, FL. EPA-600/3-84-077.

Di Toro, D.M. and J.P. Connolly, 1980. Mathematical model of water quality in large lakes. Part 2: Lake Erie. EPA-600/3-80-065, USEPA ERL, Duluth, MN.

Di Toro, D.M., J.J. Fitzpatrick, and R.V. Thomann, 1981, rev. 1983. Water Quality Analysis Simulation Program (WASP) and Model Verification Program (MVP) - Documentation. Hydroscience, Inc., Westwood, NJ. for USEPA, Duluth, MN, Contract No. 68-01-3872.

Di Toro, D.M. and J.J. Fitzpatrick, 1993. Chesapeake Bay sediment nutrient flux model. Contract report EL-93-2. Prepared for U.S. Army Engineer Waterways Experiment Station, Vicksburg, MS.

Di Toro, D.M., D.J. O'Connor, J.L. Mancini and R.V.Thomann, 1973. A preliminary phytoplankton-zooplankton-nutrient model of Western Lake Erie. In Systems Analysis & Simulation in Ecology. Volume 3. Academic Press.

Di Toro, D.M., D.J. O'Connor and R.V.Thomann, 1971. A dynamic model of the phytoplankton population in the Sacramento-San Joaquin Delta. In Non Equilibrium Systems in Natural Water Chemistry, Advances in Chemistry Series, No. 106. American Chemical Society.

Di Toro, D.M., N.A. Thomas, C.E. Herdendorf, R.P. Winfield, and J.P. Connolly, 1987. A post audit of a Lake Erie eutrophication model. J. Great Lakes Res. 13:801-825.

Gedney, R.T., 1971. Numerical calculation of the wind-driven currents in Lake Erie. Ph.D. Thesis, Department of Fluid, Thermal and Aerospace Sciences, Case Western Reserve University.

Herdendorf, C.E., 1967. Lake Erie physical limnology cruise, midsummer 1967. State of Ohio, Department of Natural Resources, Division of Geological Survey, Columbus, Ohio.

HydroQual, Inc., 1995. A water quality model for Massachusetts and Cape Cod Bays: Calibration of the Bays Eutrophication Model (BEM). Prepared for the Massachusetts Water Resources Authority by HydroQual, Inc., Mahwah, NJ and Normandeau Associates, Bedford, NH.

Hydroscience, Inc., 1973. Limnological systems analysis of the Great Lakes: Part I. Prepared for the Great Lakes Basin Commission. Hydroscience, Inc., Westwood, NJ.

Kreiger, K.A., D.W. Schloesser, B.A. Manny, C.E. Trisler, S.E. Heady, J.J.H. Ciborowski, and K.M. Muth, 1996. Recovery of burrowing mayflies (Ephemeroptera: Ephemeridae: Hexagenia) in western Lake Erie. J. Great Lakes Res. 22:254-263.

Laws, E.A. and M.S. Chalup, 1990. A microalgal growth model. Limnol. Oceanogr. 35:597-608.

Markarewicz, J.C., 1993. Phytoplankton biomass and species composition in Lake Erie, 1970 to 1987. J. Great Lakes Res. 19:275-290.

Rast, W. and G.F. Lee, 1978. Summary analysis of the North American (U.S. portion) OECD Eutrophication Project: Nutrient loading-lake response relationships and trophic state indices. EPA-600/3-78-008. USEPA ERL, Corvallis, OR.

Schloesser, D.W., T.B. Reynoldson and B.A. Manny. Oligochaete fauna of western Lake Erie 1961 and 1982: signs of sediment quality recovery. J. Great Lakes Res. 21:294-306.

Thomann, R.V., D.M. Di Toro, R.P. Winfield and D.J. O'Connor, 1975. Mathematical modeling of phytoplankton in Lake Ontario. 1. Development and Verification. EPA-660/3-75-005, USEPA ERL, Corvallis, OR.

Thomann, R.V. and J.J. Fitzpatrick, 1982. Calibration and verification of a mathematical model of the eutrophication of the Potomac Estuary. Prepared for the Department of Environmental Services, Government of the District of Columbia, Washington, D.C.

Thomann, R.V., R.P. Winfield , D.M. Di Toro and D.J. O'Connor, 1976. Mathematical modeling of phytoplankton in Lake Ontario. 2. Simulations using LAKE 1 Model. EPA-660/3-76-065, USEPA ERL, Duluth, MN.

Vollenweider, R.A., 1966. Advances in defining critical loading levels for phosphorus in lake eutrophication. Mem. 1st Ital. Idrobiol. 33:53-83.

Vollenweider, R.A., W. Rast, and J. Kerekes, 1980. The phosphorus loading concept and Great Lakes eutrophication. In Phosphorus Management Strategies for Lakes. Ed. R.C. Loehr, C.S. Martin, W. Rast. Ann Arbor Science Publishers Inc., Ann Arbor, MI.

Wright, S., L.H. Tiffany, W.M. Tidd, 1955. Limnological survey of Western Lake Erie. U.S. Department of the Interior, Fish and Wildlife Service, Special Scientific Report: Fisheries No. 139, Washington D.C.