A Kinetic Model of Copper Cycling in San Francisco Bay

A two-dimensional, depth-averaged kinetic model of copper cycling was developed for the San Francisco Bay estuary. Adsorption and desorption reaction rate constants were determined from experimental sorption experiments. To calibrate the model, processes related to aqueous speciation were included. The model was used to predict spatial and seasonal trends in the adsorption and desorption of copper. Model predictions show that copper is continually being re-partitioned between sediment and water. Re-partitioning is prevalent near tributary and anthropogenic sources. It also occurs between segments of the bay, in response to differences in salinity and the availability of organic ligands dissolved in the water. In areas of restricted circulation such as the South Bay, copper adsorbed onto settling particles during wet season storm events acts as a source to the water column during the dry season. The relative contribution of resuspended benthic sediment to dissolved copper concentrations is highly variable in the bay. In the North Bay, dissolved copper is principally introduced from the San Joaquin-Sacramento Delta. In the South and lower South bays, desorption from sediment during the dry season may contribute as much as 20% of the total mass input of dissolved copper. Improvement of water quality can be achieved by reducing loads; however, changes are predicted to take years.


INTRODUCTION
San Francisco Bay is a highly urbanized estuary that has undergone significant changes in the past 150 years. Of particular importance to the health of aquatic and benthic organisms living in the bay has been the introduction of trace metals accompanying industrialization. Anthropogenic sources of metals have been revealed in the past twenty years during a concerted effort by research institutions, environmentalists, regulators, local governments, and industries to monitor bayarea sediment, water, and biota.
One of the key findings of estuarine investigations conducted in San Francisco Bay, particularly those spearheaded by the San Francisco Estuary Institute (SFEI), and represented by its annual regional monitoring program reports, is that dissolved copper concentrations are elevated above background levels. In fact, dissolved copper concentrations have been found to often exceed water quality objectives in an area below the Dumbarton Bridge known as the lower South Bay. These exceedances have prompted the San Francisco Bay Regional Water Quality Control Board to begin developing and implementing pollution prevention strategies and source control measures to improve water and sediment quality.
Understanding the processes leading to copper enrichment and depletion in the bay is critical to developing successful remediation strategies. To date, extensive data have been collected on copper concentrations in water and sediment (Flegal and others 1991;SFEI 1994SFEI -1999Hornberger and others 1999), porewater (Rivera-Duarte and Flegal 1997), and watershed tributary sources (BASMAA 1996). In addition, copper speciation studies have demonstrated the importance of dissolved organic material on copper geochemistry (Donat and others 1994). Results from field studies have also been used to show that the two largest sources of total (dissolved plus adsorbed) copper are resuspended sediment and anthropogenic inputs, respectively (Rivera-Duarte and Flegal 1997). Finally, Flegal and others (2004) reviewed previous studies and identified interannual differences in flushing and episodic storm events as important variables affecting dissolved copper concentrations.
Powerful tools for synthesizing and interpreting data from these investigations include hydrodynamic, sediment-transport, and chemical partitioning numerical models, two of which have previously been employed to describe copper cycling in San Francisco Bay (Chen and others 1996;Wood and others 1995). Chen and others (1996) developed a model calibrated to hydrodynamic, advection, dispersion, and sediment transport conditions; however, they did not include time-dependent chemical partitioning of copper between adsorbed and dissolved forms. By contrast, Wood and others (1995) used a kinetic-based formulation to demonstrate that equilibrium and non-equilibrium partitioning of copper between dissolved and adsorbed phases may affect basin-scale concentrations. The model was not, however, calibrated with respect to either hydrodynamic or transport conditions.
The objective of the current study was to combine time-dependent chemical processes (Wood and others 1995) with representative hydrodynamic and sediment transport conditions (Chen and others 1996), to predict the relative importance of time-dependent sorption to copper cycling in San Francisco Bay. This task was accomplished by first calibrating the hydrodynamic, advection, dispersion, and sediment transport components of the model (see Appendix A). Next, site-specific adsorption and desorption reaction rate constants were experimentally derived and used to attune the kinetic-based chemical partitioning submodel. Predictions were analyzed for spatial and seasonal trends in copper concentrations and distribution coefficients. The relative importance of processes such as source loading and nonequilibrium sorption on dissolved copper concentrations were determined using mass balance calculations.
hydraulics-related phenomena in estuaries, coastal waters, and seas where stratification can be neglected. The software consists of two linked modules. The first is a hydrodynamic module (MIKE 21 HD) that solves the timedependent, vertically-integrated equations of continuity and conservation of momentum in two horizontal directions. The second is a heavy metals module (MIKE 21 ME) that uses output from the hydrodynamic module to settle and resuspend sediment and sorbed metals. It also advects and disperses sediment and dissolved metals in the water column.
MIKE 21 ME calculates the distribution of adsorbed and dissolved metals by integrating derivative expression for the change in dissolved and adsorbed concentrations with time using a fourth-order Runge-Kutta technique in an integrated two-step procedure with the advection-dispersion equation. The primary assumption inherent in this approach is that adsorption and desorption can be empirically expressed in terms of a secondorder reaction of the form: (1) where M i = the i th dissolved metal, S j = the j th surface site located on a suspended grain of sediment, and M i S j = the adsorbed metal complex (Luoma 1990).
Using this representation, the rate of adsorption to particulate matter is equal to the product of the concentrations on the left-hand side of equation 1 and a factor k + , an experimentally measured adsorption rate constant: (2) where brackets indicate a chemical concentration (or activity). The desorption rate is similarly given by the product of the concentrations on the right-hand side of equation 1 and a desorption rate constant k -: (3) Using equations 1 and 2, the change in the dissolved metal concentration Mi as a function of time is simply the difference between desorption and adsorption rates: (4) and the change in the adsorbed metal concentration is equal to the adsorption rate minus the desorption rate:

