ABSTRACT

Predictive modelling of the fate of two persistent toxic organic chemicals, PCBs and Mirex, is discussed in light of the results from oceanographic scale investigations from the Niagara River to Lake Ontario, to the St. Lawrence River Estuary aquatic ecosystem. A mathematical model, TOXFATE, is used to run simulations of the fate of Mirex in Lake Ontario, a relatively small part of the total system, using a mass balance approach. TOXFATE features simulations of the fate of Mirex in lake water, plankton, benthos, suspended and bottom sediments, small and large fish (sculpins and salmonids). A friendly user interface (TOXSHELL) facilitates running the program on microcomputers. Concentrations of "dissolved" (the fraction not removed by high speed centrifugation) persistent toxic organic chemicals in the Niagara River are in the ng L^-1 or ng m^-3 range, yet the total load transported into Lake Ontario can be considerable given the high discharge of some 6000 m³ sec^-1. The river draining Lake Ontario is the St. Lawrence, and PCB loads actually double due to the various sources along the river. The insecticide and flame retardant, Mirex was essentially introduced from only two point sources, the Niagara and Oswego Rivers. The chemical is still detectable some 1000 kilometres downstream of the main site near Niagara Falls of its original introduction to this river-lake-estuary system. The system's recovery from Mirex pollution is related to the major natural aquatic processes of the system and is not compounded by continuing point and non-point source inputs. Simulations show the fast response of Mirex concentration in water following a reduction in loadings in the early 1960s and a much slower reaction of bottom sediments and fish to the same loadings reduction.

The National Water Research Institute has been involved for several years in the study of the fate of persistent toxic organic chemicals (PTOCs) in the waterways of the Niagara River to the St. Lawrence River Estuary system (Fig. 1). This study has taken the form of data collection as well as modelling. The system under study is made up of four types of units. The first is the high discharge (6,000 to 13,000 m³ sec^-1), rapid flow rivers, namely the Niagara and St. Lawrence. Second are the shallow, short residence time (days to weeks) riverine lakes on the St. Lawrence River, namely St. Francois, St. Louis and St. Pierre. Third is the deep, long residence time (10 years), Lake Ontario. Lastly, there is the freshwater-saltwater St. Lawrence inner and outer estuary. The rivers are essentially PTOC transport systems on a grand scale. The riverine lakes provide only temporary storage or sinks even for contaminants associated with sediments because these are eventually resuspended and moved on downstream. The major sinks where long-term effects become evident are Lake Ontario and the St. Lawrence Estuary. These sites are also where sediment associated contaminants can be permanently removed by deep burial in bottom sediments. Even so, some of the PTOCs transported still pass out of the system to the Gulf of St. Lawrence and beyond.

The fate of PTOCs in aquatic ecosystems involves many processes including dilution, sorption/desorption, volatilization, burial, and degradation by chemical and biological processes. PTOCs fate is also governed by the physicochemical characteristics of the chemicals, i.e., solubility, volatility, hydrophobicity, degree of partitioning (into or onto sediment and other particulates), lipophilicity, and resistance to degradation (by hydrolytic, photochemical, and biochemical routes) and by the characteristics of the receiving waterbodies, namely suspended particulate concentrations and type, sedimentation rates, trophic level, and food webs (Allan, 1989). The available information on PTOC contamination of the Niagara River to St. Lawrence Estuary system is extensive (Allan et al., 1990; 1991). In this paper, two chemicals, PCBs and Mirex, are used as examples of the aquatic transport and fate of PTOCs as this relates to predictive models.

TOXFATE, an organic TOXic contaminants FATE model, was developed (Halfon, 1984) in 1983 (Version 1.0) at the National Water Research Institute for the specific purpose of modelling the fate of toxic contaminants in large lakes. Halfon and Oliver (1990) continued the development of the model (version 2.0) while the present version (3.6) is the most recent. The present configuration includes benthic and water food chains. The model can handle a variety of chemicals from the very volatile (Halfon et al., 1990) to the very persistent.

The TOXFATE model is made user-friendly by the addition of a graphic user interface (TOXSHELL). The interface program lets users interact easily with the model, modify existing data files or add new ones, display concentrations in seven compartments and compare mass balances of different chemicals graphically. TOXSHELL has been written to aid users to run TOXFATE on MS-DOS machines. The format of this interface is drop-down menus and windows. Also, the display of data as XY graphs, for the presentation of computer simulations, and of bar graphs, for display of mass balances, is automatic. As an extra feature for advanced users, output files are created in such a way that they can be entered into any spreadsheet program for easy analysis of numerical results or creation of publication-grade graphs.

The extensive data available on Mirex distribution in various media of the Niagara River, Lake Ontario, St. Lawrence River estuary system has allowed design and testing of predictive models (Halfon, 1984; Halfon and Oliver, 1990; Halfon et al., 1990). Note that there are no longer any significant point sources of Mirex to the system. The fate model TOXFATE (Fig. 2) has evaluated the response of lakes to contaminant loadings using the mass balance approach. The main purpose is to predict concentrations in several biotic and abiotic compartments, given a known, or estimated, loadings. The model uses loading data on Mirex and concentration data in lake and river media, as well as laboratory derived physical-chemical characteristics to arrive at equations to predict such processes as sorption, volatilization and bioaccumulation. characteristics that need to be known for an aquatic ecosystem, include flow, residence times, suspended sediment concentrations, and sedimentation rates.

