Enhanced Two Dimensional Hydrodynamic and Water Quality Model (CE-QUAL-W2) for Simulating Mercury Transport and Cycling in Water Bodies

CE-QUAL-W2 (W2) is a widely-used two-dimensional, laterally averaged, longitudinal/ vertical, hydrodynamic and water quality model. This model was modified and enhanced to include a mercury (Hg) simulation module for simulating Hg transport and cycling in water bodies. The Hg simulation module in W2 is able to model the physical and biochemical processes including adsorption and desorption of Hg species on multi-solids, settling and resuspension, sediment burial of adsorbed Hg, diffusive exchange between water column and sediment layer, volatilization, and biogeochemical transformations among Hg species. This paper describes the Hg simulation module, W2 model validation and its application to the Xiaxi River, China, a historical Hg contaminated water body. The W2 model was evaluated using the Xiaxi River data collected in 2007 and 2008. Model results show that W2 was able to predict the total Hg and methylmercury concentrations observed for the Xiaxi River. The Xiaxi River W2 model also provides a foundation for the future investigations of Hg contamination in the Xiaxi River. This application demonstrated the W2 model capability in predicting complex transport and cycling of Hg species in water bodies.


Introduction
CE-QUAL-W2 (W2) is a two-dimensional (2D), laterally averaged, longitudinal/vertical, hydrodynamic and water quality model for simulating surface water bodies. Lateral averaging assumes lateral variations in velocities, temperatures, and constituents are negligible. This model has been continuously developed and maintained by the U.S. Army Corps of Engineers (USACE) and Portland State University since late 80s [1,2]. Since its first release in 1986, the W2 model has been successfully applied to hundreds of rivers, lakes, reservoirs and estuaries throughout the world for thermal and water quality investigations [3][4][5][6][7][8]. U.S. governmental agencies including USACE, U.S. Bureau of Reclamation (USBR), U.S. Geological Survey (USGS), U.S. Environmental Protection Agency (USEPA), and Tennessee Valley Authority (TVA) as well as numerous states, county, and local agencies have used this model as a management tool to evaluate the impacts of point and nonpoint pollution sources including organic wastes, nutrients, and temperature in water bodies. The W2 model has also been applied in China [6]. Cole and Wells [2] presented the model theory, fundamental equations, and application guidelines. The model complexity has increased significantly with the increasing complexity of water quality and management issues since its release. The latest version of the W2 model can simulate the complete eutrophication process in water bodies, including carbon, nitrogen and phosphorous cycles, multiple algae, epiphyte and macrophyte species, dissolved oxygen, etc.
Additionally, a benthic sediment diagenesis module was integrated into the model to compute the sediment-water exchange fluxes of nutrients and sediment oxygen demand [9]. Mercury (Hg) is a global, persistent, and bioaccumulative toxic contaminant. Its organic form, methylmercury (MeHg), is the most biologically active and toxic form of Hg. It can accumulate in organisms within the food chain, such as fish, posing a risk to wildlife and humans [10][11][12]. Several models have been developed to simulate the transport and fate of Hg in aquatic environments. These models include the MCM [13], MERC4 [14], SERAFM [15], and WASP [16,17]. Although W2 is a powerful tool and can simulate up to 32 water quality constituents, the current W2 model does not have capability of simulating Hg transport and cycling in water bodies.
To effectively evaluate both the potential for harmful bioaccumulation and proposed measures to mitigate that harm, the spatial and temporal distribution of Hg within a water body should be characterized as completely as possible. There has been a need to bring the W2 to the state-of-art and adding Hg modeling capability. In this study, a Hg simulation module was developed and incorporated into the W2 model. With this enhancement, W2 model can directly simulate the transport and cycling of three Hg species in water bodies. After a series of model testing and validation studies, the enhanced W2 model was applied to the Xiaxi River for further evaluation of the model performance. The Xiaxi River, which drains the Hg mining area in Wanshan County, has been contaminated by numerous Hg tailings. Wanshan County, located in Guizhou Province, south-western China, was identified as a heavy Hg contamination site. The largest conglomeration of Hg mines and refining plants in China is located in this area [18]. Enhanced W2 model with a newly developed Hg simulation module was applied to the Xiaxi River. Besides evaluating the W2 performance, model results of this application were intended to provide insight into the long-term fate of MeHg and Total Hg (THg) along the river and the relative importance of different processes in Hg cycling. This application can also be used as a management tool to quantify the sources of Hg contamination, the distribution of Hg, the influence of changing environmental conditions, and long-term mitigation of Hg contamination within the watershed.

