Modeling hydro-biogeochemical transformation of chromium in hyporheic zone: Effects of spatial and temporal resolutions

Journal of Hydrology 579 (2019) 124152

Contents lists available at ScienceDirect

Journal of Hydrology journal homepage: www.elsevier.com/locate/jhydrol

Research papers

Modeling hydro-biogeochemical transformation of chromium in hyporheic zone: Effects of spatial and temporal resolutions

T

Chen Yanga,b, Fei Zhengc, Yuanyuan Liud, You-Kuan Zhanga,b, Wei Liue, Qiang Zhanga,b, ⁎ Xiaofan Yangf, a

Guangdong Provincial Key Laboratory of Soil and Groundwater Pollution Control, School of Environmental Science and Engineering, Southern University of Science and Technology, Shenzhen 518055, China b State Environmental Protection Key Laboratory of Integrated Surface Water-Groundwater Pollution Control, School of Environmental Science and Engineering, Southern University of Science and Technology, Shenzhen 518055, China c MCC Water Environment Technology Research Institute, MCC Huatian Engineering & Technology Corporation, Nanjing 210019, China d School of Earth Sciences and Engineering, Nanjing University, Nanjing 210023, China e Institute of Geological Survey, China University of Geosciences, Wuhan 430074, China f State Key Laboratory of Earth Surface Processes and Resource Ecology, Faculty of Geographical Science, Beijing Normal University, Beijing 100875, China

A R T I C LE I N FO

A B S T R A C T

This manuscript was handled by Huaming Guo, Editor-in-Chief

Effects of spatial and temporal resolutions (SR and TR) on modeling hydro-biogeochemical transformation of chromium (Cr) are important in simulating reactive transport processes. The current study was conducted in the hyporheic zone (HZ) at the Hanford Site of the U.S. Department of Energy, which has been known for its highly heterogeneous sediments and transient hydrodynamics. Distributions of hydraulic conductivity and sedimentassociated Fe concentration were averaged at a group of SRs, while measured hourly water levels were further moving averaged at daily and monthly TRs. Fe concentration is selected for assembling geochemical heterogeneity due to its important role in redox transformation of Cr at the site. Three Fe distributions, with Fe concentrated on small, medium, and large sediment grains, respectively, were also considered. A series of flow and reactive transport simulations configured with different combinations of SRs, TRs, and Fe distributions were conducted. Simulated results revealed that Cr(VI) discharged to river is underestimated if SR decreases. For both the hydrodynamics and the discharge of Cr(VI) to river, difference caused by SR can be amplified with the decrease of TR. Biogeochemical transformation of Cr is more dependent on SR while hydrodynamics is on TR. The stronger control of SR than TR on biogeochemical transformation of Cr is resulted from more sensitive increase of Fe(II) with decreasing SR than with decreasing TR. Effect of SR is highly sensitive to variations of SR with Fe on medium-size grains while is persistent to a much smaller SR with Fe on small-size grains. Results from the current study are expected to benefit modelers on the selections of the spatial-temporal resolutions for inhouse modeling and field sampling, which may also have implications for upscaling.

Keywords: Discretization scale Chromium Reactive transport modeling Hydro-biogeochemical processes Hyporheic zone

1. Introduction Physical and geochemical heterogeneities in the subsurface, such as spatial variations of the hydraulic conductivities, chemical species and concentrations, significantly regulate the reactive transport of contaminants (Li et al., 2006; Liu et al., 2014; Perujo et al., 2017; Salamon et al., 2007; Wen and Li, 2017). For instance, the reaction rates, which are largely dependent on properties of sediments (Liu et al., 2014; Li et al., 2010), in relatively homogeneous batch reactors can be several orders of magnitude larger than those in highly heterogeneous natural systems (Li et al., 2008; Salehikhoo et al., 2013; White and Brantley,

⁎

2003). Moreover, with the rapid development of hydraulic facilities (e.g., dams) and regulations, temporal variabilities of regional hydrodynamics become more complex, which further affect the migration and the transformation of contaminants in the riverine system (Shuai et al., 2019; Zachara et al., 2016). As such, the spatial and temporal resolutions (scales in discretization) used in the reactive transport modeling are critical because small-scale heterogeneities (i.e., smaller than the discretization scale) may have huge impact on macroscopic system behaviors, nevertheless, are unavoidably neglected (Li et al., 2006). Selected studies have focused on the effect of model resolution on

Corresponding author. E-mail address: [email protected] (X. Yang).

https://doi.org/10.1016/j.jhydrol.2019.124152 Received 17 July 2019; Received in revised form 12 September 2019; Accepted 14 September 2019 Available online 20 September 2019 0022-1694/ © 2019 Elsevier B.V. All rights reserved.

Journal of Hydrology 579 (2019) 124152

C. Yang, et al.