In summary, two kinds of data are necessary to run a fate model:

a) Physical-chemical characteristics of the contaminants, and

b) Hydrometeorlogical data and physical characteristics of the lake.

The physical-chemical characteristics are:

K_oc = organic carbon partition coefficient [(mg kg^-1)/(mg L^-1)]

K_ow = octanol water partition coefficient

H = Henry's Law constant [m³-atm mol^-1]

M_r = Molecular Weight

k_p = photolysis rate [m h^-1]

The model was calibrated with data from 20 chemicals present in the lake. The range of the log K_ow in the data file is from 4 to 7. The relevant data for Mirex are as follows: log K_ow, 7.13; log K_oc, 7.13; Molecular weight, 546; Henry Law's Constant, 5.41 x 10^-2 [m2 atm mol^-1]; and photolysis 6.28 x 10^-3 m day^-1. An important requirement is the mass balance of suspended sediments. Suspended sediments enter the lake from the surrounding areas, they are deposited to the bottom and leave the lake from the outflowing river. The mass balance of suspended sediments influence the value of the following parameters:

k₃ = sedimentation velocity of suspended sediments [1 m d^-1]

v_r = solids resuspension velocity 0.741 [mm year^-1]

v_p = solids permanent sedimentation velocity 0.889 [mm year^-1]

S = concentration of suspended sediments [1.2 x 10^-3 kg m^-3]

The values in the square brackets are computed for Lake Ontario. A separate mass balance calculation is required for these four parameters. Simons and Lam (1980) concluded that the uncertainty surrounding the formulation of sedimentation and resuspension may undermine the predictive ability of water quality models. The hydrometeorlogical data used to describe a lake include river inflows, for Lake Ontario the Niagara River, wind speed and water temperature. They are included in Table 1.

TOXFATE is formalized as a system of ordinary differential equations and the equations can be parameterized to represent a variety of contaminants. The equations are solved numerically. The model state variables are organic contaminant concentrations in suspended sediments, water, plankton, fish (alewife and salmonid), bottom sediments and benthos compartments.

The transport processes of advection in water are handled by a box model. The pollutant is assumed to be immediately and completely mixed within each spatial cell. In each spatial cell several mass balance equations are solved. These equations describe the movement of the contaminants across different phases. For example Eq. 1 has the following meaning: The change in time of the amount of toxic contaminants in the water column is a function of several processes including loadings to the lake, losses to the St. Lawrence River and to the atmosphere, and internal redistribution of the toxic contaminant in different phases, thus

dx

V -- = loading - river outflow - volatilization - sorption

dt to suspended sediments + desorption from suspended

sediments - sorption to plankton + desorption from

plankton - uptake from small fish - uptake from large

fish + resuspension from bottom sediments

In mathematical terms the fate model has the following form:

V dx₁/dt = L - Qx₁ - V k_v x₁ - V S k₁ K_oc foc_ss x₁ + V S k₂ x₂

- V P k₁ K_oc foc_p x₁ + V P k₂ f(T) x₃ - V F_s c₁ x₁ (1)

- V F_L c₁ x₁ + v_r A x₇ - V k_p x₁

where:

x₁ = toxicant concentration in water [mg m^-3]

x₂ = concentration of toxicant in suspended sediments [mg kg^-1]

x₃ = concentration of contaminants in plankton [mg kg^-1]

x₆ = concentration of contaminants in benthos [mg kg^-1]

k₁ = sorption from water [1.45 x 10^-6 m³ h^-1 L^-1)]

k₂ = desorption from suspended sediments [1.14 x 10^-4 h^-1]

k₃ = sedimentation velocity of suspended sediments [1 m d^-1]

k₄ = sedimentation velocity of plankton [1 m d^-1]

c₁ = ventilation rate of fish [m³ kg^-1 h^-1)]

c₂ = feeding rate of prey by predator [kg kg^-1 h^-1]

c₃ = excretion coefficient of fish [h^-1]

c₄ = feeding rate of benthos [m³ kg^-1 h^-1]

k_v = volatilization rate [h^-1]

f(T) = 1.024^T, where T is the water temperature [^oC]

Computed constants related to fish feeding:

K_oc = organic partition coefficient [(mg kg^-1)/(mg L^-1)]

K_ow = octanol water partition coefficient

H = Henry's Law constant [m³-atm mol^-1]

M_r = Molecular Weight

k_p = photolysis rate [m h-1]

The following two parameters are needed only if the Henry's Law constant is not known:

σ = Water solubility of chemical [mol m^-3]

v = Vapour pressure [mm Hg]

foc_b = fraction organic carbon of benthos [0.2 dimensionless]

foc_p = fraction organic carbon of plankton [0.18, dimensionless]

v_r = solids resuspension velocity 0.741 [mm/year]

v_p = solids permanent sedimentation velocity 0.889 [mm/year]

w_s= weight of each small fish [ 5 grams]

w_l= weight of each large fish [3100 grams]

A = area of Lake Ontario [1.95 10¹⁰ m²]

B = concentration of benthos [.075 kg dry weight m^-3]

D1 = active sediment layer thickness [0.018 m]

D2 = thickness of bottom sediment layers under D1 [0.005 m]

F_s = concentration of small fish [3.01 x 10^-5 kg m^-3]

F_L = concentration of large fish [1.78 x 10^-5 kg m^-3]

L = toxicant load [mg h^-1]

Q = river flow [2.556 x 10⁷ m³ h^-1]