Hg Simulation Module within CE-QUAL-W2
One of the major concerns of Hg and its impact on the environments is its transformation to MeHg. Speciation plays a major role in determining the bioavailability of Hg species for MeHg transformation and toxicity of Hg to living organisms. Speciation also influences transport of Hg within and between environmental media [19][20][21][22]. Hg0 and HgII are the major chemical forms of Hg input into the environment from anthropogenic or natural sources [20]. During the biogeochemical cycling of Hg, MeHg can be produced. Hg biogeochemical cycle involves three species and multi-phases in nature. The Hg simulation module in the W2 model was designed to simulate three primary Hg species of environmental concern: Hg0, HgII, and MeHg. Hg0 is assumed to exist in the dissolved phase and only in the water column. HgII and MeHg are simulated for the water column and sediment layer. The three Hg species and their kinetic equations implemented in the Hg simulation module are described in detail in Zhang and Johnson [23], and briefly summarized here.
Hg0 is slightly water soluble and has a high Henry's Law constant [24]. Hg0 constitutes very little of the total mercury in the surface water but may provide a significant pathway for the volatilization of Hg from surface waters. In the Hg simulation module, Hg0 is assumed to exist in the dissolved phase and only in the water column. MeHg and HgII exist in both the water column and the sediment layer and may distribute in dissolved phase and sorbed with dissolved organic carbon (DOC), particulate organic matter (POM), and inorganic solids. Physical and biochemical rate processes simulated in the Hg simulation module include: (1) adsorption and desorption of Hg species; (2) volatilization; (3) atmospheric deposition; (4) diffusive exchange between the water column and sediment layer; (5) deposition and re-suspension; (6) sediment burial of adsorbed Hg; (7) biogeochemical transformations among the three species. Figure 1 provides an overview of the model representation of the Hg species and their processes. Table 1 lists the symbols and definitions of concentrations of Hg species computed in the Hg simulation module.  The sediment layer represents the upper mixed layer and may be on the order of several centimeters thick due to bioturbation and mixing. Suspended sediments, with Hg bound to them,  The sediment layer represents the upper mixed layer and may be on the order of several centimeters thick due to bioturbation and mixing. Suspended sediments, with Hg bound to them, are deposited to and eroded from the sediment bed. Depending upon ambient concentrations of Hg species, the suspended sediments lost or gained to the water column become a source or sink for Hg. Once deposited, sediments coalesce to form an unconsolidated layer of variable, but often shallow depth. This active layer of sediments may undergo continual bioturbation, or mixing by organisms, so that sediments and Hg are distributed evenly throughout. Deeper sediment is considered biologically inactive, and a static repository for sediments and Hg as they are buried.
As represented in Figure 1, the transport and cycling of Hg species in a water body is a complex interaction of various physical, chemical, and biochemical processes. The physical and biochemical processes modeled in the Hg simulation module are reasonably well documented in Zhang and Johnson [23]. The three Hg species and their kinetic source and sink equations implemented in the W2 model are briefly summarized here.

Hg Partitioning
Hg species may be present in different phases in aquatic environments. In the Hg simulation module, Hg0 is assumed in the dissolved phase, HgII and MeHg are allowed to be partitioned among the dissolved phase and DOC, and solids adsorbed phases. Solids can be assigned to abiotic mineral solids, detrital, or various classes or size categories of solids. Under equilibrium partitioning, HgII and MeHg are assumed to be instantly exchangeable between the dissolved and adsorbed phases. The fractions of each phase are functions of the partition coefficients and the concentrations of DOC, POM and inorganic solids. The linear equilibrium partitioning fractions of HgII and MeHg based on a total volume basis in the water column are calculated as follows.
Sediment HgII and MeHg are also allowed to be partitioned among dissolved, DOC, POM and inorganic solids sorbed phases. The partitioning fractions of HgII and MeHg in the sediment layer are identically calculated as that for the water column.
The symbols in Equations (1)-(10) are defined in Table 2. Hg0 oxidation yield coefficient settling velocity of inorganic solid "n" m·d −1 v rn sediment resuspension velocity of solid "n" m·d −1 v som settling velocity of POM m·d −1