the simulation results. Vieux and Needham (1993) indicated that the sediment yield increased as much as 32% due to the shortened stream lengths caused by large cell size used in the Agriculture Nonpoint Source Pollution Model. Molnar and Julien (2000) tested a group of grid-cell sizes varying from 127 m to 914 m using the CASC2D rainfallrunoff model, and suggested that a proper application of a coarser grid under certain circumstances requires subsequent compensation, such as using more sophisticated descriptions of the parameters in the model. Yu et al. (2014) examined the results of a hydrological model using different spatial (15–200 m) and temporal (5–60 min) resolutions, and found that the spatio-temporal resolution significantly affected the hydrological response of the watershed because the resolutions were correlated to the elevation and precipitation which were important in the simulations. Stenta et al. (2017) pointed out that in order to obtain the same hydrological response when the grid-cell size was increased, other parameters in the hydrological model needed to be re-scaled. Unfortunately, these studies are limited in hydrological modeling without reactive transport. To our best knowledge, the combined effects of spatial and temporal resolutions on modeling biogeochemical transformation of the redox sensitive chemicals were rarely studied. In addition, due to the excessive computational and fieldwork costs in reactive transport modeling (Li et al., 2010, 2011), the conduction of such a study is also expected to help modelers to choose the most appropriate resolution in their simulations to balance the model precision and the computational cost. Last, but not the least, although previous studies have tried to establish upscaling laws of the reaction rates to bridge the gaps between scales (Wen and Li, 2017), either theoretical or numerical upscaling (coarse graining) is difficult to generalize for modeling in practice. An intuitive knowledge of the effects of spatiotemporal resolutions on reactive transport modeling may provide implications for such efforts on upscaling. In this study, we investigated the role of spatial and temporal resolutions (SR and TR) in reactive transport modeling based on a representative hyporheic zone (HZ) located in the 100-BC Area at the Hanford Site of the U.S. Department of Energy (US DOE). The subsurface environment of the Hanford 100-BC Area is composed of heterogeneous aquifers with transient hydrodynamics due to the poorlysorted sediments (Hammond et al., 2011) and the upstream dam operations (Yabusaki et al., 2008). Chromium (Cr) contamination is serious at the site due to the release of radioactive material (Williams et al., 2008), and the HZ is the last barrier determining the fate of Cr: discharged to river as Cr(VI) or retarded in the HZ as Cr(III) (Yang et al., 2018). A reaction network for biogeochemical transformation of Cr (Liu et al., 2017; Yang et al., 2018), specially established for Cr retardation in the HZ at the Hanford Site, was selected for reactive transport modeling. In the reaction network, Cr(VI) can be reduced to immobilized Cr(III) in the HZ by sediment-associated Fe(II) while the dissolved oxygen (DO) intrusion from the river inhibits the immobilization of Cr(VI). Distribution of sediment-associated Fe(II) is largely dependent on the heterogeneities of subsurface materials (Li et al., 2011) while DO intrusion is mainly controlled by hydrodynamics (Yang et al., 2018). In addition, previous study (Yang et al., 2018) on the effect of TR promoted the further exploration on SR and combined TR and SR. In the following sections, the study area is briefly introduced in Section 2. Details of the reactive transport modeling of Cr are shown in Section 3. Recent numerical models in Yang et al. (2018) were adopted in this study with additional configurations as follows. Spatial distributions of the physical and geochemical parameters were obtained with a stochastic approach based on measurements in previous studies (e.g., Williams et al., 2008). Then the stochastically generated or calculated spatial parameters and the measured water levels were successively upscaled and moving averaged from baseline scales to larger scales, respectively. A group of scenarios with different SRs and TRs were simulated using a community groundwater solver for water flow coupled with an in-house simulator for reactive transport. In Section 4,

Fig. 1. Plan view of the Hanford Site and location of the 2D transect.

effects of TRs and SRs on modeling Cr transport and transformation in the HZ are discussed based on the simulation results from Section 3. Finally, conclusions are summarized in Section 5. 2. Study area The study site is located at the Hanford 100-BC Area in the southeastern Washington State of the U.S. Detailed descriptions of this site were given in Yang et al. (2018), and here only a brief introduction is summarized. A two-dimensional (2D) transect is selected as the modeling domain (indicated by the yellow dash line in Fig. 1). Cr(VI) contamination is critical in the Hanford 100-BC Area, which has potential environmental risks due to the discharge of Cr(VI) to the nearby Columbia River and its surrounding ecosystems (Truex et al., 2015). Particularly, river water discharged from upstream dam operations induces extra diurnal and weekly fluctuations of the river stage (Yabusaki et al., 2008). Thus, the enhanced exchange between river water (RW) and groundwater (GW) in the HZ could cause additional environmental risks (Yang et al., 2018). The selected transect starts at Well 199-B2-14 and ends at the center of the river with a height of 45 m and a length of 314 m (Fig. 2). Sediments in the transect are categorized as the Hanford formation with hydraulic conductivities over 2000 m/d and Ringold E formation with typical values of 40–120 m/d (Ma et al., 2010). Specifically, a 0–3 m thick alluvium layer with low permeability (~4 m/d) along the river bank is found to be an active biogeochemical zone (hot spots) where sediments have the redox potential to reduce Cr (VI) to Cr(III) (Herzog et al., 2015; Liu et al., 2017; Moser et al., 2003; Williams et al., 2008). 3. Models setup 3.1. Stochastic simulation and analysis of the soil properties Stochastic simulation of heterogeneities is prerequisite in flow and reactive transport modeling in this study. For physical heterogeneity, hydraulic conductivity (K) was selected (Li et al., 2010; Salamon et al., 2007; Zhang, 2002). For geochemical heterogeneity, concentration of sediment-associated Fe (CFe) was selected due to its important role in the immobilization of Cr(VI) in the HZ (Liu et al., 2017; Yang et al., 2018). Distributions of both K and CFe are closely related to grain size 2

Journal of Hydrology 579 (2019) 124152

C. Yang, et al.

Fig. 2. Geometry and stratigraphy of the model domain. Dashed blue line represents the average water table, red box the contaminant plume, and yellow dots the observation points. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.) 2 ρ g n3 ⎤ ⎛ dg ⎞ Ks = ⎜⎛ w ⎟⎞ ⎡ ⎟ ⎜ 2 ⎢ ⎥ ⎝ μ ⎠ ⎣ (1 − n) ⎦ ⎝ 180 ⎠

distribution of the sediments. For example, sediments composed of larger size grains usually have higher permeability and less CFe than smaller size grains (Li et al., 2010). Hence, in this study, the grain size distribution was determined based on the stochastically generated gamma log data (gg, borehole geophysical logging data), and then the K and CFe were calculated based on grain size distribution using their corresponding relationships.

where ρw and μ are the density and viscosity of the water, respectively, g is the gravitational constant, and n is the porosity of the Ringold formation. 3.1.3. Fe concentration CFe was calculated based on the grain size, which is considered as an improvement from previous studies that the chemical content was determined by a negative correlation with K (Li et al., 2010, 2011). Three Fe distributions were considered in this study to cover a wide range of grain size, which are listed in Table 2: (1) Distribution 1: Fe was concentrated on small-size grains, which is commonly found in previous studies (Brook and Moore, 1988; Lin et al., 2003; Martincic et al., 1990; Okweye et al., 2016; Salomons and Förstner, 1984; Yao et al., 2015). Thus, a 100% measured concentration (MC) (Liu et al., 2017) was set for dg < 0.08 mm and 1% MC for dg > 0.08 mm. (2) Distribution 2: Fe was concentrated on medium-size grains (Shang et al., 2011), which is a common distribution of redox sensitive Uranium at the Hanford Site. Thus, 100% MC was set for 0.08 mm < dg < 2.36 mm, and 1% MC for dg < 0.08 mm and dg > 2.36 mm. (3) Distribution 3: Fe was concentrated on small- and large-size grains (i.e., two sides) since previous studies reported similar or even higher heavy metal concentrations on coarser grains (Brook and Moore, 1988; Lin et al., 2003; Singh et al., 1999; Tessier et al., 1982; Tsai et al., 2003). Thus, 100% MC was set for dg < 0.08 mm and dg > 2.36 mm, and 1% MC for 0.08 mm < dg < 2.36 mm. There was no Fe distributed on grains with size larger than 8.25 mm, since only a < 8 mm size fraction was sieve-collected in Liu et al. (2017).