(5)
Because rate constants that describe changes between the dissolved and adsorbed phases are dependent on a unique chemical environment where temperature, pressure, and concentration of all other components in the system are invariant, the actual equations utilized by MIKE 21 ME to calculate changes in the dissolved and adsorbed phases are: where the additional terms following the rate constants are forcing functions that allow adsorption and desorption to vary with pH, redox potential ( ), temperature (T ), and salinity (SAL), respectively, and the final terms in equations 6 and 7 are diffusive exchange with the porewater and sediment exchange, respectively. The forcing functions for temperature (T ) and salinity (SAL) are: where a, b, c = forcing function coefficients determined during calibration using timevarying temperature data provided by the U.S. Geological Survey (USGS 2000) and salinity data from a calibrated salinity model for the selected simulation period. (The pH and redox forcing functions were not used in this study because there was insufficient spatial and temporal resolution in existing data to distinguish such effects.)

Experimental Rate Constants
Sediment and water samples necessary to determine the experimental rate constants for MIKE 21 ME were collected from the area designated as the South Bay station in Figure 1. Surface sediment grab samples were taken with a 0.1 m 2 , Kynar®-coated stainless Gray-O'Hara box core. Sediment was sealed into three-gallon plastic buckets fitted with lowdensity polyethylene (LDPE) liners and placed in walk-in coolers maintained at 4°C. Because it was necessary to keep the sediment in its ambient oxidation state, as much care as possible was taken to avoid exposure of the sediments to air during handling both in the field and laboratory. Ambient water was collected into pre-cleaned containers from the bay sites by use of a peristaltic pump equipped with pre-cleaned Teflon® tubing. Sediment was first centrifuged at 4,000 rpm to remove porewater. Ambient baywater was similarly centrifuged at 4,000 rpm to remove suspended material. Prior to the experiment, 40 L of centrifuged water were placed in a covered 60-L stainless steel cockpot coated with Halar®, a commercial Teflon® coating. (Initial investigations into different coatings available indicated that Halar® had the best properties to prevent chemical adsorption.) Water temperature was 16°C, which was maintained throughout the experiment by placing the vessel in an isolated room maintained at a temperature of 16 ± 2°C. The contents of the vessel were stirred continuously using a Halar® coated stainless steel propeller-type stirrer.
At the start of the experiment, water was fortified to achieve an initial dissolved copper concentration of 46.8 µg L -1 (micrograms of copper per liter of water, equivalent to parts per billion [ppb]). Twenty grams of de-watered sediment were subsequently introduced into the vessel and allowed to react with the fortified water for a period of 1 h (subsampling at 1 h showed that the dissolved copper had achieved its lowest concentration due to adsorption on the sediment). Additional subsamples were taken at 3, 6, 9, 12, and 24 h, plus 5 d and 10 d (a total of 9 times including 0 and 1 h). At each subsampling interval, water temperature was determined by measuring the temperature of a separate container of water to avoid contaminating the experiment. Enough water was removed from the vessel to analyze for total and dissolved copper and total suspended sediment concentration. For the total metals analyses, the samples were immediately acidified with trace-metal-grade Depth (meters) nitric acid. For the dissolved metals subsamples, water was immediately centrifuged at 4,000 rpm prior to filtration to ensure that the 0.45-µm filters did not clog and that the waters would be removed from the sediments as soon as possible. These subsamples were immediately filtered and the filtrates acidified with nitric acid. Analyses were performed using standard EPA procedures (EPA Method 7211 or 6020). A review of metals data from the first and last experimental sampling events (Table 1) shows that the total copper concentrations were similar to the spike amounts. This indicates that loss due to adsorption on the vessel walls was minimal.
The measured distribution coefficient as a function of time is shown on Figure 2. A decrease in the distribution coefficient in the figure indicates net desorption and an increase indicates adsorption. As shown on the figure, adsorption occurred during the first hour of the experiment, but was followed by net desorption for the remaining 239 hours. These results are consistent with a state of metastable equilibrium being quickly established between the ions added to the solution and the suspended sediment. Reactions such as aqueous metal complexing and replacement of metal surface species by competing cations likely resulted in a continuous increase in the amount of dissolved metal as a state of equilibrium was approached. To obtain adsorption and desorption rate constants from these experimental results, it was first assumed that metal binding sites were present in large excess over the dissolved metals (Jannasch and others 1988). Sorption could then be expressed as a first-order reaction of the form: where M i and M s = the dissolved and adsorbed phases, respectively. The change in the equilibrium distribution coefficient with time was subsequently written in terms of the adsorption k + and desorption krate constants: where K D = the distribution coefficient between the dissolved and adsorbed phase (L kg -1 ), t = time, SSC = suspended sediment concentration, M T = total metal concentration, M i,o = initial dissolved concentration. Regression on the desorption phase of the experiment produced the predicted timevarying distribution coefficient shown on Figure 2, which corresponds to adsorption and desorption rate constants of 0.010 L mg -1 d -1 and 0.47 d -1 , respectively.
The initial adsorption process shown on Figure 2 is consistent with the hypothesis that adsorption between copper ions and surface sites is faster than subsequent desorption and aqueous complexation. Prior to the sorption experiment, the unspiked sample had a measured dissolved concentration of 3.9 µg L -1 . After the spike and particulate matter had reacted in the sample for a period of 1 h, this concentration was 3.4 µg L -1 . The fact that the concentration did not increase implies that the copper spike was adsorbed onto sorption sites on the suspended particulate matter, without forming additional aqueous complexes. The occurrence of higher dissolved concentrations during subsequent sampling is consistent with desorption and complexation of copper by organic and/or inorganic ligands through slower chemical reaction mechanisms.
The desorption rate constant calculated in this study is ten times faster than the rate calculated by Glegg and others (1988) for sediment from the Tamar estuary (0.47 d -1 vs. 0.05 d -1 , respectively). This difference is consistent with the opposing processes believed to be occurring in the two respective studies. In the desorption experiments of Glegg and others (1988), an initial decrease was observed in the distribution coefficient at short time scales (i.e., desorption was occurring). This was followed by an increase in the distribution coefficient (i.e., adsorption was occurring) during later subsampling. The initial metastable equilibrium in their experiment was attributed to rapid equilibration between the solution and surface-bound copper. At longer time scales, true equilibrium between solution and matrix-bound copper was believed to be occurring. In this study, no distinction between surface-and matrix-bound copper could be made because observed distribution coefficients continually decreased with time.
The desorption rate constant calculated in this study for copper (0.47 d -1 ) is also larger than the upper and lower estimated values of Gee and Bruland (2002) (0.43 d -1 and 0.12 d -1 ). Differences between the upper and lower bounded values of Gee and Bruland (2002) were attributed to the inclusion and exclusion of the first datum at 20 min, respectively. They hypothesized that during the first 20 min of their experiments, the copper spike overwhelmed the capacity for organic ligands to chelate the metal. Consequently, elevated concentrations of ionic copper rapidly adsorbed to the particulate matter, producing higher adsorption rates (and correspondingly higher desorption rates) than would be expected to naturally occur in San Francisco Bay.
Although the initial adsorption process shown on Figure 2 could be explained by a limited amount of organic matter available for complexing, the rates calculated from the desorption phase of the experiment should be representative of natural conditions. First of all, the distribution of copper on particulate surfaces is likely controlled by surface reactions involving copper ions. Secondly, the desorption rate constant in this study (0.47 d -1 ) is less than that obtained by Millward and others (1992) for the uptake of dissolved copper by solid particulate matter from the Tamar estuary (about 0.7 d -1 ). This indicates that the rates used in this study fall within the range of reported values.