Hg Transformations and Kinetic Equations
The oxidation, reduction, methylation and demethylation transformations of Hg species are assumed to be widespread in aquatic environments. Since some of the reaction processes have not been well understood and defined, Hg transformation reactions are simulated using first order rate kinetics. These reactions may act on Hg depending on the speciation and partitioning of Hg. To account for this, the base rate of reaction is modified by the fraction of Hg dissolved in the aqueous phase, complexed with DOC, and adsorbed to abiotic and biotic particles. As represented in Figure 1, major processes simulated in the Hg simulation module also include partitioning of HgII and MeHg, particulate settling, resuspension and burial, sediment-water diffusion, and volatilization. Meanwhile, Hg species can be transported together with water, DOC, and suspended solids. The total internal source and sink terms of three species in the water column are calculated as follows.
Hg0 oxidation +1.33 HgI I reduction +1.33 MeHg photolytic degradation (11) Hg0 oxidation +1.33 MeHg sediment−water trans f er HgI I methylation −1.33 Hg settles through water column and deposits to or erodes from benthic sediments. Within the sediment bed, biotic and abiotic reactions can result in methylation, demethylation, and reduction, just as in the water column. However, because of varying environmental conditions, the rate of these reactions can vary dramatically from those in the water column. The sediment layer in the Hg simulation module is assumed to have constant properties including the thickness, volume, porosity, and bulk density. Only HgII and MeHg are simulated for the sediment layer. Concentrations of Hg species are expressed in terms of mass per unit volume of total sediments (ng·L −1 ). The total internal source and sink terms of HgII and MeHg in the sediment layer are calculated as follows.
Within the sediment layer, dissolved species migrate downward or upward through percolation and pore water diffusion. The sediment-water flux of dissolved Hg species included in Equations (12) to (15) is the product of effective mass transfer velocity and the difference between bulk fluid concentrations on either side of the interface. The mass transfer velocity between each modeled species in the sediment layer and the water column is calculated in the W2 model. The following four equations are included in the Hg simulation module for calculating the sediment-water transfer of dissolved Hg species: (1) Thibodeaux et al. [25]; (2) Di Toro et al. [26]; (3) Boyer et al. [27]; and (4) Schink and Guinasso [28].
where D b is biodiffusion coefficient representing particle diffusivity (m 2 ·d −1 ), ρ b is sediment bulk density (g·cm −3 ), v m is user-defined sediment-water transfer velocity (m·d −1 ), z 2 is sediment bioturbation layer thickness (m), MW is molecular weight of Hg species (g·mol −1 ), D m is molecular diffusion coefficient (m 2 ·d −1 ), z l is pore water diffusion layer thickness (0.01 m), u * is flow shear velocity along the bed, which is approximately 10 percent of the mean velocity of flow (m·s −1 ), v is kinematic viscosity of water (m 2 ·d −1 ).

Enhanced CE-QUAL-W2 Model
The W2 model uses fixed computation grids with a static bathymetric surface onto which longitudinal segments and vertical layers are mapped; the hydrodynamic and water quality computations are performed at the intersections of these segments and layers. A direct coupling between hydrodynamic and water quality simulations is implemented in W2 model. The 2D water quality component computes the transport and fate of constituents with their kinetic reaction rates expressed in source and sink terms. The general 2D water quality transport and mass balance equation of W2 model is written as follows [2].

∂(BΦ) ∂t
where Φ is laterally averaged constituent concentration (g·m −3 ), B is control volume width (m), U and W are longitudinal and vertical velocities respectively (m·s −1 ), D x is longitudinal temperature and constituent dispersion coefficient (m 2 ·s −1 ), D z is vertical temperature and constituent dispersion coefficient (m 2 ·s −1 ), q Φ is lateral inflow or outflow mass flow rate of constituent per unit volume (g·m −3 ·s −1 ), S Φ is laterally averaged source/sink term (g·m −3 ·s −1 ). The sources/sinks for each water quality state variable in the W2 model are separated into two arrays, one contains boundary sources/sinks, and the other one contains internal sources/sinks due to kinetic interactions. Figure 2 provides the overall flow chart of the W2 model solving the transport and mass balance of three Hg species.

Enhanced CE-QUAL-W2 Model Validation and Evaluation
A comprehensive validation of the newly developed Hg simulation module in the W2 model would require an extensive testing and validation. In this study, enhanced W2 model was validated through a series of case studies and a real world application. The Hg simulation module developed here includes similar Hg cycling processes as WASP model. The WASP model is also capable of simulating the transport and cycling of three Hg species (Hg0, HgII and MeHg) in water bodies. The model results generated from WASP were first used to validate individual Hg processes simulated in the W2 model. This model validation was designed to assure that the algorithms and formulations are implemented correctly in the enhanced W2 model. The enhanced W2 model was also applied to the Xiaxi River in Wanshan County, China. Wanshan County in east Guizhou

Enhanced CE-QUAL-W2 Model Validation and Evaluation
A comprehensive validation of the newly developed Hg simulation module in the W2 model would require an extensive testing and validation. In this study, enhanced W2 model was validated through a series of case studies and a real world application. The Hg simulation module developed here includes similar Hg cycling processes as WASP model. The WASP model is also capable of simulating the transport and cycling of three Hg species (Hg0, HgII and MeHg) in water bodies. The model results generated from WASP were first used to validate individual Hg processes simulated in the W2 model. This model validation was designed to assure that the algorithms and formulations are implemented correctly in the enhanced W2 model. The enhanced W2 model was also applied to the Xiaxi River in Wanshan County, China. Wanshan County in east Guizhou Province has a long history of Hg mining and processing. Multi million tons of calcine from the roasted ore was deposited in this area from 1949 to early 1990s [29]. Approximately 4.4 g Hg per kg of calcines were detected in Wanshan County [30]. Although Hg mining and retorting activities have been banned since 2001, the field investigation found the Hg leak from historic deposits of calcines and Hg contamination in the Xiaxi River was detected. The application of enhanced W2 model to the Xiaxi River was intended, in part, to evaluate the merit of the W2 model for simulating the transport and cycling of Hg species in a water body.