P = concentration of plankton [10^-4 kg m^-3]

Φ = porosity of bottom sediment [.899 dimensionless]

ρ = bulk density [2400 kg m^-3]

S = concentration of suspended sediments [1.2 x 10^-3 kg m^-3]

T = water temperature [^oC]

V = volume of Lake Ontario [1.68 10¹² m³]

V_s = volume of the top bottom sediment layer [D1 A m3]

Model calibration

Parameters k1 and k₂ describe the affinity of toxic contaminants for suspended sediments and plankton. These parameters were calibrated during model development and should not be modified: Their values are as follows:

k₁ = sorption from water [1.45 x 10^-6 m³ h^-1 L^-1)]

k₂ = desorption from suspended sediments [1.14 x 10^-4 h^-1]

The rate constants associated with these parameters are regulated by other measured quantities, such as K_oc and f_oc. Thus k₁ and k₂ are only scaling constants and should not be modified.

Volatilization

Volatilization is an important removal process of contaminants from the lake surface to the air. Most models, including TOXFATE, represent this process using the standard two-layer representation of the water surface (Liss, 1973); the volatilization rate is modelled using the well known two-film representation of the water surface (Liss, 1973). Recently, Achman et al. (1993) and Hornbuckle et al. (1993) have measured the volatilization of polychlorinated biphenyls (PCBs) in Green Bay, Lake Michigan; they found that measured volatilization rates usually exceed predicted values. In TOXFATE, the volatilization parameter, k_v [h^-1] is computed as

k_v = 1 /(K_l + K_g) - A/V (2)

where, K_l [h m^-1] is the inverse of the liquid phase mass transfer coefficient, K_g [h m^-1] is the inverse of the gas-phase mass transfer coefficient, and A [m²] and V [m³] are the areas and volumes of the lake (lake compartments if the lake is divided into spatial compartments). The parameter K_l is computed as:

K_l = 1/ [k_O2 (32/M_r)^.5] (3)

where k_O2 is the oxygen (molecular weight 32) exchange constant [m h^-1] computed according to measured wind speeds (u) at 10 m over the water surface (Banks, 1975).

k_O2 = 1.51 x 10^-2 u^.5 for u < 5.5 m/s

(4)

k_O2 = 1.15 x 10^-3 u² for u > 5.5 m/s.

In this model the effects of intermittent turbulent and advective transport events are not included since, as Burns et al. (1981) noted, Whitman models usually differ very little from the predictions of more complex (e.g., surface renewal) models (Danckwerts, 1970). In Lake Ontario the average wind velocity is 6.8 metres per second at a height of 1 m and therefore k_O2 has a value of 0.053 [m h^-1]. The k_O2 factor is adjusted for water temperature with the following formula,

k_O2 = k_O2 1.024 ^(T-20) (5)

The gas resistance, K_g, is computed as follows:

K_g = 1 / {[W H (18/M_r)^0.5] / R (T+273.15)} (6)

where W [cm h^-1] = 0.1857 + 11.36 u [m/s] at 10 m over the water surface (Liss, 1973), H [m³-atm/mol] is Henry's Law constant, 18 is the molecular weight (M_r) of water and R is the gas constant. If H is not known, it can be computed using information on the vapour pressure v (in mm Hg) and the contaminant solubility σ (mol m^-3), then

H = v / (760 σ) (7)

The factor 760 in Eq. 7 converts the vapour pressure from mm Hg to atmospheres.

The sorption rate on suspended sediments and plankton depends on the K_oc of the chemical and on the organic content of the compartment. Sorption is modelled assuming nonequilibrium between water and suspended sediments.

V S dx₂/dt = V S k₁ K_oc foc_ss x₁ - V S k₂ x₂ - S k₃ A x₂

(8)

+ V F_s c₃ x₄ + V F_L c₃ x₅

where x₂ is the concentration of contaminants in suspended sediments [mg/kg], k₁, k₂ and k₃ are parameters whose respective value of 4.56x10^-7 [m³ L^-1 h^-1], 115/K_ow [h^-1] and 1 [m/day] were obtained by calibration (k₃ was provided by Endicott et al. 1990) when the model was developed for Lake Ontario (Halfon and Oliver, 1990), K_oc is the affinity of the contaminant for organic matter [(mg/kg) / (mg/L)], foc_ss is the fraction of organic carbon in suspended sediments, about 0.12 in Lake Ontario, x₁ is the concentration of the contaminant in water [mg/m³], V is the volume of water [m³] and S is the concentration of suspended sediments [kg m^-3].

Equation 8 is also used to model the uptake of toxic contaminants by plankton (phytoplankton and zooplankton are combined in one compartment) from water, the main difference between Eqs. 8 and 9 is that excretion is temperature dependent or:

V P dx₃/dt = V P k₁ K_oc foc_p x₁ - V P k₂ f(T) x₃ - P k₄ A x₃

(9)

- V F_s p_s c₂ x₃ - V F_L p_l c₂ x₃

where x₃ is the concentration of contaminants in plankton [mg/kg], foc_p is the fraction of organic carbon in plankton, about 0.18, k₄ is the sedimentation velocity, 1 [m/day], p_s is the proportion of plankton eaten by small fish in relation to the total food eaten (plankton + benthos), p_l is the proportion of plankton eaten by large fish in relation to the total food eaten (plankton + small fish) and all the other parameters are the same as in Eq.s 1 and 8. P is the concentration of plankton [kg m^-3].