3.1.1. Gamma log data Geo-statistical analysis has not been performed/documented in Hanford 100-BC Area, so the stochastic distribution of gg in the Hanford 300 Area was utilized because both sites have similar geological properties (Shuai et al., 2019; Williams et al., 2008). The variogram of gg is described by a spherical model with three nested structures as follows (Williams et al., 2008): 3

γ (h) = C0 +

h

h

3

∑ Ci ⎡⎢1.5 ⎛ a ⎞ − 0.5 ⎛ a ⎞ ⎤⎥ i=1

⎣

⎜

⎟

⎜

⎟

⎝

i⎠

⎝

i⎠

(1)

⎦

where γ is the variogram value, h is the lag separation distance, C0 is the nugget, and Ci and ai are the sills and actual ranges, respectively. Values of these parameters are listed in Table 1 (Williams et al., 2008). The gg field was generated in the modeling domain using the SGSIM program (Sequential Gaussian Simulation program) in the geostatistical software library (GSLIB) (Deutsch and Journel, 1998). Generated gg was further scaled to the mean and the standard deviation of the Ringold formation, which are 215.832 and of 33.418, respectively (Williams et al., 2008). 3.1.2. Grain size and saturated hydraulic conductivity A correlation between the grain size and the gg was established as Eq. (2) (Williams et al., 2008):

dg = p2 + (p1 − p2 ) × [1 + (p3 × gg ) p4]−p5

3.2. Spatial and temporal resolutions

(2) Gamma log data and grain size were generated and calculated, respectively, on a grid with 5 cm × 5 cm resolution. Then the Ks and CFe calculated on 5 cm × 5 cm grid-cells were up-scaled to grid-cells of 10 cm × 10 cm, 50 cm × 50 cm, and 200 cm (horizontal) × 100 cm (vertical) respectively. They are referred to as 10, 50, and 200 cm gridcells in the rest of the paper, representing SRs from the laboratory to the field. Ks were up-scaled using the geometric mean (Wen and GomezHernandez, 1996), while CFe were up-scaled using the volume averaging (Williams et al., 2008):

where dg (mm) is the geometric mean diameter as the representative grain size metrics, and p1–p5 are constants as 28, 0.00001, 0.006, 20.9152, and 0.9522. Thus, dg was calculated based on the gg generated in section 3.1.1. The 20%, 60%, and 77% quantiles of the cumulative distribution function for the calculated diameters are 0.08 mm, 2.3566 mm and 8.2494 mm, respectively, which are comparable to the results of grain size analysis at the Hanford Site in Liu et al. (2017). Saturated hydraulic conductivity (Ks) of the sediment was estimated using the Kozeny-Carmen equation (Bear, 1972; Williams et al., 2008) as follows:

Table 2 Three CFe distribution patterns with grain size.

Table 1 Parameters for variogram of the gamma log data. Direction / / Horizontal Vertical

C0 0 / / /

C1 0.35 a1 30 4

(3)

Patterns C2 0.3 a2 200 14

C3 0.35 a3 1400 25

Small Medium Small and large

3

Diameters < 0.08 mm

0.08–2.36 mm

2.36–8.25 mm

> 8.25 mm

100% MC 1% MC 100% MC

1% MC 100% MC 1% MC

1% MC 1% MC 100% MC

0 0 0

Journal of Hydrology 579 (2019) 124152

C. Yang, et al.

temporal scales (Fig. 5). These TRs represent scales from the level dominated by human activities to the level controlled mainly by natural regulations. Similarly, a hypothetical time scale was used, meaning that the hourly time step was actually used in all scenarios. Hereinafter, increasing TR means the increase of temporal resolution, i.e., the increase of temporal variability, and vice versa. Therefore, a total of 9 velocity fields (10 cm-hourly, 10 cm-daily, 10 cm-monthly, 50 cmhourly, 50 cm-daily, 50 cm-monthly, 200 cm-hourly, 200 cm-daily, and 200 cm-monthly) and 27 reactive transport fields (each velocity field was used to calculate reactive transport with three CFe distributions) were simulated in this study. 3.3. Flow and reactive transport modeling 3.3.1. Modeling flow A massively parallel subsurface flow and reactive transport model (PFLOTRAN) (Hammond et al., 2014) was used to simulate variably saturated flow in the modeling domain. The water levels in Well 199B2-14 and at the river gauge were used as the specified heads on the left side and on the right slope under the river stage, respectively. The river bank above the river stage was treated as the seepage boundary. The top and bottom boundaries were no flow due to the low precipitation and the impermeable underlying sediments, respectively. The center of the river at the right boundary was assumed to be no flow as well. The characteristic curve of the soil moisture was based on the Brooks-Corey function, and the related parameters are listed in Table S1. In each scenario, steady-state flow field based on annual means of the well level and the river stage was first calculated as the initial condition for the corresponding transient flow field. After spin-up of the transient flow fields to reach quasi-steady state, one more year was simulated, the results of which were repeatedly used for 10 years of simulation in the reactive transport.

Fig. 3. Distributions of Ks at (a) 10 cm, (b) 50 cm, and (c) 200 cm SRs. n

CFeC =

∑ j = 1 CFeF , j VF , j n

∑ j = 1 VF , j

(4)

where V is the volume of a fine grid-cell, the subscripts C, F and j refer to the coarse grid-cell (10, 50, or 200 cm), the fine grid-cell (5 cm), and the index of the fine grid-cell, respectively. The summations were taken over all the fine grid-cells (n) belonging to the same coarse grid-cell. Distributions of Ks and CFe at different SRs are shown in Figs. 3 and 4, and S1 and S2 in the Supporting Information. Ks shown in Fig. 3 are those adjusted in the flow simulations. It should be noted that a 50 cm grid-cell was further divided into a group of 10 cm sub-grid cells, and the parameters applied on these sub-grid cells were kept the same to those on that 50 cm grid-cell (Fig. S3). Such treatment was also done for the grid of 200 cm resolution. Therefore, all simulations were actually conducted on the same grid of 10 cm resolution. As such, difference induced by SR could be retained while additional numerical deviations could be avoided. For all three SRs, the model domain was divided into 207,480 grid cells with 62,218 inactive cells to accommodate the irregular river boundary. Hereinafter, increasing SR means the increase of spatial resolution, i.e., the increase of parameter resolution, and vice versa. Time series of hourly water levels in Well 199-B2-14 and at the river gauge were moving averaged at daily (24 h) and monthly (720 h)