Model Validation through Comparing WASP Model Results
In this study, an equivalent segment was set up for W2 and WASP. The two models were run with, to the extent possible, a consistent set of kinetic coefficients and parameters, initial conditions, fluxes from the water column and water column concentrations. A user defined sediment-water transfer velocity was used in both models because different approaches are included. The same Hg partitioning coefficients on solids were specified for the water column and sediment layer due to the limitation of the WASP model. In the W2 model, the impact of sulfate reduction on methylation rate in the sediment layer is simulated, which has been reported in many literatures [31][32][33][34]. However, the WASP model does not simulate the impact of sulfate on the methylation. In the test case, sulfate reduction rate in W2 was specified to ensure that the same methylation rate (k d23_2 = 0.01 d −1 in Table 3) was used in both models. The WASP model simulates the demethylation from MeHg to HgII in the water column process through a bacterial degradation process, but this process is simulated in the W2 model through a photolytic degradation. Therefore, light impacts on the demethylation in this test case were ignored in the W2 model. The photolytic degradation of MeHg has been reported in many literatures [35][36][37]. Meanwhile, WASP model always needs to be coupled with other hydrodynamic models for water quality investigations, such as EFDC [38,39].
In this test case, the water column was set as 100 m long, 10 m wide with an average of 2.5 m deep. The water temperature was set as 25 • C. The solar radiation was set as 500 W·m −2 , and light extinction coefficient was 1.0 m −1 . Both inorganic and organic solid classes were simulated in both models. The suspended concentrations of silt and fines, sand and organic matter were defined as 20, 50 and 15 mg·L −1 respectively. The bottom sediment layer was made up of four solid classes with a fixed porosity of 0.7, on average, 83.33% sand, 10.71% silt and fines, and 5.96% organic matter. The sediment temperature was set as 20.0 • C. The settling velocities for silts and fines, and organic matter were respectively set as 0.5, 1.2 and 0.3 m·d −1 . The resuspension velocities for silts and fines, and organic matter were set as 1.1 × 10 −4 , 8.0 × 10 −5 and 9.0 × 10 −4 m·d −1 . Burial velocity was set as 5.03145 × 10 −6 m·d −1 . The sediment-water transfer velocity was specified as 0.0864 m·d −1 . The initial concentrations of Hg0, HgII and MeHg in the water column were defined as 1.0, 10.0 and 0.0 ng·L −1 , and the initial concentrations of HgII and MeHg in the sediment layer were specified as 50.0 and 0.0 ng·g −1 . Most of data in this test case was based on observed data in Lake Taihu in China [40]. The other water quality parameters are given in Table 3. These parameters were estimated from literature values provided by Knightes et al. [17], Allison and Allison [41] and Muresan et al. [42]. Both W2 and WASP models were ran with a fixed time step of 0.1 days for a year simulation period. Table 3. List of water quality parameters used in the test cases.

Variable
Units Value  Table 3. Cont.