Two classes of fish are included in the model, small fish of 32 g (w_s) in wet weight and large fish with a wet weight w_l of 2.2 Kg. These weight classes were chosen given the fact that some data were available from Lake Ontario (Halfon and Oliver, 1990). The individual fish weights are model parameters (w_s and w_l) and they can be changed. The contaminant absorption from water (Eqs. 16-22) is conceptualized according to a model by Neely (1979). Uptake from food (Eq. 23) follows a formulation introduced by Elliott (1975). In TOXFATE, the feeding preference of small fish depends upon the relative concentrations of plankton and benthos and the feeding preference of large fish depends upon the relative concentrations of plankton and small fish. This preference is quantified by the parameters b_s, p_s, p_l and s_l. These parameters are:

b_s = 2150. [kg/m³ / kg/m³] (10)

p_s = 0.5 [kg/kg] (11)

p_l = 0.1 [kg/kg] (12)

s_l = 0.9 [kg/kg] (13)

Model users might want to change the values of these four parameters rather than the defaults if feeding preferences of fish are known from experiments.

Equation 25 describes the excretion of contaminants by fish; here the assumption is that excretion is inversely proportional to K_ow and to the weight of the fish. This novel formulation is necessary to explain the observed increasing bioconcentration factors as the degree of chlorine substitution in the aromatic ring becomes greater; for example in fish the bioconcentration factor for 1,2,4-trichlorobenzene (1,2,4-TCB) is 5,800, for 1,2,3,4-tetrachlorobenzene (1,2,3,4-TeCB) 28,500, for pentachlorobenzene (QCB) 260,000 and for HCB, 2,400,000. Clearly these bioconcentration factors are not a linear function of log K_ow. The fish bioconcentration factors in Lake Ontario are much higher than those observed in laboratory microcosms where experiments are performed to a maximum exposure time of three to four months.

Fish absorb organic contaminants both from water through the gills and from their food, i.e., plankton, benthos, or smaller fish. To follow the movement of organic contaminants into two fish size classes, small fish like alewife and large fish like trout, the following nonequilibrium mass balance model is used (Leidy and Ploskey, 1980):

small fish,

V F_s dx₄/dt = V F_s c₁ x₁ + V F_s p_s c₂ x₃ + V_s F_s b_s c₂ x₆

(14)

- V F_s c₃ x₄ - V F_L s_l c₂ x₄

large fish,

V F_L dx₅/dt = V F_L c₁ x₁ + V F_L p_l c₂ x₃ + V F_L s_l c₂ x₄

(15)

- V F_L c₃ x₅

where x₄ and x₅ are the concentrations of the contaminant in small fish (F_s) and large fish (F_L) [mg/kg (dry weight)], p_s is the proportion of plankton eaten by small fish in relation to the total food eaten (plankton + benthos), b_s is the proportion of benthos eaten by small fish in relation to the total food eaten (plankton + benthos), s_l is the proportion of small fish eaten by large fish in relation to the total food eaten (plankton + small fish), p_l is the proportion of plankton eaten by large fish in relation to the total food eaten (plankton + small fish); the first term is the accumulation rate of contaminant from water [m³/(kg h)], c₂ (computed in Eq. 23 below) is the assimilation rate from food [kg/(kg h)], where predator is the fish itself and prey might be plankton, benthos, and for large fish, small fish. c₃ [h^-1] is the excretion rate computed with Eq. 25 below. F_s and F_L are the concentrations of small and large fish, respectively [kg m^-3]. This model formulation assumes that the fish has not necessarily achieved equilibrium with the current toxicant concentration in the natural water from which it was taken. For each fish class a different specific formulation is used for c₂.

More complex models have been presented in the literature (Jensen et al., 1982) to account for blood circulation and oxygen in the water. Given the lack of data or blood water partition coefficients and fat blood partition coefficient for most contaminants and fish, an unnecessarily complex model including such factors is deemed unreasonable until a definite need is shown.

Intake of organic contaminants from water is parameterized as a function of K_ow (Neely, 1979):

c₁ = 1000 A RV / W [m³ kg^-1 h^-1] (16)

where W, the individual fish weight [either w_s or w_l], is 5 g for the small fish compartment and 3180 g for the large fish compartment and 1000 is a conversion factor [g/kg]

RV = Q / (E_ox C_ox) [ m³ h^-1 ] (17)

and

Q = (32 1000 α W^.8) / [82 (T+273.15)] [mg O₂ h^-1] (18)

where, 32 and 82 are constants used to convert millilitres of O₂ to grams of O₂ and 1000 is to convert grams of O₂ to milligrams of O₂ (see Neely, 1979); for small fish,