3.3.2. Reaction network of Cr in the HZ The biogeochemical reaction network in this study was first established by Liu et al. (2017) and further improved by Yang et al. (2018). Details (Table S2) of the reaction network included: (1) multi-rate Cr (VI) sorption and reduction by sediment-associated Fe(II) (not in aqueous phase) with the corresponding Fe(II) oxidation to Fe(III), (2) microbially-mediated regeneration of Fe(II) with dissolved organic carbon (DOC) as the electron donor and carbon source, (3) abiotic oxidation of Fe(II) by DO that competes with Cr(VI), (4) DO consumption by microorganisms as the electron acceptor with DOC as the electron donor and carbon source, and (5) DOC production from the particulate organic carbon (POC) transformation. The reactive Fe(II) for Cr reduction was assumed to be only present in the alluvium, so that in other regions Cr was only subject to sorption (shown in the blue and purple areas in Fig. 2). Specific parameters of the reactions are provided in Table S3. 3.3.3. Modeling reactive transport Reactive transport simulations of 10 years for Cr(VI) with other relevant chemical species (i.e., organic carbon, Fe(II)/Fe(III), and DO) in the alluvium and Cr(VI) sorption to sediments outside the alluvium were performed based on the flow fields obtained a priori. Initial conditions based on measurements (Liu et al., 2017) for these reactions are provided in Table S4. Cr(VI) with a concentration of 2000 μg/L (Dresel et al., 2008; Johnson, 2016; Liu et al., 2017; Smoot et al., 2011; Truex et al., 2015; Yang et al., 2018) was injected as a plug at the fifth year, shown as the 10 m × 5 m red box in Fig. 2. For all the mobile species, the top, bottom, and right sides of the model domain were assumed to be no-flux boundaries, and the left side was treated as constant concentration when the water flow was into the model domain and as zero concentration gradient when the water flow was out of the domain. Species concentrations in the RW were assumed to be constant. The reactive transport was solved by an in-house developed simulator as an alternative of the reaction sandbox provided in PFLOTRAN, which has

Fig. 4. Distributions of CFe at (a) 10 cm, (b) 50 cm, and (c) 200 cm SRs. Fe is on small and large size grains. 4

Journal of Hydrology 579 (2019) 124152

C. Yang, et al.

River

Elevation (m)

124

Well

122 120 118 116 0

(a) Hourly 4380 t (hour)

(b) Daily 8760

0

(c) Monthly 4380 t (hour)

8760

0

4380 t (hour)

8760

Fig. 5. Groundwater levels in Well 199-B2-14 and stages of river at (a) hourly, (b) daily, and (c) monthly TRs. Grey lines in each subplot are of hourly TR used for comparison.

In the first column of Fig. 7, Cr(VI) mass to river is minimum at 200 cm SR for all three Fe distributions. With a lower SR, Fe(II) is more averaged to larger size grains, where there is no or less Fe(II) at a higher SR. Thus, the retardation capability of the HZ to Cr(VI) becomes stronger, and Cr(VI) is more prevented to river. This overestimation of retardation capability in the HZ means the underestimation of environmental risk about the Cr(VI)-discharge to river. This phenomenon was drawn in Fig. 8 by taking a scenario with Fe distributed on medium size grains. Cr(VI) enter the river through corridors with a small amount of Fe(II) at 10 cm SR, while Cr(VI) are retarded in the HZ at 200 cm SR. Similar phenomena have been reported in previous studies. Elfeki et al. (2012) showed that higher heterogeneity, which corresponds to higher SR here, might lead to larger variance of the contaminant plume (i.e., more dispersion of the plume). Zheng and Jiao (1998) described the small patches connected to the main body of the plume by a narrow neck, and the fringy feature on the iso-surface of the plume at the MADE site, where the sediments are extremely heterogenous. These phenomena are mainly caused by the preferential flow paths due to high heterogeneity (Elfeki et al., 2012; Zheng and Jiao, 1998). In summary, with higher SR, Cr(VI) is more easily discharged to river along paths with larger grain sizes, where K is larger and CFe is lower. When Fe is on small size grains, Cr(VI) enters the river at a maximum value and Cr(III) precipitated in the HZ is the minimum, regardless of the SR. This can be observed by that the grey line-group is always the highest in the first column while the lowest in the second column in Fig. 7. This indicates that the highest risk to river occurs when Fe is on small size grains. When Fe is on small and large size grains, higher Cr(III) mass is observed (Fig. 7j), meaning more retarded Cr(VI) in the HZ, since large-size grains are always the preferential flow paths for Cr(VI). This is consistent with Li et al. (2011), in which the largest amount of U(IV) was generated at the edge of the preferential flow paths interfacing with Fe(III)-rich areas. Dissolved and adsorbed Cr(VI) in the inland area is independent of Fe distribution and is of slight difference at different SRs (the third and fourth columns in Fig. 7). This difference is due to the hydrodynamic discrepancies induced by different physical heterogeneities at different SRs. Moderate hydrodynamics at a lower SR result in more dissolved and adsorbed Cr (VI) in the inland area. This is shown by the fact that red lines are always higher than blue and black lines with a maximum discrepancy of ~1 g. This also indicates that the effect of geochemical heterogeneity is more important than that of the physical heterogeneity, which was rarely discussed in previous studies. For two sides Fe distribution, Cr(VI) mass to river is not sensitive to SR (Fig. 7i). For medium size distribution, it is the most sensitive case to SR due to the large discrepancy among three SRs. Hence, with the same decrease of SR, retardation capability of the HZ is the most overestimated when Fe is on the medium size grains, and thus the environmental risk to river is the most underestimated. With medium size distribution, effect of upscaling lasts to about 200 cm SR, since Cr(VI)discharge to river is totally prevented at this SR (Fig. 7e). This means