Variable Units Value
The test case study was used to validate the performances of the W2 model with respect to various aspects such as multi-phase partitioning, methylation, demethylation transformations, and sediment-water transfer. Figure 3a,b show the comparison of W2 and WASP computed concentrations of HgII in the water column and sediment layer. Figure 4a,b show W2 computed HgII internal flux by each process. Initially, total concentration of HgII in the water column increased mainly by the sediment resuspension, since this flux remained at a high level (nearly 1.8 ng·L −1 ·d −1 ), even though HgII settling increased gradually during the simulation period. Total concentration of HgII in the water column reached the equilibrium after 20 days. In this test case, HgII settling and resuspension were the two main internal sink and source in the water column. Hg0 oxidation, HgII reduction, HgII methylation and sediment-water transfer rates were relatively small. Total concentration of HgII in the sediment layer was more stable for the simulation period. The sediment layer become a large storage reservoir for Hg.
Water 2017, 9,643 12 of 23 The test case study was used to validate the performances of the W2 model with respect to various aspects such as multi-phase partitioning, methylation, demethylation transformations, and sediment-water transfer. Figure 3a,b show the comparison of W2 and WASP computed concentrations of HgII in the water column and sediment layer. Figure 4a,b show W2 computed HgII internal flux by each process. Initially, total concentration of HgII in the water column increased mainly by the sediment resuspension, since this flux remained at a high level (nearly 1.8 ng·L −1 ·d −1 ), even though HgII settling increased gradually during the simulation period. Total concentration of HgII in the water column reached the equilibrium after 20 days. In this test case, HgII settling and resuspension were the two main internal sink and source in the water column. Hg0 oxidation, HgII reduction, HgII methylation and sediment-water transfer rates were relatively small. Total concentration of HgII in the sediment layer was more stable for the simulation period. The sediment layer become a large storage reservoir for Hg.
The test case study was used to validate the performances of the W2 model with respect to various aspects such as multi-phase partitioning, methylation, demethylation transformations, and sediment-water transfer. Figure 3a,b show the comparison of W2 and WASP computed concentrations of HgII in the water column and sediment layer. Figure 4a,b show W2 computed HgII internal flux by each process. Initially, total concentration of HgII in the water column increased mainly by the sediment resuspension, since this flux remained at a high level (nearly 1.8 ng·L −1 ·d −1 ), even though HgII settling increased gradually during the simulation period. Total concentration of HgII in the water column reached the equilibrium after 20 days. In this test case, HgII settling and resuspension were the two main internal sink and source in the water column. Hg0 oxidation, HgII reduction, HgII methylation and sediment-water transfer rates were relatively small. Total concentration of HgII in the sediment layer was more stable for the simulation period. The sediment layer become a large storage reservoir for Hg.     MeHg main source was from the methylation of HgII in the water column. Settling and resuspension were not the dominant processes affect the fate of MeHg in the water column. In the sediment layer, main sources of MeHg come from MeHg deposition from water column and HgII methylation. It is important to understand the sources of MeHg in order to control the production of MeHg in the water column and sediment layer and reduce the environmental risk. To investigate individual Hg transformations, solid sorption, settling, resuspension, and sediment-water transfer processes were turned off. Figure 7a,b show dissolved concentrations of HgII and MeHg in the water column. Figure 8a,b show Hg transformation fluxes respectively. W2 predicted concentrations of HgII and MeHg are consistent with WASP model results. Dissolved HgII was lost by its reduction and methylation in the water column. In the meantime, MeHg was gained and reached the peak around 50 days, then decreased by its demethylation and photoreduction. The above case studies indicate that the algorithms and formulations were correctly implemented in the enhanced W2 model.  To investigate individual Hg transformations, solid sorption, settling, resuspension, and sediment-water transfer processes were turned off. Figure 7a,b show dissolved concentrations of HgII and MeHg in the water column. Figure 8a,b show Hg transformation fluxes respectively. W2 predicted concentrations of HgII and MeHg are consistent with WASP model results. Dissolved HgII was lost by its reduction and methylation in the water column. In the meantime, MeHg was gained and reached the peak around 50 days, then decreased by its demethylation and photoreduction. The above case studies indicate that the algorithms and formulations were correctly implemented in the enhanced W2 model.  To investigate individual Hg transformations, solid sorption, settling, resuspension, and sediment-water transfer processes were turned off. Figure 7a,b show dissolved concentrations of HgII and MeHg in the water column. Figure 8a,b show Hg transformation fluxes respectively. W2 predicted concentrations of HgII and MeHg are consistent with WASP model results. Dissolved HgII was lost by its reduction and methylation in the water column. In the meantime, MeHg was gained and reached the peak around 50 days, then decreased by its demethylation and photoreduction. The above case studies indicate that the algorithms and formulations were correctly implemented in the enhanced W2 model. To investigate individual Hg transformations, solid sorption, settling, resuspension, and sediment-water transfer processes were turned off. Figure 7a,b show dissolved concentrations of HgII and MeHg in the water column. Figure 8a,b show Hg transformation fluxes respectively. W2 predicted concentrations of HgII and MeHg are consistent with WASP model results. Dissolved HgII was lost by its reduction and methylation in the water column. In the meantime, MeHg was gained and reached the peak around 50 days, then decreased by its demethylation and photoreduction. The above case studies indicate that the algorithms and formulations were correctly implemented in the enhanced W2 model.

Study Area
The Xiaxi River is located in Wanshan County, Guizhou Province, China. This river has been influenced by the world-famous Wanshan Hg mining, which is one of the thirteen large and super large-scale Hg mines located within the Guizhou Province, China. Figure 9 shows the extent of the Xiaxi River, and its location in Guizhou Province. The Hg tailings are distributed throughout the headwater reaches of the Xiaxi River. Due to the high geochemical background of Hg in this area and high Hg emissions from coal combustion using Hg enriched coals as well as from mining activities, Hg contamination in this area has become a concern [43][44][45][46].