Model Calibration
The first step in calibrating the numerical model was to properly define the system to be modeled and adequately identify the important processes to be included. For this investigation, the model domain was the San Francisco Bay estuary from the delta to the Pacific Ocean, divided into 200 x 200-m rectangular grid cells ( Figure 1 shows the final bathymetry). The processes included in the model were tides, wind, waves, erosion, deposition, diffusion, adsorption, and desorption. In addition, loading from major watersheds draining to the bay and delta was important for sediment, salt, and copper. Calibration procedures for the hydrodynamic, advection-dispersion, and sediment transport components are described in Appendix A, and only the metals component is described below. Because one modeling objective was to understand seasonal differences in copper cycling, model scenarios included a wet season during a wet water year and a dry season where sufficient dissolved copper calibration data were available. The 1996 wet season and 1997 dry season were subsequently selected due to their very different hydrographs ( Figure 3).
Boundary conditions and loads for MIKE 21 ME are listed in Table 2. At each boundary, both the dissolved and the adsorbed copper concentrations were specified. Dissolved and adsorbed metal concentrations at the Pacific Ocean and delta boundaries were calculated from equilibrium distribution coefficients, either adsorbed metal concentrations on the suspended sediment or total water column concentrations, and total suspended sediment concentrations (set at 20 mg L -1 at the Pacific Ocean boundary [SFEI 1994[SFEI -1999 and as a time series at the delta [USGS 2000]). Distribution coefficients and metal concentrations on suspended sediment for the Delta boundary were derived from the average of all measurements from the San Joaquin and Sacramento Regional Monitoring Program stations for 1993-1998 (SFEI 1994(SFEI -1999. Distribution coefficients for the Pacific Ocean boundary were similarly obtained from the Golden Gate station. Finally, adsorbed and dissolved concentrations for the Pacific Ocean boundary were calculated from total watercolumn metal concentrations measured in the northern Pacific Ocean (Burton and Statham 1990).
Thirty-seven industrial and wastewater treatment plants (publicly owned treatment works [POTW] sources) were also included in the model. Flows, suspended sediment, and total metal loads for POTWs were obtained from 1996 and 1997 National Pollutant Discharge Elimination System (NPDES) selfmonitoring reports. Sixty-nine tributary flows were calculated for each simulated water year by first delineating individual watershed drainage areas in San Francisco Bay. USGS stream gauges in a number of creeks were then used to estimate flows in nearby streams by normalizing flows by watershed area. Because suspended sediment and metal concentrations are generally higher during storm events, different suspended sediment and metal concentrations were used in this study during high and low tributary flows. High flows were defined as flow rates greater than twice the July and August base flow for a given year. Concentrations of suspended sediment  (BASMAA 1996). Concentrations of suspended sediment and copper for storm events were similarly obtained for each watershed from a land-use summary of these data .
Because the MIKE 21 ME module explicitly accounts for resuspension and transport of benthic sediment with specified chemical concentrations, initial benthic sediment concentrations had to be specified to the model. Considering that copper concentrations on suspended sediment are relatively constant in the bay (Figure 4), an initial benthic copper concentration of 52 mg kg -1 was used in most grid cells. In cells where toxic "hotspots" have been identified (SFBRWQCB 1999), measured values were substituted. Initial porewater concentrations were set to 8.9 µg L -1 , based on the average of concentrations measured at the PW stations shown on Figure 1. This concentration falls within the range of measurements of Rivera-Duarte and Flegal (1997). . Correlation between suspended sediment concentrations and adsorbed copper concentration in water for different regional monitoring stations (SFEI 1994(SFEI -1998. The slope of the correlation line is equal to the copper concentration on sediment (52.4 mg kg -1 )  The desorption rate constant determined in this study was used as direct input to the model based on the observation that desorption rate constants typically vary within a restricted range of 0.1 to 1.0 d -1 (Nyffeler and others 1984), but adsorption rate constants have been shown to vary several orders of magnitude (Glegg and others 1988). Values for adsorption rate constants for different segments of the bay were obtained through an iterative process where temperature and salinity forcing function changes necessitated modifications to the adsorption rate constants and subsequent revisions to the forcing functions. Calibration progressed first by noting that salinity is relatively constant in the South Bay (20 to 25 ppt), but temperature varies both spatially and temporally in this region. This observation supported calibrating the temperature function to achieve seasonal dependence of metal concentrations in the lower South Bay. Salinity was calibrated by isolating the section of the bay with the greatest salinity gradient (i.e., Suisun-San Pablo Bay). The final adsorption coefficient used for the lower South Bay in the wet season was the same value that was determined from the laboratory desorption experiment.

NORTH BAY
The lowest adsorption rate constants used in this study were for the lower South Bay (Table 2). This is consistent with higher concentrations of organic ligands increasing copper solubility (Kuwabara and others 1989). The lower adsorption rate constant required to calibrate the dry season simulation relative to the wet season in the lower South Bay is evidence that when POTW inflows are the largest source of freshwater, organic copper complexing becomes increasingly important (Sedlak and others 1997).

Model Performance
Measured concentrations used for comparison to model predictions included grab samples from the San Francisco Bay Estuary Institute (SFEI) Regional Monitoring Program (RMP) (SFEI 1999) and the City of San Jose monitoring data for the lower South Bay (City of San Jose 1999). Statistics on predicted and measured concentrations are compared for the calibration years (1996 and 1997) and verification years (1993 and 1998) in Table 3. In general, results of the calibration show good agreement with the measured data. Laboratory precision reported by the RMP was ±15%, and was measured at 15%. Average instantaneous values in the North and Central bays agree to within 0.1 µg L -1 for all calibration and verification years, and average instantaneous differences are within 0.3 µg L -1 (Table A-3). For the South Bay, average instantaneous values are within 0.2 µg L -1 for all water years, and average instantaneous differences are within 0.2 µg L -1 . Finally, average instanta- neous values in the lower South Bay agree to within 0.2 µg L -1 , and average instantaneous differences are within 0.5 µg L -1 for all water years except the 1993 dry season. The lower South Bay datum for the 1993 period is 0.8 µg L -1 higher than predicted, although the range of the data is comparable. This difference may be related either to uncertainty in measurement time (±1 h) or the sampling interval of the MIKE 21 ME calculations (output was only saved every 3 h). There is considerable variation in the lower South Bay, with the model being able to capture the general trend and average value, but not the absolute variation.

Dissolved Copper Concentrations
Predicted and measured dissolved copper concentrations at 4 of the 18 monitoring stations shown on Figure 1 are within 0.5 µg L -1 for all data ( Figure 5). The highest average dissolved copper concentrations during the wet season occur in the lower South Bay (2.2 µg L -1 ) (Table 4), where there is a strong correlation between concentration and tributary loading. For example, the large storm event around February 25 (Figure 3) produces an increase in concentration of 1.5 to 2.0 µg L -1 that remains in the lower South Bay throughout the simulation ( Figure  5). This increase in copper concentration is not noticeably propagated to the South Bay monitoring station. Average dissolved concentrations in the South Bay and North Bay are comparable to one another (1.7 to 1.8 µg L -1 ) ( Table 4), except the North Bay station is more heavily influenced by changes in San Joaquin-Sacramento Delta loads (e.g., the February 25 storm event produces an increase in concentration of approximately 0.5 µg L -1 ). Finally, the lowest average dissolved concentrations of dissolved copper occur in the Central Bay. Short-term tidal oscillations in copper concentrations in this region of the bay cause variations as high as 0.5 µg L -1 . By contrast, oscillations of 1 µg L -1 are predicted in the lower South Bay.
Average dissolved concentrations during the dry season simulation are lower than the wet season for the North and Central bays, but are higher for the South and lower South bays ( Table 4). Part of the trend in the North and Central bays can be attributed to declining delta inflow and greater intrusion of lowerconcentration Pacific Ocean water. In the South Bay, increased amounts of organic material available for complexing may result in greater partitioning into the dissolved state (Sedlak and others 1997).
Two plots of dissolved copper concentration are displayed on Figure 6 for the wet and dry season to illustrate the effect of tributary loads and wind-wave resuspension of benthic sediment, respectively. The plot on the top of the figure is the concentration during the February 25 storm event. Concentrations are generally between 1.0 and 1.5 µg L -1 outside of the Golden Gate, 1.5 to 2.0 µg L -1 in much of the bay, and greater than 5.0 µg L -1 in the lower South Bay and in the immediate vicinity of watershed tributary sources. The plot on the bottom shows concentrations during a strong wind-wave event causing resuspension of sediment in the shoals of the bay. Of importance in this plot are the smaller contribution of tributary sources and the greater intrusion of the Pacific Ocean into the Central Bay. Also, a net translation of dissolved copper from the lower South Bay along the main channel of the bay is illustrated in the figure. Finally, there is an increase in dissolved concentration on the order of 0.5 µg L -1 in the area of high resuspension at the eastern shoals of the South Bay.

Copper Distribution Coefficients
Hypothetical equilibrium distribution coefficients reported in Table 4 were calculated from adsorption and desorption rate constants (Table 2), assuming a salinity of 20 parts per thousand (ppt) and temperature of 17°C. The highest distribution coefficients are predicted to be in the Central Bay, a result consistent with the lower concentrations of dissolved organic matter available for complexing (Sedlak and others 1997). Assuming the difference in average distribution coefficients between each adjoining segment of the bay is approximately 5,000 L kg -1 , the time for equilibrium to be reached between water from different segments is about two days, with adsorption occurring at a rate that is several hours faster than desorption.
Calculated equilibrium distribution coefficients and model-predicted distribution coefficients for four monitoring stations are displayed on Figure 7 as dark and light lines, respectively. Where predicted distribution coefficients are larger than equilibrium coefficients, net desorption is predicted to occur. Where these values are smaller, net adsorption is predicted. In the North Bay, equilibrium distribution coefficients are near 40,000 L kg -1 during the two months of the wet season simulation, and drop as much as 10,000 L kg -1 during the end of the simulation. Predicted distribution coefficients oscillate near equilibrium for both simulations, although there is net desorption at the beginning of the wet season and throughout the dry season. In

WET SEASON 1997 DRY SEASON
the Central Bay, equilibrium distribution coefficients are greater than 40,000 L kg -1 in the wet season and 30,000 L kg -1 in the dry season. In contrast to the North Bay, net adsorption is predicted to be occurring. In the South Bay, predicted distribution coefficients are greater than 30,000 L kg -1 in the wet season and 20,000 L kg -1 in the dry season, and are close to equilibrium in both. Finally, in the lower South Bay, predicted distribution coefficients are greater than 20,000 L kg -1 in the wet season and are near 10,000 L kg -1 during the dry season. The drop in distribution coefficients during February 1996 in North and lower South bays is consistent with the large influx of lower salinity freshwater tributary loads and the salinity dependence determined from calibration.
Two plots of predicted distribution coefficients are displayed on Figure 8 for the same storm and wind-wave events of Figure 6. Predicted distribution coefficients generally decrease from the Pacific Ocean into the bay during the wet season simulation, with the smallest values found within several kilometers of tributary sources. For the dry season simulation, the lowest distribution coefficients are found throughout lower South Bay. Because the plots represent ebb tides, there is a net transfer of dissolved copper from the North and South bays to the Central Bay, and from the lower South Bay to the South Bay. These plumes of dissolved copper produce predicted distribution coefficients that are lower along the main channel of the estuary because there is insufficient time for equilibration to occur (equilibration rates are on the order of days [ Table 4], but tidal oscillations are on the order of hours).
Ratios between predicted and equilibrium distribution coefficients are shown on Figure 9. Values greater than 1.0 indicate net desorption, and values less than 1.0 net adsorption. During the February 25 storm event (upper plot), the area of the bay undergoing net adsorption is larger than the area undergoing desorption. Areas of adsorption are found along the main channel of the Central and South bays, where dissolved copper is transferred from adjoining segments. Adsorption also occurs near the mouths of tributaries with lower equilibrium distribution coefficients than the receiving water of the bay. This is most noticeable in the figure at the mouths of the Petaluma, Napa, Sonoma, and Alameda tributaries. Although net desorption is evident in the lower South Bay grid cells where tributary loads are located, net adsorption is observed in surrounding cells. This result is consistent with the relatively lower equilibrium distribution coefficients of the lower South Bay compared to the watershed sources, but the greater potential of dissolved copper introduced by the tributaries to be advected.
In contrast to the wet season storm event, the area undergoing desorption during the dry season is relatively larger (lower plot). Desorption during this period is greatest at the boundary between the South and lower South bays, consistent with increased organic complexing that is expected to occur. The magnitude of the change shown on the figure, however, is likely an artifact of delineating sharp adsorption rate constant boundaries during calibration (actual desorption processes in the bay are more gradual than modeled). Of equal importance on the figure is that net desorption occurs in the shoals and net adsorption in the channel. This suggests that resuspended sediment acts as a source of dissolved copper and that there is a mass flux of dissolved copper toward the main channel of the bay. 5.5e+004 5e+004 -5.5e+004 4.5e+004 -5e+004 4e+004 -4.5e+004 3.5e+004 -4e+004 3e+004 -3.5e+004 2.5e+004 -3e+004 2e+004 -2.5e+004 1.5e+004 -2e+004 1e+004 -1.5e+004 5000 -1e+004 1000 -

Copper Mass Fluxes
Mass fluxes of copper were calculated between dissolved and adsorbed phases within and between each of the bay segments depicted on Figure 1. This process consisted of first calculating desorptive and diffusive fluxes for each grid cell using model output and the respective components of equation (6). Values were subsequently summed over each bay segment for each timestep to obtain net desorption and diffusion mass transfer rates. Finally, mass fluxes between segments were determined from mass balance calculations using the following: (1) model output of dissolved and adsorbed copper in the water column and benthos in each grid cell for each timestep; (2) calculated sorption-related and diffusive mass exchange; and (3) copper mass loading due to tributary and POTW sources.
Net copper mass flux rates are reported in Table 5 in units of 10 3 kg yr -1 . In the North Bay, there is no net change in dissolved copper concentration during the wet or dry seasons. The largest source of dissolved copper is the San Joaquin/Sacramento Delta (104 x 10 3 kg yr -1 ), followed by tributary loads during the wet season (22 x 10 3 kg/year) and diffusive loads during the dry season (0.8 x 10 3 kg yr -1 ). Most dissolved copper introduced into the North Bay is advected to the Central Bay, with desorption being a minor process. In the Central Bay, there is predicted to be a net loss of dissolved copper during the wet season (3 x 10 3 kg yr -1 ), and no net change during the dry season. Most copper exchanged in this region of the bay is between the Pacific Ocean (-142 x 10 3 kg yr -1 ) and North Bay (129 x 10 3 kg yr -1 ). Net adsorption occurs (3 x 10 3 kg yr -1 ); however, its relative importance compared to other processes is minor. In the South and lower South bays, there is a net accumulation of dissolved copper for the seasons modeled (1 to 20 x 10 3 kg yr -1 ). During the wet season, tributary inflows are the largest source of dissolved copper (4 to 10 x 10 3 kg yr -1 ). During the dry season, diffusion is most important in the South Bay (8 x 10 3 kg yr -1 ), and POTW loads are most important in the lower South Bay (1 x 10 3 kg yr -1 ). A distinguishing feature of the latter two regions is that net adsorption occurs during the wet season (0.4 to 3 x 10 3 kg yr -1 ) and desorption during the dry season (0.1 to 0.9 x 10 3 kg yr -1 ). The relative input of dissolved copper from desorption during the dry season is approximately 7% of the total for the South and lower South bays.

Copper Residence Times
The residence time of copper in different bay segments was calculated using a simple box model from average seasonal mass inputs and fluxes. This method assumes mass inputs are constant in hydraulic flow rate, dissolved concentration, suspended sediment concentration, and adsorbed concentration. In addition, the mass input method assumes mass within the system is constant while calculating the time required to replenish the existing reservoir. A simplified expression for residence time using the mass input method is shown below: where Total Mass = mass in the active sediment layer plus mass in the water column, and Mass Input = source loading from POTWs, major creeks, tributaries, and adjacent bay segments. Chemical mass in the water column and benthic sediment was calculated within every model grid cell at every time step. In order to perform the calculations, an active sediment layer, defined as the approximate depth of resuspendable sediment within the basin of interest, had to be specified. The active layer was initially defined as the top 10 cm of bay sediment. This value was selected based on the upper limits of erosional depths in the South Bay obtained from

Dissolved Copper Cycling
Dissolved copper predominantly enters the North Bay from the San Joaquin-Sacramento Delta system and exits through the Central Bay without accumulation (Table 5). Because the average distribution coefficient assigned to the delta is sometimes above and sometimes below calculated North Bay equilibrium coefficients during the wet season, sorption processes are predicted to be relatively unimportant to the mass balance of dissolved copper. The primary effect of adsorption and desorption processes in the North Bay is to redistribute copper introduced by tributaries. Because the equilibrium distribution coefficient of the delta was modeled as being larger than the average North Bay coefficient during the dry season, it might be anticipated that net desorption would be observed in the channel. As shown on Figure 9, the opposite trend is exhibited. This indicates that there is a net source of dissolved copper to the channel during times of high sediment resuspension, either from North Bay tributaries or copper desorbed from resuspended sediment in the adjacent shoals.
The primary source and sink of dissolved copper in the Central Bay is the North Bay and Pacific Ocean, respectively (Table 5). Although there is predicted to be a net loss of dissolved copper during the wet season, the decrease is roughly equivalent to the amount of copper adsorbed. Consequently, even though adsorption only represents a small fraction of the total copper mass during a given season, it may be important to the redistribution of copper between bay segments. The relatively constant adsorbed copper concentrations shown on Figure 4 are empirical evidence that some redistribution occurs.
The importance of adsorption and desorption to the mass balance of dissolved copper is greatest in the South and lower South bays (Table 5). The primary source of dissolved copper during the wet season is tributary loads, and the primary sink is adsorption onto suspended sediment. This dependence is evident on Figure 9, where despite the fact that there is some desorption of copper near segment boundaries and lower South Bay tributaries, adsorption occurs over a relatively larger area. During the dry season, the primary source of dissolved copper is predicted to be diffusion and POTW flows for South and lower South bays, respectively. The primary sink for the two regions is exchange with the Central Bay and South Bay, respectively. The importance of adsorption and desorption during the dry season is shown on Figure 9, where there is an apparent mass transfer of dissolved copper to suspended sediment in the South Bay channel from the South Bay shoals and lower South Bay.

Implications for Water Quality
Because adsorbed copper concentrations are relatively uniform in the bay (Figure 4), dissolved copper concentrations are predominantly controlled by seasonally and spatially varying changes in water chemistry related to changes in salinity and the availability of organic ligands in different segments of the bay. Sorption processes related to disequilibria continuously occur between bay segments and episodically occur during wet season storm events or dry season wind-wave events (it should be noted that disequilibria between bay segments is more gradual than can currently be modeled; some of the disequilibria is caused by applying different rate constants to different bay segments). The net result of episodic events is to redistribute copper that is initially sedimentbound (i.e., adsorbed). This redistribution is in the direction of the main channel of the bay, away from tributary and POTW sources.
If mass loading rates were relatively constant, the bay would eventually reach a steady state condition where the total mass of copper is unchanged. During an "average" year, lower equilibrium distribution coefficients in the dry season would result in greater partitioning and transport of dissolved copper. Dissolved copper removed during this season would later be replenished during wet season storm events. Under steady state conditions, the immediate effect of reducing dissolved copper concentrations in tributary or POTW loads would be to raise the equilibrium distribution coefficient in the source water. This would subsequently cause greater desorption from the sediment during initial mixing with the receiving water. Because dissolved copper is more mobile, it may be expected that immediate improvements in water quality would be possible; however, results from the lower South Bay indicate that if dissolved concentrations are higher than ambient water, net adsorption occurs close to the source. In other words, suspended sediment provides a buffering effect against changes in dissolved copper concentrations. Given that residence times of copper adsorbed to sediment are on the order of years, water quality improvements may be gradual.

Model Uncertainty
Sensitivity runs reported in Appendix B indicate that the magnitude of tributary source loads affects predicted dissolved copper concentrations. This is particularly true during the wet season, when tributary sources are a larger fraction of the total mass of copper introduced to the water column. Because this effect cannot be separated from the calibration, trends and copper cycling processes predicted in this study may not accurately describe copper partitioning near tributaries.
Another source of uncertainty in the model is the adsorption and desorption rate constants selected during calibration. Initial estimates were based on a single experiment and are only valid for the set of chemical conditions present during the experiment. The greatest difference in predictions due to changes in rate constants would be near tributary sources and bay segment boundaries. If the desorption rate was 0.05 d -1 (Glegg and others 1988), then adsorption and desorption processes would take longer than residence times of dissolved copper in some parts of the bay during the wet season and all segments of the bay in the summer. Consequently, the bay would be further from equilibrium between the dissolved and adsorbed phases than predicted in this study.
Although the initial porewater concentration does not greatly affect the predicted dissolved concentration (Appendix B), porewater concentrations can affect the predicted importance of sorption processes discussed above. Assuming that porewater copper concentrations are approximately equal to water column concentrations over long time periods (data from Rivera-Duarte and Flegal [1997] suggest that benthic mass fluxes are sometimes positive and sometimes negative), and assuming copper partitioning is predominantly governed by an equilibrium distribution coefficient (Figure 7), the importance of desorption would be greater than predicted in this study. This is because a smaller (or negative) diffusive flux would require more copper desorption from suspended sediment to maintain equilibrium. The contributions of desorption to the total copper flux for the case where diffusive fluxes are replaced by desorption was determined from a sensitivity run that excluded diffusion. The resulting net contribution to the total mass input of dissolved copper was found to be 20% and 8% for the South and lower South bays, respectively.

CONCLUSIONS
Model calibration shows that dissolved copper concentrations in the bay generally oscillate near concentrations predicted from seasonally varying equilibrium distribution coefficients. Using adsorption and desorption rates calculated in this study as initial estimates, calibration was also only found to reproduce measured dissolved concentrations when different adsorption rates were specified for distinct segments of the bay. This result is consistent with the relatively larger amounts of organic ligands available for complexing in the South and lower South bays (Kuwabara and others 1989), although differences between Bay segments is more gradual than can be modeled. During wet-season storm events or dry-season wind-wave events, adsorption or desorption may locally alter dissolved concentrations. This effect is particularly pronounced in the lower South Bay, where the combination of high sediment resuspension and organic complexing during the dry season produces dissolved concentrations that episodically exceed federal and state water quality standards.
Model results show that the primary sources of dissolved copper to the bay are the tributaries that feed it, particularly those entering through the San Joaquin-Sacramento Delta. The primary sink is the Pacific Ocean, with net fluxes of dissolved copper occurring in this direction. There is some internal cycling of copper in different segments of the bay that redistributes copper initially adsorbed in particulate matter. In the South and lower South bays, net adsorption of copper may periodically occur due to high tributary loading. Some of the copper irreversibly adsorbed during the wet season is later released during wind-wave resuspension events in the dry season. The net result of episodic events on total copper mass balance is to repartition sediment-bound adsorbed copper to the dissolved phase. This redistribution occurs in the direction of the main channel of the bay, away from tributary and POTW sources.