been proved to be feasible for the current study site (Yang et al., 2018). Then the Cr(VI) mass dissolved in the inland area, adsorbed to the sediments, and discharged to river, the Cr(III) mass precipitated in the HZ, and the DO mass transported into the model domain and dissolved in the inland area (Yang et al., 2018) were calculated for the following analysis. 4. Results and discussion 4.1. Physical and chemical distributions at different SRs and TRs In Fig. 3, with decreasing SR, preferential flow paths are gradually reduced, which can be observed by the shrinkage of the red area. The reason is that large Ks are removed to compensate for extremely low values during upscaling. At places where the color becomes warmer, hydraulic connections are strengthened due to the replacement of small Ks with larger values. Flow fields in the domain would become smoother with fewer sudden changes of velocity, which will be discussed in Section 4.2. Similar patterns are observed for CFe in Fig. 4. CFe are more evenly distributed at lower SR, since some parts with lower CFe are supplemented in up-scaling, which could sacrifice the strong redox ability of the nearby places. In Fig. 5, with the moving average from hourly to monthly TR, crossovers between GW level and river stage in the middle of a year are gradually filtered out. This can cause decrease of the flow reversal frequency in the modeling area (Yang et al., 2018). Flow reversal here means RW intrusion to subsurface instead of the general GW discharge to river. It has importance on redox transformation of Cr in the HZ (Yang et al., 2018). This effect might be magnified when superimposed on the spatial scaling effect in this study. 4.2. Hydrodynamics and Cr transformation at different SRs Fig. 6 shows horizontal velocities (U) at three observation points (P1–P3) located from the well side (Well 199-B2-14) to riverside. First of all, simulated horizontal velocities are consistent with results from previous studies conducted at the Hanford Site, e.g., −25 to 25 m/d in Yabusaki et al. (2008), indicating the accuracy of the model in this study. In this section, the first and second rows are focused since here we only discuss the effect of SR. At P1, the water is almost stagnant at 10 cm SR (~0 m/d), and becomes dynamic with decreasing SR (Fig. 6a and d), e.g., over 0.2 m/d at 200 cm SR. In contrast, U at P2, where it is extremely high, decreases to about one third (~15 m/d to ~5 m/d) at peak values from 10 cm to 200 cm SR, leading to a more moderate state (Fig. 6b and e). This indicates that upscaling tends to result in a more uniform flow field in space. Hence, variations of hydrodynamics with SR discussed here are consistent with the prior judgement in Section 4.1. What’s more, the effect of spatial upscaling is not obvious in the middle of the year when interactions between GW and RW are especially intense, e.g., little difference in three SRs at about 4000 h in Fig. 6a. 5

Journal of Hydrology 579 (2019) 124152

C. Yang, et al.

Fig. 6. Simulated horizontal velocities (U) at three selected observation points, P1 (531.80 m, 120.45 m), P2 (571.8 m, 120.45 m) and P3 (679.95 m, 120.45 m). (a)–(f), (g)–(l), and (m)–(r) represent hourly, daily, and monthly TRs, respectively. The first, third, and fifth rows represent the whole year, while the second, fourth, and sixth rows represent 6500–7000 h interval.

hydrodynamics of different SRs become more obvious with decreasing TR. In general, for hydrodynamics, either the discrepancy induced by TR is much more obvious than that induced by SR, or hydrodynamics are more controlled by TR. Thus, at places with intense anthropogenic activities inducing extra water flow fluctuations at small time-scales, it is necessary to use higher TR to capture such hydrodynamics while the discrepancy induced by SRs is insignificant under such transient conditions. With lower TR, more Cr(III) precipitated in the HZ (the second column in Fig. 9) regardless of the SR. This means an overestimated retardation capability of the HZ. Therefore, the risks seem to be further underestimated when compared to that caused by a lower SR. For example, from 10 cm-hourly to 200 cm-hourly then to 200 cm-monthly (Fig. 9j), Cr(III) precipitation continuously increases. However, Cr(VI) to river (the first column in Fig. 9) also increase with decreasing TR, which in fact neutralize a small portion of the underestimated risks from the individual decrease of SR, e.g., from 10 cm-hourly to 50 cmhourly then to 50 cm-monthly (Fig. 9i). The higher the SR, the more the neutralization at lower TR. Taking Fig. 9i for example, difference between black line (10 cm SR) and corresponding grey line is larger than that for red line (200 cm SR). Then, increase of both Cr(VI) discharge and Cr(III) precipitation are balanced by decrease of dissolved and adsorbed Cr(VI) in the inland area at lower TR. This can be observed by

that a maximum underestimation of environmental risk to river has been achieved at 200 cm SR, and the underestimation would not increase any more with the continued decrease of SR. For small size distribution (Fig. 7a), Cr(VI) mass to river is less sensitive to SR than that with medium size distribution. However, upscaling has the most persistent effect on Cr(VI) mass to river since there is still a large amount of Cr(VI) discharged to river at 200 cm SR. Effect of upscaling should disappear at a SR lower than 200 cm at which Cr-discharge to river becomes almost zero as that in Fig. 7e. As such, the selection of SR deserves particular attention when environmental risk assessment is conducted at places where redox sensitive chemicals are distributed on small size grains due to the persistent effect of upscaling. Attention should also be paid when Fe is distributed on medium size grains, which might be the dominating distribution pattern of Fe at the Hanford Site.

4.3. Hydrodynamics and Cr transformation at different combinations of SR and TR With decreasing TR, fluctuations of U at higher TR are smoothed out and hydrodynamics become gentle. More importantly, both frequency and magnitude of the reversed flow decrease with decreasing TR, which can be observed, for instance, by comparing Fig. 6f with l. Difference in 6

Journal of Hydrology 579 (2019) 124152

C. Yang, et al.

Fig. 7. Calculated masses of Cr(VI) to river (first column), Cr(III) precipitated in the HZ (second column), and Cr(VI) (dissolved + adsorbed) in the inland domain (third and fourth columns). The first to third rows represent Fe distributed on small, medium, and small and large grains, respectively. They were based on hourly water levels. Dashed grey lines in the second and third rows represent corresponding lines in the first row, and are shown for easy comparison.