Study Area
The Xiaxi River is located in Wanshan County, Guizhou Province, China. This river has been influenced by the world-famous Wanshan Hg mining, which is one of the thirteen large and super large-scale Hg mines located within the Guizhou Province, China. Figure 9 shows the extent of the Xiaxi River, and its location in Guizhou Province. The Hg tailings are distributed throughout the headwater reaches of the Xiaxi River. Due to the high geochemical background of Hg in this area and high Hg emissions from coal combustion using Hg enriched coals as well as from mining activities, Hg contamination in this area has become a concern [43][44][45][46].

Study Area
The Xiaxi River is located in Wanshan County, Guizhou Province, China. This river has been influenced by the world-famous Wanshan Hg mining, which is one of the thirteen large and super large-scale Hg mines located within the Guizhou Province, China. Figure 9 shows the extent of the Xiaxi River, and its location in Guizhou Province. The Hg tailings are distributed throughout the headwater reaches of the Xiaxi River. Due to the high geochemical background of Hg in this area and high Hg emissions from coal combustion using Hg enriched coals as well as from mining activities, Hg contamination in this area has become a concern [43][44][45][46].  Water samples were respectively collected at the eight sampling sites in the watercourse (Figure 9) during moderate and high flow periods in 2007 and 2008 [47,48]. Concentrations of THg, dissolved mercury (DHg), MeHg and total suspended solids (TSS) were measured in the laboratory for all collected samples. The collection, storage, and preservation techniques of samples complied with the USEPA Method 1631. Figure 10 shows the distribution of 2008's observed THg concentrations in the watercourse of the Xiaxi River system. The THg concentrations were much higher along the Xiaxi River and then decreased rapidly downstream from the source. Based on observed data collected in 2007 and 2008, the THg concentrations ranged in the surface water from 4.5 to 2100 ng·L −1 [47,49]. The THg concentrations in the benthic sediments ranged from 1.1 to 360 mg·kg −1 [49].
Water samples were respectively collected at the eight sampling sites in the watercourse (Figure 9) during moderate and high flow periods in 2007 and 2008 [47,48]. Concentrations of THg, dissolved mercury (DHg), MeHg and total suspended solids (TSS) were measured in the laboratory for all collected samples. The collection, storage, and preservation techniques of samples complied with the USEPA Method 1631. Figure 10 shows the distribution of 2008's observed THg concentrations in the watercourse of the Xiaxi River system. The THg concentrations were much higher along the Xiaxi River and then decreased rapidly downstream from the source. Based on observed data collected in 2007 and 2008, the THg concentrations ranged in the surface water from 4.5 to 2100 ng·L −1 [47,49]. The THg concentrations in the benthic sediments ranged from 1.1 to 360 mg·kg −1 [49].

Model Development and Calibration
The W2 model was applied to the mainstem Xiaxi River from the headwater until the outlet of the watershed (Figure 9). The mainstem river is about 17 km long. The Xiaxi River is across the whole watershed and has been the main concern of Hg contamination. Table 4 listed the geometry of the Xiaxi River reaches and major tributaries. The Xiaxi River W2 model mainly consists of reaches R320, R410, R180 and R20 (Figure 9). River reaches R450 and R500 were treated as Hg contaminated sources, and R50 and R440 were treated as uncontaminated sources. The plan view of the Xiaxi River W2 model grid is shown in Figure 11. The Xiaxi River was discretized with fifty-three horizontal segments with Δx ranging from 338.8 to 468.7 m. The averaged water depth was 0.5 m. The vertical layer thickness was set up with 0.25 m.