α = 0.005 T - 0.018 [ml O₂ h^-1 g^-.8] (computed from data in Table 1 in Neely's paper) (19)

and, for large fish,

α = 0.017 T + 0.058 [ml O₂ h^-1 g^-.8] (computed from data in Table 1 in Neely's paper) (20)

E_ox is the ventilation efficiency (Neely, 1979) [0.75 dimensionless),

C_ox = 1000 [14.45 - 0.413 T + 0.00556 T²] [mg O₂ m^-3], (21)

where there are 1000 [L m^-3], and

A = (0.07 log K_ow - 0.02), [dimensionless] (22)

The weight of the large fish was derived from available trout data. The function A (the efficiency transfer coefficient) was derived by Neely to account for the different assimilation efficiency of different contaminants. Note that the absorption efficiency is a function of log K_ow, not K_ow itself, therefore small adsorption efficiency differences between two contaminants may not result in large differences in observed bioconcentration factors (Halfon and Oliver, 1990). This issue has been accounted for in TOXFATE with a new formulation of the excretion process, Eq. 25 below.

The second major source of contaminants is food. In TOXFATE small fish feed on plankton and benthos while the big fish feed on small fish and plankton. The feeding rate c₂ is defined as:

c₂ = D / W [kg /(kg h)] (23)

where W can be either w_s or w_l, D is computed using the model by Elliott (1975) assuming 70% growth (Elliott provides model formulation for growth rates from 10-90%); the general formulation of the Elliott model is:

b₁

D = a W exp (b₂T) [kg h^-1] (24)

where a [kg g^-1], b₁ [dimensionless] and b₂ [^oC^-1] have the relative values of [a: 3.58x10^-5; 1.33x10^-4; 3.16x10^-3], [b₁: 0.769; 0.781, 0.764], [b₂: 0.335; 0.129; -0.76] for the three temperature ranges (0-6.7 ^oC, 6.7-14.9 ^oC, and higher than 14.9 ^oC).

To solve the problem of computing the very different bioconcentrations observed for contaminants with K_ows ranging over several orders of magnitude the excretion rate process in TOXFATE was formalized as:

c₃ = 8.3 10³ / K_ow^.72 [h^-1] for small fish and (25)

c₃ = 0.15 10³ / K_ow^.56 [h^-1] for large fish, (26)

where 2.28x10³ (g^0.75 h^-1) is an empirical scaling constant. The scaling constant is very large since K_ow can range from 10,000 to 1,000,000. The justification for this formulation (Eq. 25) is that contaminants with a high K_ow are retained in the fish much more than less hydrophobic compounds and bigger fish retain more contaminants than smaller fish since their metabolism is relatively slower. Before developing Eq. 25, many other formulations published in the literature were tested, but none were successful since they are all based on log K_ow rather than on K_ow itself. Thomann (1989) has also proposed a slightly more complex formulation for Eq. 25. His rationale is that this formulation implies the same mechanisms that hinder or enhance transport into the organism are operative in the transport out of lipid pools across lipoprotein membranes and into the excretory systems.

Benthos accumulates toxic contaminants from the bottom sediments. Losses occur by excretion and predation to the bottom feeding fish. Equation 26 describes the flow of toxic contaminants:

V_s B dx₆/dt = V_s B c₄ x₇ - Vsed F_s b_s c₂ x₆ (27)

where x₆ is the concentration of contaminants in benthos [mg/kg], B is the concentration of benthos [kg/m³], c₄ [h^-1] is computed according to a modified feeding a model by Leidy and Ploskey (1980). The feeding rate c₄ is defined as:

c₄ = G1 AS [m³ kg^-1 h^-1] (28)

where B is the concentration of benthos [kg per cubic metre of bottom sediments] and AS is an assimilation factor = K_ow / ( K_ow + k₅) [dimensionless, k₅ = 2 10⁶].

G1 = (G2 foc_b D1 1000)/24 [m3 kg-1 h-1] (29)

where foc_b is the fraction organic carbon of benthos [0.5 dimensionless]and (using a notation similar to that of Leidy and Ploskey, 1981),

G2 = Z [1 - e^{(-Y BB)}] [m³ /(mg C day)] (30)

(where Y = 10 ^{[-2.9664 - 0.9787 LOG10 (GMAX)]} (31)

and

GMAX = Z BB, (32)

where

Z = 10 ^{[-3.2295 -.0678 LOG10 (BB)]} (33)

and

BB = ρ (1-Φ) foc_s D1 1000 [mg C / m²] (34)

and foc_s is the fraction organic carbon of bottom sediments [0.03, dimensionless], Φ is the bottom sediments porosity and ρ is the bulk density.

Bottom sediments (top layer)

Bottom sediments are divided into 34 layers. The top layer has a depth of 1.8. The other 33 layers have a depth of 5 mm each. The main sources of toxic contaminants to bottom sediments are the suspended sediments and plankton. Another source is the release from benthos. This source is relatively minor. Bottom sediments lose toxic contaminants to the overlaying water, by resuspension, to benthos and to the deep sediment layers. Equation 35 describes these processes:

V_s dx₇/dt = S k₃ A x₂ + P k₄ A x₃- v_r A x₇ - v_p A x₇

(35)

- V_s B c₄ x₇

where x₇ is the concentration of contaminant in the bottom sediments [mg/m³], V_s is the volume of the bottom sediments computed as the area of the lake, A, multiplied by the active sediment layer, D1 [2 cm] and v_p are 0.741 [mm/year] and 0.889 [mm/year], respectively (Endicott et al., 1990)

Bottom Sediments (deeper layers)

The top layer of the bottom sediments is continuously covered by new sediments sinking from the water column. From a mathematical point of view toxic contaminants seem to move downward at a rate v_p. Thus, the equations representing the concentration of toxic contaminants in the deep layers of the bottoms sediments are:

V_s dx_i/dt = v_p A x_i-1 - v_p A x_i (36)

where i ranges between 8 and 40.

In the calibration stage TOXFATE has quantified the behaviour of four chlorobenzenes, 1,2,4-TCB, 1,2,3,4-TeCB, QCB and HCB, in Lake Ontario. The purpose of this procedure is to verify the correct formulation of the model using four chemicals which belong to the same family of chlorobenzenes but which have different physico-chemical characteristics. The predicted and observed concentrations in the different components of the ecosystem agree within 50% for averages and are completely indistinguishable if standard deviations are considered (Halfon and Oliver,1990). This fact is encouraging since the four chlorobenzenes have very different chemical properties, the estimated past loadings vary widely between 15,000 kg over the whole 75 years for HCB and over 300,000 for 1,2,4-TCB during the same period. Furthermore, the main encouraging factor is that the model was able to predict concentrations in large fish, here trout, quite accurately even if the low chlorinated chlorobenzenes, 1,2,4-TCB, 1,2,3,4-TeCB and QCB have fish concentrations about 10 times lower than those of HCB (Halfon and Oliver, 1990) while the log K_ow's range only over a factor of two.

TOXSHELL has been written to aid users to run TOXFATE on MS-DOS machines. The format of this interface is drop down menus and windows (Fig. 3). The program controls storage and retrieval of data in three subdirectories, LAKES, CHEM and RESULTS. The user, however, is not aware of the automatic internal mechanism of data handling. Also automatic is the display of data as XY graphs for the presentation of computer simulations, and of bar graphs for display of mass balances. As an extra feature for advanced users, the files in the RESULTS subdirectory are created in such a way that they can be entered into spreadsheet programs for easy analysis of numerical results or creation of graphs. When TOXSHELL is used to run TOXFATE, the user must choose a chemical filename and a lake file name. TOXFATE can be run from TOXSHELL only if these two input files have been loaded. The input files can also be created or edited in a text editor, for example BRIEF, but the editor within TOXSHELL has been designed for the specific purpose to enter data for TOXFATE. TOXFATE's output is written to four different files PARAM.DAT, WAT.DAT, SED.DAT and CHART.DAT. TOXSHELL then renames the files and saves them in the appropriate directories with new extensions .PRM, .WAT, .SED and .CHT.

Viewing the results using a spreadsheet program.

The shell program TOXSHELL copies the results to four files, three of which can be imported into a spreadsheet (e.g. Quattro Pro) for analysis. It is important that the user specifies different names for each chemical output. The output files are stored in the subdirectory RESULTS.

As output TOXSHELL creates four files: --------

*.PRM (containing the relevant parameters used) \ Can be

*.WAT (water concentrations) = Imported to

*.SED (sediment concentrations) / Quattro Pro

or Lotus

-------

*OUT.CHT (file used for checking the mass balance;

used mostly for debugging the software)

These four files contain information about the results of the simulations. These files (except *.CHT) can be imported individually in a spreadsheet program to look at the numbers or to plot the results. The capability of assigning input file names allows the user to run TOXFATE with different chemicals. The capability of assigning output file names allows the user to modify the input file of a chemical, for example change the K_oc, and save the results in different output files. Different scenarios can then be compared since all output files, one for each scenario, are stored on disk.

Concentrations of "dissolved" (the fraction not removed by high speed centrifugation) PTOCs in the Niagara River are in the ng L^-1 or ng m^-3 range, yet the total load transported into Lake Ontario can be considerable given the high discharge of some 6000 m³ sec^-1. Although PTOCs with large partitioning coefficients such as the higher chlorinated PCBs, are concentrated in the suspended solids phase (ng g^-1 range), much of the chemical load is still transported in the operationally defined dissolved phase because of the low suspended solids concentrations in the Niagara River (1 to 10 mg L^-1). At its source in Lake Ontario, the flow of the St. Lawrence River is some 7700 m³ sec^-1 and rises to 12,700 m³ sec^-1 at Quebec City. A given water mass traverses the river in five to seven days. Concentrations of suspended solids rise from some 1 mg L^-1 at the river source to some 10 mg L^-1 Quebec City, at the upper end of the St. Lawrence Estuary. Water samples are processed by high speed centrifuge followed by large volume extraction of the effluent to determine changes in the fraction of the PTOC load transported in the operationally defined particulate and dissolved phases. PCB concentrations in the suspended solids decreased downstream from around 1000 ng g^-1 to 200 ng g^-1, but because there is a ten times increase in suspended solids, the PCB load actually doubled due to the various sources along the river. These sources include major aluminum smelters and automobile plants and other industries and several large urban centres including Montreal. The dissolved phase PCBs were determined in the centrifugate from the Westfalia centrifuge. The particulate phase was detained in the extract from the Westfalia. However, the centrifugate from the Westfalia was then centrifuged using a high speed Sorval centrifuge. As can be seen in Fig. 1, just downstream of Montreal, a significant portion of the Westfalia so called 'dissolved' phase (647 kg yr^-1) can be removed as recaptured fine particles or colloids by the Sorval (153 kg yr^-1.) Concentrations of PCBs in the super-centrifuged colloids (particles approximately less than 0.1 μ size extracted by the Sorval centrifuge, Fig. 4) were one to two orders of magnitude higher than in the particulates normally extracted by only high speed centrifugation (particles approximately less than 1 μ size extracted by the Westfalia centrifuge). The PCBs introduced to the river are mainly in the dissolved phase as are those introduced along its course but they become more fractionated into or onto the particulate phase during passage downstream.

The insecticide and flame retardant, Mirex, was first found in Lake Ontario in the early 1970s. By the late 1980s, sampling the St. Lawrence had revealed that this PTOC had been transported in the water column as far as the St. Lawrence Estuary and beyond. The chemical was still detectable some 1000 kilometres downstream of the main site near Niagara Falls of its original introduction to this river-lake-estuary system. The river draining Lake Ontario is the St. Lawrence. At its point of outflow from Lake Ontario, the concentration of Mirex in suspended solids is around 5 ng g^-1 and decreases to some 1 ng g^-1 near Quebec City. The highest value for Mirex in the operationally defined dissolved phase was 13 pg L^-1. Further evidence of the long-range aquatic transport of Mirex comes from bottom sediment cores from Lake Ontario and the Laurentian Trough (lower St. Lawrence Estuary) (Fig. 1). The peak concentration in the Lake Ontario core was 65 ng g-1 in 1967, and in the Laurentian Trough 0.19 ng g^-1 in 1981 (Fig. 5).

The Mirex in the system was essentially introduced from only two point sources, the Niagara and Oswego Rivers at the southwest and southeast ends of Lake Ontario, respectively (Fig. 1). Because of this, the system's recovery from Mirex pollution is related to the major natural aquatic processes of the system and is not compounded by continuing point and non-point source inputs. The results of the Mirex model would also apply to similar PTOCs with high partition coefficients, for example some of the more chlorinated PCB isomers. Some simulation results from the model are shown in Fig. 6. These simulations show the fast response of Mirex concentration in water following a reduction in loadings in the early 1960s and a much slower reaction of bottom sediments and fish to the same loadings reduction.

The best way to resolve the transport and fate of PTOCs within aquatic ecosystems and to have confidence in models predicting future levels of exposure or rates of ecosystem recovery is by correct collection and accurate analysis of multi-media samples. Interpretation of this and related limnological data provides hypotheses for further testing either by subsequent sampling; laboratory microcosm simulations; physical-chemical testing; or by theoretical modelling based on the properties of the contaminants and specific polluted waterbodies. Only on the basis of such knowledge can meaningful and realistic control actions, regulations, and remedial action plans be decided upon and implemented.

TOXSHELL has been written in C for the specific purpose of letting users use the TOXFATE model in a friendly way. TOXFATE can also run as a FORTRAN program on micros and mainframe computers. The fact that the output can be viewed both in graphical mode, using TOXSHELL, or in both a numerical or graphical way, using any spreadsheet, allows much flexibility. Often models take much effort in data collection, thought and programming, only to be left on the shelf because of poor user interface. TOXSHELL is a first effort to eliminate this drawback and let interested users employ TOXFATE to analyze their problems.

A simulation model, TOXFATE, has been employed to quantify the behaviour of PCBs and Mirex in Lake Ontario. Except for PCBs, the concentrations of toxic contaminants in Lake Ontario are much below Ontario Water Quality Objectives but given the large volume of the lake, a considerable amount of toxic contaminants are stored in this ecosystem. The main problem of the prediction of toxic contaminants fate, once a model has been satisfactorily validated for a few compounds, is the lack of loading data; this problem is particularly evident in a large lake system. Presently, the Niagara River is regularly monitored for compounds in the water and in suspended sediments (Durham and Oliver, 1983), however, such monitoring has taken place only since the late 1970's and earlier loadings data can only be inferred from bottom sediments data. For new or recently released compounds, or for compounds not yet identified in the lake, these data are not available and may not be available for a foreseeable future; for contaminants which are already entering the lake, the bottom sediments are a good source of past loadings data provided that a comprehensive data base with concentrations and dating is available, as shown by Oliver et al. (1989).

Achman, D.R., K.C. Hornbuckle, and S.J. Eisenreich. Volatilization of Polychlorinated Biphenyls from Green Bay, Lake Michigan. Env. Sci. Technol., 27: 75-87; 1993.

Allan, R.J. Factors affecting sources and fate of persistent toxic organic chemicals: examples from the Laurentian Great Lakes. Chapter 9.3 in Aquatic Toxicology, Boudou A. and Ribeyre F. (Eds), Pub CRC Press, pp. 219-248; 1989.

Allan, R.J., Campbell, P.G.C., Forstner, U. and Lum, K.R. (Eds.). The fate and effects of toxic chemicals in large rivers and their estuaries. Sci. Total Environment 97/98, 871 pp.; 1990.

Allan, R.J., Ball, A.J., Cairns, V., Fox, G.A., Gilman, A.P., Peakall, D.B., Piekarz, D., Van Oostdam, J.C., Villeneuve, D.C. and Williams, D.C. Toxic chemicals in the Great Lakes and associated effects. Volumes 1 and 2. Pub. Govt. of Canada, Dept. Supply and Services 1991, Ottawa, Cat. No. En 37-95/1990-1E, ISBN 9-662-18317-7, 748 pp; 1991.

Banks, R.B. Some features of wind action on shallow lakes. J. Environ. Eng. Div., Proc. ASCE 101EE5): 813-827; 1975.

Burns, L., Cline, D.M. and Lassiter, R.R., 1981. Exposure Analysis Modeling System EXAMS): User manual and system documentation. U.S.EPA, Athens, Georgia, 440 pp.

Comba, M.E., Palabrica, V.S., Wasslen and Kaiser, K.L.E. St. Lawrence River trace organic contaminants study (Part III), 1985. National Water Research Institute, Burlington, Ontario, Contribution No. 90-1, 39 pp.; 1990.

Comba, M.E., Norstrom, R.J., MacDonald, C.R. and Kaiser, K.L.E. A Lake Ontario-Gulf of St. Lawrence dynamic mass balance for Mirex. Env. Sci. Technol., 27: 2198-2206; 1993.

Danckwerts, P.V. Gas-liquid reactions. McGraw-Hill Report Co., New York, 276 pp; 1970.

Durham, R.W. and B.G. Oliver. History of Lake Ontario contamination from the Niagara River by sediment radiodating and chlorinated hydrocarbon analysis. J. Great Lakes Res., 9:160-168; 1983.

Elliott, J.M. The growth rate of brown trout Salmon trutta L.) fed on reduced ratios. J. Animal Ecol., 44:823-842; 1975.