of organic carbon (OC) is generally lower at 200 cm SR (Fig. 10b and f) because it is more consumed by biological activities at the corresponding higher DO level. In contrast, Fe(II) is higher at a lower SR, triple concentration at 200 cm SR to that at 50 cm SR (Fig. 10h), and that is why more Cr(III) is generated at 200 cm SR (Fig. 9b, f and j). Therefore, more OC, particularly the POC, is consumed to maintain this relatively high Fe(II) level. It is clearly shown in Fig. 10g that POC concentration is about 10 μmol/g lower at 200 cm SR than that at 10 cm SR. In Fig. 10a, flow reversals at daily and monthly TRs carry less DO from RW to GW, which can be observed due to smaller peaks (~80 μmol/L) relative to those at hourly TR (~300 μmol/L). Thus, less OC is consumed, as shown in the second and third columns in Fig. 10. Therefore, more Fe(II) is left and more Cr(III) is precipitated in the HZ at lower TRs, which is the explanation for the results in Fig. 9. However, Fe(II) content increases moderately with the decrease of TR when compared to that with the decrease of SR. In Fig. 10h, concentration of Fe(II) increase from 0 to 0.01 then to ~0.04 μmol/g with decreasing SR while from 0.01 to ~0.02 μmol/g from hourly to daily and monthly. Therefore, the redox transformation of Fe(II) and thus Cr is more sensitive to SR than TR. Besides, concentrations of OC (Fig. 10b, f, c and g) decrease with time as the gradual intrusion of DO from RW to GW. When a large RW intrusion occurs in the middle and at the end of the year, DO has obvious concentration peaks (Fig. 10e) while OC has clear valleys (Fig. 10f and g). Fe(II) also gradually decreases with time due to the DO intrusion. The above findings are also consistent with those in Yang et al. (2018).

Fig. 8. Simulated concentration distribution of Cr(VI) at 0.027 and 0.082 years for 10 (a and b), 50 (c and d), and 200 cm (e and f) SRs. Fe is on medium size grains, and hourly water levels were used.

comparing grey line-group with other line-groups in the third and fourth columns of Fig. 9. Decrease of dissolved and adsorbed Cr(VI) are observed at lower TR, which is consistent with Yang et al. (2018). In addition, difference of Cr(VI) mass to river among different SRs are increased with decreasing TR (Fig. 9i). This can be observed by comparing the grey line-group with other line-groups in the first column of Fig. 9. In contrast to hydrodynamics, redox transformation of Cr is more dependent on the SR than the TR (e.g., 10 cm-hourly/200 cm-hourly vs. 10 cm-hourly/10 cm-monthly).

5. Conclusions

4.4. Mechanisms for effects of SR and combined SR and TR on Cr transformation

In this study, taking the redox sensitive Cr as an example, numerical simulations were conducted to illustrate effect of spatio-temporal resolution on modeling biogeochemical transformation of contaminants in a hyporheic zone (HZ) at the US DOE Hanford Site. The SRs are 10, 50, and 200 cm, representing scales from the laboratory to the field, while TRs are hourly, daily, and monthly, representing scales from the

In Fig. 10a and 10e, DO concentration has insignificant dependence on SR at P4. This is shown as the slight increase of DO with decreasing SR. However, the difference of DO among different SRs is much more obvious at the regional scale (i.e., the whole modeling domain) in Fig. 11, and its content is the highest at 200 cm SR. The concentration 7

Journal of Hydrology 579 (2019) 124152

C. Yang, et al.

Fig. 9. Calculated masses of Cr(VI) to river (first column), Cr(III) precipitated in the HZ (second column), and Cr(VI) (dissolved + adsorbed) in the inland domain (third and fourth columns). First to third rows represent hourly, daily, and monthly TRs. Fe is on small size grains. Dashed grey lines in the second and third rows correspond to lines in the first row and are shown for easy comparison.

Hydrodynamics are more controlled by TR, and decrease of TR could amplify the effect of SR. Thus, at places where the human activities are intensive, a higher TR should be preferentially considered. 2. Environmental risk of Cr(VI) to river is underestimated with decreasing SR. Decrease of TR could neutralize a small portion of this underestimation. Difference of Cr(VI) mass to river among different SRs are amplified with the decrease of TR. Redox transformation of the Cr is more dependent on SR than TR. 3. With the same decrease of SR, risk is underestimated the most for medium size distribution. For small size distribution, upscaling has the most persistent effect on Cr(VI) mass to river since there is still a large amount of Cr(VI) discharged to river at 200 cm SR in this study.

100

400

5

1000 500 0 0 800

(e)

200 100 0 4

1500

10

300

4.5 t (year)

5

(b)

200cm-hourly 400 (c) 300 C (POC) (μmol/g)

200

2000

50cm-hourly

5

400 200 0 4

0 0 30

(f)

4.5 t (year)

5

0.08

0 0

10

0.06

(g)

25 20 15 10 4

5

50cm-monthly 0.16 (d)

100

10

600

50cm-daily

200

C (POC) (μmol/g)

C (DOC) (μmol/L)

(a)

300

0 0 C (DO) (μmol/L)

10cm-hourly 2500

C (DOC) (μmol/L)

C (DO) (μmol/L)

400

C (Fe(II)) (μmol/g)

1. Lower SRs and TRs tend to generate more uniform hydrodynamics.

C (Fe(II)) (μmol/g)

level dominated by human activities to that controlled mainly by natural regulations. Generated stochastic physical and geochemical fields and measured water levels are successively upscaled to larger spatial and temporal scales. Simulation results are of particular importance for help choosing the proper SRs and TRs when conducting monitoring and/or modeling of reactive contaminants and for improving accuracy of environmental risk assessment on contaminants discharge to river or other water bodies. Results also have general implications for riverine areas that have been polluted by other redox-sensitive pollutants such as As or Sb, which have high mobility at a low valence state and low mobility at a high valence state. Several important conclusions of this study can be drawn as follows:

4.5 t (year)

5

5

10

4.5 t (year)

5

(h)

0.04 0.02 0 4

Fig. 10. Calculated concentrations of DO (a, e), DOC (b, f), POC (c, g) and Fe(II) (d, h) at observation point P4 (690.05 m, 121.05 m) when Fe is on small size grains. The first row is for 10-year simulation while the second for the fifth-year simulation. 8

Journal of Hydrology 579 (2019) 124152

C. Yang, et al.