Model Development and Calibration
The W2 model was applied to the mainstem Xiaxi River from the headwater until the outlet of the watershed (Figure 9). The mainstem river is about 17 km long. The Xiaxi River is across the whole watershed and has been the main concern of Hg contamination. Table 4 listed the geometry of the Xiaxi River reaches and major tributaries. The Xiaxi River W2 model mainly consists of reaches R320, R410, R180 and R20 (Figure 9). River reaches R450 and R500 were treated as Hg contaminated sources, and R50 and R440 were treated as uncontaminated sources. The plan view of the Xiaxi River W2 model grid is shown in Figure 11. The Xiaxi River was discretized with fifty-three horizontal segments with ∆x ranging from 338.8 to 468.7 m. The averaged water depth was 0.5 m. The vertical layer thickness was set up with 0.25 m.  The inputs used for the Xiaxi River W2 model include external loads, initial and boundary conditions, and water quality parameters. The input data were assembled based on field data collected during the project period (2007)(2008), outputs generated from other models, and literatures. Because Xiaxi River lacks routine hydrological monitoring, the HEC-HMS model was used to compute the stream flow [18] and fed into the W2 model. Xiaxi River inflow discharges for reaches R450, R500, R50 and R440 are shown in Figure 12a The inputs used for the Xiaxi River W2 model include external loads, initial and boundary conditions, and water quality parameters. The input data were assembled based on field data collected during the project period (2007)(2008), outputs generated from other models, and literatures. Because Xiaxi River lacks routine hydrological monitoring, the HEC-HMS model was used to compute the stream flow [18] and fed into the W2 model. Xiaxi River inflow discharges for reaches R450, R500, R50 and R440 are shown in Figure 12a  The inputs used for the Xiaxi River W2 model include external loads, initial and boundary conditions, and water quality parameters. The input data were assembled based on field data collected during the project period (2007)(2008), outputs generated from other models, and literatures. Because Xiaxi River lacks routine hydrological monitoring, the HEC-HMS model was used to compute the stream flow [18] and fed into the W2 model. Xiaxi River inflow discharges for reaches R450, R500, R50 and R440 are shown in Figure 12a Two suspended solids (fine and coarse solids) were simulated to characterize the partitioning and distribution of HgII and MeHg along the river. The fine solids represent particles with corresponding particle diameter of 0.01 mm, and the coarse solids represent particles with corresponding particle diameter of 0.1 mm. Figure 13 shows the solids grading curves of four sampling sites in the river. Fine solids account for 88.5% of the TSS. The boundary conditions of TSS and THg were calculated based on the following linear regression relationships developed by Lin et al. [ where Q is the stream flow (m 3 ·s −1 ), TSS is total suspended solids (mg·L −1 ), THg is the total concentration (ng·L −1 ). Two suspended solids (fine and coarse solids) were simulated to characterize the partitioning and distribution of HgII and MeHg along the river. The fine solids represent particles with corresponding particle diameter of 0.01 mm, and the coarse solids represent particles with corresponding particle diameter of 0.1 mm. Figure 13 shows the solids grading curves of four sampling sites in the river. Fine solids account for 88.5% of the TSS. The boundary conditions of TSS and THg were calculated based on the following linear regression relationships developed by Lin et al. [18]. [THg] R450 = 1435Q + 71 (21c) [THg] R500 = 2496Q + 592 (21d) where Q is the stream flow (m 3 ·s −1 ), TSS is total suspended solids (mg·L −1 ), THg is the total concentration (ng·L −1 ). Two suspended solids (fine and coarse solids) were simulated to characterize the partitioning and distribution of HgII and MeHg along the river. The fine solids represent particles with corresponding particle diameter of 0.01 mm, and the coarse solids represent particles with corresponding particle diameter of 0.1 mm. Figure 13 shows the solids grading curves of four sampling sites in the river. Fine solids account for 88.5% of the TSS. The boundary conditions of TSS and THg were calculated based on the following linear regression relationships developed by Lin et al. [ where Q is the stream flow (m 3 ·s −1 ), TSS is total suspended solids (mg·L −1 ), THg is the total concentration (ng·L −1 ). Based on available observed data, the Xiaxi River W2 model was calibrated in several steps. First, the manning's roughness coefficient 'n' was calibrated so that simulated discharge approximated observed flow discharge [18]. The final manning's roughness coefficient 'n' was set as 0.035. Next, the settling velocities of fine and coarse particles were calibrated, and simulated TSS approximated observed TSS concentrations along the river. In the Hg simulation module, the partitioning coefficients and transformation rates for HgII and MeHg were calibrated through a sensitivity analysis. Finally, the W2 model performance was evaluated by root mean square error (RMSE) and relative error (RE).
where C M is the modeled concentration (ng·L −1 for Hg and mg·L −1 for TSS), C O is the observed concentration (ng·L −1 for Hg and mg·L −1 for TSS). Figure 14 presents the comparison of W2 simulated flow discharge and simulated data [18] at the mouth of the Xiaxi River. The average error is 0.12 m 3 /s. High flood events mainly took place during the summer time (June to August), which were consistent with observed rainfall events. The W2 model was able to reproduce observed hydrodynamic results in the river.

Results and Discussion
Based on available observed data, the Xiaxi River W2 model was calibrated in several steps. First, the manning's roughness coefficient 'n' was calibrated so that simulated discharge approximated observed flow discharge [18]. The final manning's roughness coefficient 'n' was set as 0.035. Next, the settling velocities of fine and coarse particles were calibrated, and simulated TSS approximated observed TSS concentrations along the river. In the Hg simulation module, the partitioning coefficients and transformation rates for HgII and MeHg were calibrated through a sensitivity analysis. Finally, the W2 model performance was evaluated by root mean square error (RMSE) and relative error (RE).
where CM is the modeled concentration (ng·L −1 for Hg and mg·L −1 for TSS), CO is the observed concentration (ng·L −1 for Hg and mg·L −1 for TSS). Figure 14 presents the comparison of W2 simulated flow discharge and simulated data [18] at the mouth of the Xiaxi River. The average error is 0.12 m 3 /s. High flood events mainly took place during the summer time (June to August), which were consistent with observed rainfall events. The W2 model was able to reproduce observed hydrodynamic results in the river. After the hydrodynamics was calibrated, the W2 model was further calibrated for suspended solids. The modeled TSS concentrations were compared with observed data in Table 5. The results show that the RMSE for TSS is 0.10 mg·L −1 in 4 September 2007, and the relative error is 10%. For the TSS in 8 August 2008, the RMSE is 0.081 mg L −1 , and the relative error is 2.62%. The RMSE and RE values in high flow period are much lower than that in moderate flow period. The RE value increased when TSS concentration decreased, especially when TSS concentration is lower than 1.0 mg·L −1 . The settling velocities of fine and coarse particles with 5 m·d −1 and 25 m d −1 were used in the model. After the hydrodynamics was calibrated, the W2 model was further calibrated for suspended solids. The modeled TSS concentrations were compared with observed data in Table 5. The results show that the RMSE for TSS is 0.10 mg·L −1 in 4 September 2007, and the relative error is 10%. For the TSS in 8 August 2008, the RMSE is 0.081 mg L −1 , and the relative error is 2.62%. The RMSE and RE values in high flow period are much lower than that in moderate flow period. The RE value increased when TSS concentration decreased, especially when TSS concentration is lower than 1.0 mg·L −1 . The settling velocities of fine and coarse particles with 5 m·d −1 and 25 m d −1 were used in the model. The partitioning coefficient of Hg on inorganic solids has a wide range in literatures. Allison and Allison [41] provided a range of 10 2 -10 6 L·kg −1 . From a study of several Wisconsin lakes conducted by Rice et al. [50], the partitioning coefficient of Hg ranged from 8 × 10 4 to 3.4 × 10 5 . Suggested values for both methylation and demethylation rates are 10 −5 -10 −2 d −1 [41]. Methylation rate ranged between 0.01 to 0.09 d −1 in a clear lake in Wisconsin [51]. Studies by Carroll et al. [52] and Oremland et al. [53] indicated that methylation rates in Carson River ranged between 0.0025 and 0.011 d −1 and demethylation rates ranged between 0.12 and 0.55 d −1 . A sensitivity analysis was performed using a range of partitioning coefficient (K d ) and methylation/demethylation rate (M/D). W2 modeled DHg and MeHg concentrations with a range of values of K d and M/D are presented in Tables 6 and 7. Based on the sensitivity analysis, the partitioning coefficients of Hg on fine and coarse solids were set to 10 6 and 10 4 L·kg −1 . The methylation and demethylation rates were set to 0.01 and 0.1 d −1 in this study.

Summary and Conclusions
The Hg simulation module has been developed and incorporated into the 2D hydrodynamic and water quality model CE-QUAL-W2. The enhanced W2 model was tested and validated against WASP mercury model. The model results indicated that the enhanced W2 model was well designed, and the algorithms and formulations were correctly implemented. Application of the enhanced W2 model to the Xiaxi River, a river course in Wanshan Hg mining area of Guizhou Province, China, further evaluated the capability of the enhanced W2 model to simulate the transport and cycling of Hg in water bodies. The enhanced W2 model was able to reproduce observed TSS, THg, DHg and MeHg. THg, DHg and MeHg concentrations decreased rapidly downstream from the source and stabilized at lower levels, which is likely due to high association of Hg to solids and the rapid settlement of these solids. The source control of Hg discharge is the key to reduce the Hg contamination along the Xiaxi River. This application demonstrated the W2 model capability in predicting complex transport and cycling of Hg species in water bodies. When more observed data become available, the robustness of the Xiaxi River W2 model will be further calibrated and

Summary and Conclusions
The Hg simulation module has been developed and incorporated into the 2D hydrodynamic and water quality model CE-QUAL-W2. The enhanced W2 model was tested and validated against WASP mercury model. The model results indicated that the enhanced W2 model was well designed, and the algorithms and formulations were correctly implemented. Application of the enhanced W2 model to the Xiaxi River, a river course in Wanshan Hg mining area of Guizhou Province, China, further evaluated the capability of the enhanced W2 model to simulate the transport and cycling of Hg in water bodies. The enhanced W2 model was able to reproduce observed TSS, THg, DHg and MeHg. THg, DHg and MeHg concentrations decreased rapidly downstream from the source and stabilized at lower levels, which is likely due to high association of Hg to solids and the rapid settlement of these solids. The source control of Hg discharge is the key to reduce the Hg contamination along the Xiaxi River. This application demonstrated the W2 model capability in predicting complex transport and cycling of Hg species in water bodies. When more observed data become available, the robustness of the Xiaxi River W2 model will be further calibrated and validated. The current model can also be used to perform scenarios analysis and quantify the best management practice for the XiaXi River.