Endicott, D.D, Richardson, W.L., Parketon, T.F. and DiToro, D.M. A steady state mass balance and bioaccumulation model for toxic chemicals in Lake Ontario. U.S. EPA.; 1990

Halfon, E. Error analysis and simulation of Mirex behaviour in Lake Ontario, Ecol. Model., 22: 213-252; 1984.

Halfon, E.H. and Oliver, B.G. Simulation and data analysis of four chlorobenzenes in a large system, Lake Ontario, with TOXFATE, a contaminant fate model, S.E. Jorgensen (Ed.) Modelling in Ecotoxicology, Elsevier, pp. 197-214; 1990.

Halfon, E., T.J. Simons and W.M. Schertzer. Modelling the spatial distribution of seven halocarbons in Lake St. Clair in June 1984 using the TOXFATE model. J. Great Lakes. Res., 16: 90-112; 1990.

Hornbuckle, K.C., Achman, D.R., and S.J. Eisenreich. Over-water and over-land Polychlorinated Biphenyls in Green Bay, Lake Michigan. Env. Sci. Technol., 27: 87-98; 1993.

Jensen, A.L., Spigarelli, S.A. and Thommes, M.M. PCB uptake by five species of fish in Lake Michigan, Green Bay of Lake Michigan and Caiuga Lake, New York. Can. J. Fish. Aquat. Sci., 39:700-709; 1982.