Fig. 11. Calculated DO mass from river to model domain and that contained in the model domain, when Fe is on small size grains. Grey dash lines in (b) and (c) are DO mass in (a) for easy comparison. Elfeki, A.M.M., Uffink, G., Lebreton, S., 2012. Influence of temporal fluctuations and spatial heterogeneity on pollution transport in porous media. Hydrogeol. J. 20 (2), 283–297. https://doi.org/10.1007/s10040-011-0796-0. Hammond, G.E., Lichtner, P.C., Rockhold, M.L., 2011. Stochastic simulation of uranium migration at the Hanford 300 Area. J. Contam. Hydrol. 120–21, 115–128. https:// doi.org/10.1016/j.jconhyd.2010.04.005. Hammond, G.E., Lichtner, P.C., Mills, R.T., 2014. Evaluating the performance of parallel subsurface simulators: an illustrative example with PFLOTRAN. Water Resour. Res. 50 (1), 208–228. https://doi.org/10.1002/2012WR013483. Herzog, S.P., Higgins, C.P., McCray, J.E., 2015. Engineered streambeds for induced hyporheic flow: enhanced removal of nutrients, pathogens, and metals from urban streams. J. Environ. Eng. 42 (1), 1–10. Johnson, K.H., 2016. Groundwater contaminant plume maps and volumes, 100-K and 100-N Areas, Hanford Site, Washington (2016-1161). Retrieved from Reston, VA: http://pubs.er.usgs.gov/publication/ofr20161161. Li, L., Peters, C.A., Celia, M.A., 2006. Upscaling geochemical reaction rates using porescale network modeling. Adv. Water Resour. 29 (9), 1351–1370. https://doi.org/10. 1016/j.advwatres.2005.10.011. Li, L., Steefel, C.I., Yang, L., 2008. Scale dependence of mineral dissolution rates within single pores and fractures. Geochim. Cosmochim. Acta 72 (2), 360–377. https://doi. org/10.1016/j.gca.2007.10.027. Li, L., Steefel, C.I., Kowalsky, M.B., Englert, A., Hubbard, S.S., 2010. Effects of physical and geochemical heterogeneities on mineral transformation and biomass accumulation during biostimulation experiments at Rifle, Colorado. J. Contam. Hydrol. 112 (1–4), 45–63. https://doi.org/10.1016/j.jconhyd.2009.10.006. Li, L., Gawande, N., Kowalsky, M.B., Steefel, C.I., Hubbard, S.S., 2011. Physicochemical heterogeneity controls on uranium bioreduction rates at the field scale. Environ. Sci. Technol. 45 (23), 9959–9966. https://doi.org/10.1021/es201111y. Lin, J.G., Chen, S.Y., Su, C.R., 2003. Assessment of sediment toxicity by metal speciation in different particle-size fractions of river sediment. Water Sci. Technol. 47 (7–8), 233–241. Liu, C.X., Shang, J.Y., Shan, H.M., Zachara, J.M., 2014. Effect of subgrid heterogeneity on scaling geochemical and biogeochemical reactions: a case of U(VI) desorption. Environ. Sci. Technol. 48 (3), 1745–1752. https://doi.org/10.1021/es404224j. Liu, Y., Xu, F., Liu, C., 2017. Coupled hydro-biogeochemical processes controlling Cr reductive immobilization in Columbia River hyporheic zone. Environ. Sci. Technol. 51 (3), 1508–1517. https://doi.org/10.1021/acs.est.6b05099. Ma, R., Zheng, C., Prommer, H., Greskowiak, J., Liu, C., Zachara, J., Rockhold, M., 2010. A field-scale reactive transport model for U(VI) migration influenced by coupled multirate mass transfer and surface complexation reactions. Water Resour. Res. 46. https://doi.org/10.1029/2009wr008168. Martincic, D., Kwokal, Z., Branica, M., 1990. Distribution of zinc, lead, cadmium and copper between different size fractions of sediments. 1. The Limski Kanal (North Adriatic Sea). Sci. Total Environ. 95, 201–215. https://doi.org/10.1016/00489697(90)90065-3. Molnar, D.K., Julien, P.Y., 2000. Grid-size effects on surface runoff modeling. J. Hydrol. Eng. 5 (1), 8–16. https://doi.org/10.1061/(Asce)1084-0699(2000) 5:1(8). Moser, D.P., et al., 2003. Biogeochemical processes and microbial characteristics across groundwater−surface water boundaries of the hanford reach of the Columbia River. Environ. Sci. Technol. 37 (22), 5127–5134. https://doi.org/10.1021/es034457v. Okweye, P.S., Garner, K.G., Overton, A.S., Moss, E.M., 2016. Factor-cluster analysis and effect of particle size on total recoverable metal concentration in sediments of the lower Tennessee River Basin. Comput. Water, Energy, Environ. Eng. 05 (01), 17.

4. DO content is higher at lower SR, and more dependent on SR at the regional scale than at the local scale such as at P4. More sensitive increase of Fe(II) with decreasing SR than with decreasing TR results in the stronger control of SR than TR on generation of Cr(III). 5. In general, for reactive transport, the effect of geochemical heterogeneity is more important than that of physical heterogeneity, which, in turn, is also more critical than that of the temporal variability of the hydrodynamics. Geochemical heterogeneity could amplify the risk caused by physical heterogeneity. Declaration of Competing Interest The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper. Acknowledgements This study was supported by the National Natural Science Foundation of China, China (Grant No. 41807198, 41877183) and by the Guangdong Provincial Key Laboratory of Soil and Groundwater Pollution Control, China (Grant No. 2017B030301012). We thank the tremendous help from Dr. Glenn Hammond and other developers for using PFLOTRAN in the hydrodynamics simulations. Appendix A. Supplementary data Supplementary data to this article can be found online at https:// doi.org/10.1016/j.jhydrol.2019.124152. References Bear, J., 1972. In: Dynamics of Fluids in Porous Media. Environmental Science Series (New York, 1972-). American Elsevier Pub. Co., New York, pp. 764. Brook, E.J., Moore, J.N., 1988. Particle-size and chemical control of as, Cd, Cu, Fe, Mn, Ni, Pb, and Zn in bed sediment from the Clark Fork River, Montana (USA). Sci. Total Environ. 76 (2-3), 247–266. https://doi.org/10.1016/0048-9697(88)90111-8. Deutsch, C.V., Journel, A.G., 1998. In: GSLIB Geostatistical Software Library and user's Guide. Oxford University Press, New York, pp. 1 computer laser optical disc 4 3/4 in. + 1 user's guide. Dresel, P.E., Qafoku, N.P., McKinley, J.P., Fruchter, J.S., Ainsworth, C.C., Liu, C., et al., 2008. Geochemical Characterization of Chromate Contamination in the 100 Area Vadose Zone at the Hanford Site (PNNL-17674). Retrieved from Richland, WA.