Leidy, G.R. and Ploskey, G.R. Simulation modelling of zooplankton and benthos in reservoirs: documentation and development of model constructs. US Army Engineer Waterways Experiment Station Environmental Laboratory, Vicksburg, Miss., Technical Report E-80-4, 300 pp.; 1980.

Liss, P.S. Processes of gas exchange across an air-water interface. Deep-Sea Res., 20:221-238; 1973.

Neely, W.B. Estimating rate constants for the uptake and clearance of chemicals by fish. Environ. Sci. Technol., 13:1506-1510; 1979.

Oliver, B.G., Charlton, M.N. and Durham, R.W. Distribution, redistribution and geochronology of polychlorinated biphenyl congeners and other chlorinated carbons in Lake Ontario sediments. Environ. Sci. Technol., 23: 200-208; 1989.

Simons, T.J. and Lam, D.C.L. Some limitations of water quality models for large lakes: a Case study of Lake Ontario. Wat. Res. Res. 16: 105-116; 1980.

Thomann, R.V. Bioaccumulation model of organic chemical distribution in aquatic food chains. Environ. Sci. Technol., 23: 699-707; 1989.

Table 1: Data that describe Lake Ontario as input to TOXFATE.

Total volume 1.68 x 10¹² m³

Total area 1.95 x 10¹⁰ m²