9

Journal of Hydrology 579 (2019) 124152

C. Yang, et al.

Water Res. Pl-Asce 119 (2), 141–157. https://doi.org/10.1061/(Asce)07339496(1993)119:2(141). Wen, X.H., GomezHernandez, J.J., 1996. Upscaling hydraulic conductivities in heterogeneous media: an overview. J. Hydrol. 183 (1–2), R9–R32. https://doi.org/10. 1016/S0022-1694(96)80030-8. Wen, H., Li, L., 2017. An upscaled rate law for magnesite dissolution in heterogeneous porous media. Geochim. Cosmochim. Acta 210, 289–305. https://doi.org/10.1016/j. gca.2017.04.019. White, A.F., Brantley, S.L., 2003. The effect of time on the weathering of silicate minerals: why do weathering rates differ in the laboratory and field? Chem. Geol. 202 (3–4), 479–506. https://doi.org/10.1016/j.chemgeo.2003.03.001. Williams, M.D., Rockhold, M.L., Thorne, P.D., Chen, Y., 2008. Three-Dimensional Groundwater Models of the 300 Area at the Hanford Site. Pacific Northwest National Laboratory, Washington State. Yabusaki, S.B., Fang, Y.L., Waichler, S.R., 2008. Building conceptual models of field-scale uranium reactive transport in a dynamic vadose zone-aquifer-river system. Water Resour. Res. 44 (12). https://doi.org/10.1029/2007wr006617. Yang, C., Zhang, Y.K., Liu, Y.Y., Yang, X.F., Liu, C.X., 2018. Model-based analysis of the effects of dam-induced river water and groundwater interactions on hydro-biogeochemical transformation of redox sensitive contaminants in a hyporheic zone. Water Resour. Res. 54 (9), 5973–5985. https://doi.org/10.1029/2018wr023286. Yao, Q.Z., Wang, X.J., Jian, H.M., Chen, H.T., Yu, Z.G., 2015. Characterization of the particle size fraction associated with heavy metals in suspended sediments of the Yellow River. Int. J. Environ. Res. Pub He 12 (6), 6725–6744. https://doi.org/10. 3390/ijerph120606725. Yu, Z.B., et al., 2014. Spatial and temporal scale effect in simulating hydrologic processes in a watershed. J. Hydrol. Eng. 19 (1), 99–107. https://doi.org/10.1061/(Asce)He. 1943-5584.0000762. Zachara, J.M., Chen, X.Y., Murray, C., Hammond, G., 2016. River stage influences on uranium transport in a hydrologically dynamic groundwater-surface water transition zone. Water Resour. Res. 52 (3), 1568–1590. https://doi.org/10.1002/ 2015WR018009. Zhang, D., 2002. Stochastic Methods for Flow in Porous Media: Coping With Uncertainties. Academic, , San Diego, Calif, London, pp. 350. Zheng, C., Jiao, J.J., 1998. Numerical simulation of tracer tests in heterogeneous aquifer. J. Environ. Eng.-Asce 124 (6), 510–516. https://doi.org/10.1061/(Asce)07339372(1998)124:6(510).

https://doi.org/10.4236/cweee.2016.51002. Perujo, N., Sanchez-Vila, X., Proia, L., Romani, A.M., 2017. Interaction between physical heterogeneity and microbial processes in subsurface sediments: a laboratory-scale column experiment. Environ. Sci. Technol. 51 (11), 6110–6119. https://doi.org/10. 1021/acs.est.6b06506. Salamon, P., Fernandez-Garcia, D., Gomez-Hernandez, J.J., 2007. Modeling tracer transport at the MADE site: the importance of heterogeneity. Water Resour. Res. 43 (8). https://doi.org/10.1029/2006wr005522. Salehikhoo, F., Li, L., Brantley, S.L., 2013. Magnesite dissolution rates at different spatial scales: the role of mineral spatial distribution and flow velocity. Geochim. Cosmochim. Ac 108, 91–106. https://doi.org/10.1016/j.gca.2013.01.010. Salomons, W., Förstner, U., 1984. In: Metals in the Hydrocycle. Springer-Verlag, Berlin; New York, pp. 349. Shang, J.Y., Liu, C.X., Wang, Z.M., Zachara, J.M., 2011. Effect of grain size on uranium (VI) surface complexation kinetics and adsorption additivity. Environ. Sci. Technol. 45 (14), 6025–6031. https://doi.org/10.1021/es200920k. Shuai, P., et al., 2019. Dam operations and subsurface hydrogeology control dynamics of hydrologic exchange flows in a regulated river reach. Water Resour. Res. https://doi. org/10.1029/2018wr024193. Singh, A.K., Hasnain, S.I., Banerjee, D.K., 1999. Grain size and geochemical partitioning of heavy metals in sediments of the Damodar River – a tributary of the lower Ganga, India. Environ. Geol. 39 (1), 90–98. https://doi.org/10.1007/s002540050439. Smoot, J.L., Biebesheimer, F.H., Eluskie, J.A., Simpkin, T., Spiliotopoulos, A., Tonkin, M. J., 2011. Groundwater Remediation at the 100-HR-3 Operable Unit, Hanford Site, Washington, USA - 11507 (CHPRC-01149-FP). Retrieved from Richland, WA: http:// pubs.er.usgs.gov/publication/ofr20161161. Stenta, H.R., Riccardi, G.A., Basile, P.A., 2017. Grid size effects analysis and hydrological similarity of surface runoff in flatland basins. Hydrol. Sci. J. 62 (11), 1736–1754. https://doi.org/10.1080/02626667.2017.1349315. Tessier, A., Campbell, P.G.C., Bisson, M., 1982. Particulate trace-metal speciation in stream sediments and relationships with grain-size – implications for geochemicalexploration. J. Geochem. Explor. 16 (2), 77–104. https://doi.org/10.1016/03756742(82)90022-X. Truex, M.J., et al., 2015. Assessment of Hexavalent Chromium Natural Attenuation for the Hanford Site 100 Area. PNNL-24705, Richland, WA. Tsai, L.J., Yu, K.C., Ho, S.T., Chang, J.S., Wu, T.S., 2003. Correlation of Particle Sizes and Metals Speciation in River Sediment. Diffuse Pollution Conference, Dublin. Vieux, B.E., Needham, S., 1993. Nonpoint-pollution model sensitivity to grid-cell size. J.

10