Flow 7100 m³ s^-1

Water temperature 9 ^oC

Wind at one metre above water 4.5 m s^-1

Suspended sediment conc. 1.1 mg L^-1

Plankton concentration 0.1 mg L^-1 dry weight

Small fish concentration 0.05 mg L^-1 dry weight

Large fish concentration 0.018 mg L^-1 dry weight

Benthos concentration 75. g m^-3 dry weight

Bottom sediments density 2400 kg m^-3

Depth of top layer 0.018 m

Depth bot. layer 0.005 m

Permanent sediment velocity 0.889 mm y^-1

Resuspension velocity 0.741 mm y^-1

Sediment. vel. suspended sed. 1 m day^-1

Sedimentation vel. plankton 1 m day^-1

foc suspended sediment 0.12 dimensionless

foc plankton 0.18 dimensionless

foc benthos 0.20 dimensionless

Michaelis-Menten assim. param. 2 x 10⁶ dimensionless

foc sediment 0.03 dimensionless

Sediment porosity 0.888 dimensionless

Sorption from water 1.45 x 10^-6 m³ h^-1 L^-1

Desorption from sus. sed 1.14 x 10^-4 h^-1

Loading kg day^-1

LIST OF FIGURES

FIGURE 1: Schematic of the main components of the Niagara River to St. Lawrence River Estuary system and their primary role in PTOC transport and fate, and location of sediment cores in Lake Ontario and the Laurentian Trough.

FIGURE 2: Structure of the fate model TOXFATE. Each box represents the concentration of toxic contaminants in a state variable.

FIGURE 3: TOXSHELL, main menu with options.

FIGURE 4: "Dissolved", colloidal and particulate PCB flux [kg yr-1] in the St. Lawrence River. (adapted from Comba et al., 1990). See text for explanation.

FIGURE 5: Comparison of Mirex profiles in age dated bottom sediment cores taken from Lake Ontario and the Laurentian Trough at sites shown in Fig. 1 (adapted from Oliver et al., 1989 and Comba et al., 1993).

FIGURE 6: Predicted Mirex concentrations (water μg L^-1, other media μg g^-1) in Lake Ontario media with time.