Diavik Waste Rock Project: Geostatistical Analysis of Sulfur, Carbon, and Hydraulic Conductivity Distribution in a Large-Scale Experimental Waste Rock Pile

Abstract

One of the large-scale field waste rock experiments (test piles) conducted as part of the Diavik Waste Rock Project was deconstructed, providing a spatially located set of geochemical, mineralogical, and particle-size distribution samples. Geostatistical analyses were conducted for sulfur and carbon content and saturated hydraulic conductivity, which affect the geochemical evolution of waste rock, to investigate the spatial dependence of these parameters. Analyses included population statistics, experimental semi-variogram estimation, and theoretical semi-variogram fitting. Population statistics were calculated for additional data sets from samples collected during the construction of the test piles. The population statistical analyses indicated that log-normal distribution provided the best fit for all investigated data sets. Experimental semi-variograms were estimated for the spatially located data set (test pile deconstruction) using the classical estimator, and theoretical semi-variograms were fitted. This investigation showed that the spatial distribution of sulfur, carbon, and hydraulic conductivity within the core of the test-pile experiments can be approximated using a log-normal distribution with a mean and standard deviation calculated using the samples collected during construction of the piles, and that little to no spatial relationship was present for these parameters at the scale of sampling. That the saturated hydraulic conductivity of the matrix material can be represented by the same statistical distribution throughout the test pile is significant because water flow, as well as mineral surface area and reactivity are dominantly controlled within the matrix portion of the test pile. Reactive transport simulations are included to demonstrate the influence of the matrix material on effluent geochemistry.

1. Introduction

Mineral and metal mining occurs throughout the world and generates large volumes of wastes estimated at a rate exceeding 100 billion tonnes per year [1]. Poor quality effluent from mine wastes is an extremely difficult problem to mitigate given the large volumes of waste involved and geochemical complexity of mineral weathering processes. Mine-waste rock stockpiled at the surface and exposed to oxygen and moisture can lead to high porewater concentrations of sulfate and transition metals and metalloids, and acidic effluent. To predict the quality of mine-waste effluent, empirical (e.g., scale factors) and mechanistic (e.g., reactive transport models) methods have been used (e.g., [2,3,4,5,6,7]). Most methods of effluent assessment or prediction are dependent on the mine-waste characteristics that contribute to the geochemical evolution of mine-waste rock, including the mineralogy of the waste and the physical characteristics controlling water flow. Investigations of waste rock weathering often rely on global parameters such as mean sulfur (S) and carbon (C) content and saturated hydraulic conductivity (K_S). This investigation was focused on these three parameters because all are directly related to the geochemical evolution of waste rock and all can be measured using simple well-defined techniques. Saturated hydraulic conductivity was included, even though most waste rock piles are unsaturated, because it can be estimated using particle-size distribution (PSD) and is a parameter that is readily incorporated into models (compared to D₁₀—the diameter of particle size for which 10% of a given sample is smaller—for example). Limited investigations have been conducted regarding the heterogeneity of these fundamental parameters. Fala et al. [8] investigated the role of heterogeneity in hydrogeological properties (K_S, suction, and volumetric water content) on water flow through a hypothetical waste rock pile and found that calculated flow patterns generally followed the distribution of K_S within a waste rock pile (i.e., higher flow corresponding to higher K_S zones). The work of Fala et al. [8] showed that coupling a heterogeneous K_S distribution with a heterogeneous S and/or C distribution in a mechanistic implementation could define further the role of mineralogical and K_S heterogeneity in the geochemical evolution of mine-waste rock.

There has been very little research on waste rock piles to indicate the spatial distribution of mineralogical and physical parameters. Stockwell et al. [9] deconstructed a 12 m high waste rock pile at Key Lake, Saskatchewan, and collected samples to characterize physical and geochemical parameters. The results of the study indicated highly variable paste pH and aqueous SO₄ concentrations that showed no spatial dependence within the pile, and a correlation between these geochemical parameters and particle-size distribution was not identified. Khalil et al. [10] collected mine waste samples at the abandoned Ketarra mine in Morocco to assess the extent of metal impacts. They suggested a normal statistical distribution was appropriate for solid samples collected at their study site, where samples were widely distributed and collected from near the waste rock pile surface, and therefore may not be indicative of distribution through the depth of the pile. Marescotti et al. [11] used geostatistics to quantify the distribution of metals in waste rock samples from a sulfidic waste rock pile. Their semi-variogram work indicated spatial relationships for most of the assessed parameters, but, similar to the work of Khalil et al. [10], samples were collected close to the surface of the waste rock pile. Recent work by Blannin et al. [12] used geostatistics to assess sampling strategies for investigation of spatial variability of selected metals in a tailings impoundment. Although this study focused on mill tailings, which fundamentally differ from waste rock, the results indicated that geostatistics can be successfully used in the assessment of anthropogenically deposited materials.

The Diavik Waste Rock Project experiments included two large-scale test piles that were constructed from 2005 to 2006 at the Diavik Diamond Mine, NT, Canada (Diavik), to investigate the geochemical evolution of low-sulfide waste rock. The densely instrumented test piles were approximately 50 m × 60 m × 15 m and were constructed with waste rock separated by S content: Type I (wt.% S < 0.04), Type II (0.04–0.08 wt.% S), and Type III (wt.% S > 0.08). One of the test pile experiments was constructed entirely of Type I material (average 0.035 wt.% S), and a second test pile experiment was constructed entirely of Type III material (average 0.053 wt.% S) [13]. The waste rock used to construct the test-pile experiments consisted of approximately 75% granite (primarily quartz, K-feldspar, and albite), 14% pegmatitic granite (primarily quartz, K-feldspar, and albite), 10% biotite schist (primarily albite, quartz, and biotite), and 1% diabase (primarily plagioclase) [13,14,15]. The metal sulfide minerals, dominated by pyrrhotite [Fe_0.852Ni_0.004Co_0.001S] [16], with small contents of chalcopyrite [CuFeS₂], sphalerite [ZnS], and pentlandite [(Fe,Ni)₉S₈], were primarily contained in the biotite schist [15]. Mineralogical analyses indicated C was present only as carbonate [16].

The Type I test pile was deconstructed in 2014 to, in part, investigate the progress of mineral weathering within the experiment. The test pile was excavated in a series of benches that allowed sample collection for geochemical parameters and PSD throughout the pile [17].

Effluent geochemistry from the Type III experiment was documented by Sinclair et al. [18] and generally consisted of increasing concentrations of solutes (e.g., SO₄, Ni, Co, Cu, Zn) and acidic pH (average annual pH ranging from 4.7 to 4.3 between 2008 and 2012) indicative of sulfide mineral weathering. Effluent geochemistry from the Type I experiment was documented by Bailey [19] and was generally similar to the Type III effluent with lower concentrations of solutes indicative of sulfide mineral weathering and circumneutral to slightly acidic pH (average annual pH ranging from 5.9 to 6.7 between 2008 and 2012).

A summary of Diavik Waste Rock Project studies, focused primarily on the Type III test pile experiment, includes initial characterization [13,20,21]; investigation of gas-transport [22,23]; temperature [24]; hydrology [25]; microbiology, geochemistry, and mineralogy [13,18,26,27,28,29,30]; and reactive transport simulation [31,32].

The purpose of this study is to characterize the spatial distribution of S, C, and K_S of the matrix material (i.e., diameter < 5 mm) based on samples collected during the construction of the Type I and Type III test piles and the deconstruction of the Type I test pile. The geostatistical analysis is used to characterize the geochemical and physical heterogeneity of a waste rock pile from the perspective of the influence of these parameters on effluent quality. Additionally, reactive transport simulations are included to demonstrate the influence of matrix material proportion on the geochemical evolution of waste rock.

2. Methodology

Geostatistical calculations for S, C, and K_S were conducted on samples collected for mineralogical analysis and PSD as part of the construction (Type I and Type III) and deconstruction (Type I) of the test-pile experiments. Samples collected during the construction of the Type I and/or Type III test piles are referred to as ‘construction’; samples collected during the deconstruction of the Type I test pile are referred to as ‘deconstruction’. Samples collected during the construction phase of the experiment generally had limited associated spatial information; samples collected during the Type I deconstruction were spatially located using a Real Time Kinetic Global Positioning System with +/−2 cm accuracy [17].

During construction of the test pile experiments, samples for analyses of S, C, and PSD were collected from the truck loads prior to being bulldozed over the tip face [13]. The majority of the construction samples were analyzed for S only. A total of 434 samples was collected for analysis of S during the construction of the test piles (247 samples from the Type I test pile and 187 from the Type III test pile). A total of 292 samples was collected for PSD analysis during construction (165 samples from the Type I test pile and 127 from the Type III test pile). Thirty of these PSD samples were analyzed for S and C in each grain size fraction.

The construction PSD samples consisted of 5–10 kg of the <50 mm fraction of the waste rock [20]. The geostatistical analysis documented here focuses on the PSD of the matrix material (<5 mm fraction) because flow through matrix material was found to dominate infiltration of water through the test pile [25]. Only PSD samples from the interior of the pile were included in the statistical analyses of K_S. Further detail regarding the collection of mineralogical and PSD samples during the construction of the test piles was documented by Smith et al. [13,20].

The test-pile experiments were constructed in phases to allow installation of instrumentation along tip faces of the piles at consistent intervals. Tip faces on which instruments were installed were referred to as instrumentation faces. Spatial information obtained for the construction samples was limited to the interval between two instrumentation faces of the test pile. As a result, approximate horizontal spatial location was obtained, but no vertical spatial information could be obtained, because each waste rock load was pushed over the tip face to construct the test pile.

During deconstruction of the Type I test pile (Figure 1), spatially located samples were collected for mineralogical and PSD analysis. The Type I test pile was excavated in benches approximately 3 m deep; subsequently, sloped trenches, approximately 20 m wide, were excavated to facilitate sample collection. Sample locations were selected to allow construction of vertical and horizontal profiles of mineralogical, physical, and geochemical parameters throughout the test-pile experiment. Samples for analysis of mineralogy and PSD were collected from approximately the same locations. A total of 562 samples for mineralogy analysis and 244 samples for PSD analysis was collected as part of the deconstruction program. A subset of 118 deconstruction mineralogy samples was selected for the S and C portion of the geostatistical analysis. Analysis of PSD was conducted for 141 spatially located samples from the Type I test pile deconstruction experiment.

Samples of approximately 70 kg were collected to assess the PSD of the <75 mm fraction of the waste rock in the Type I test pile after 8 years of operation [33]. Consistent with the construction samples, only the PSD of the matrix material (<5 mm fraction) was used for the geostatistical analysis. More details regarding the deconstruction program and specifics about the sampling procedures were documented by Atherton [17] and Barsi et al. [33].

Four methods of calculating K_S from the results of the PSD analysis were evaluated as part of this investigation. The methods of Hazen [34], Schlichter [35], Terzaghi [36], and Chapuis [37] (Equation (1)) were used to calculate K_S and compared to the results of the Neuner et al. [25] measurements.

$$K_S=2.4622{\left[\frac{D_{10}^2{e}^3}{\left(1+e\right)}\right]}^{0.7825}$$

(1)

where K_S is saturated hydraulic conductivity (cm s⁻¹), D₁₀ is the grain diameter at which 10% of the mass of a given sample is finer (mm), and e is the void ratio (m³ m⁻³). The use of the Kozeny–Carman equation [38] was rejected due to the requirement for a soil specific surface value. Neuner et al. [25] measured K_S on 18 matrix material samples from the Diavik experiments using constant-head permeameter tests.

The data analysis followed a standard geostatistical approach for characterization of parameters to quantify the spatial distribution of S, C, and K_S within the Type I test pile. The approach included calculation of sample statistics, assessment of stationarity, estimation of semi-variograms, and fitting of theoretical semi-variograms to estimated semi-variograms.

Experimental semi-variograms were estimated for the transformed data (i.e., Ln(S), Ln(C), and Ln(K_S)) from the deconstruction of the Type I test pile. Three common variogram estimation methods were applied to assess the capability of the estimators to describe the distribution of the waste rock in the Type I test pile. The classical variogram estimator [39] is:

$$\gamma \left(s\right)=\frac{1}{2N\left(s\right)}{\displaystyle \sum }_{i=1}^{N\left(s\right)}{\left(y_i-y_{i+s}\right)}^2$$

(2)

where N(s) is the number of data pairs separated by lag distance s, y_i is the log-normally transformed parameter at a given location, and y_i+s is the log-normally transformed parameter at a second location of lag distance (s) from y_i.

Two estimators generally considered to be more robust include the Cressie–Hawkins estimator ([40]; Equation (3)) and the squared median of the absolute deviations (SMAD) estimator ([41]; Equation (4)):

$$\gamma \left(s\right)=\frac{1}{2}{\left[\frac{1}{N\left(s\right)}{\displaystyle \sum }_{i=1}^{N\left(s\right)}{\left(y_i-y_{i+s}\right)}^{0.5}\right]}^4/\left(0.457+\frac{0.49}{N\left(s\right)}\right)$$

(3)

$$\gamma \left(s\right)=2.198{\left(median\left|y_i-y_{i+s}\right|\right)}^2$$

(4)

The theoretical semi-variogram models evaluated were pure nugget (Equation (5)), exponential (Equation (6)), and Gaussian ([42]; Equation (7)).

$$\gamma \left(s\right)={\sigma }_Y^2\left\{1-\delta \left(s\right)\right\}$$

(5)

where σ²_Y is the global variance δ(s) and is the Kronecker delta equaling 1 when lag distance s = 0 and 0 when s > 0.

$$\gamma \left(s\right)={\sigma }_0^2+{\sigma }_Y^2\left\{1-exp\left(-\raisebox{1ex}{$s$}\!\left/ \!\raisebox{-1ex}{$\lambda $}\right.\right)\right\}$$

(6)

where σ²₀ is the nugget and λ (m) is the correlation length.

$$\gamma \left(s\right)={\sigma }_0^2+{\sigma }_Y^2\left\{1-exp\left(-\raisebox{1ex}{$s^2$}\!\left/ \!\raisebox{-1ex}{${\lambda }^2$}\right.\right)\right\}$$

(7)

To demonstrate the influence of matrix material content on the geochemical evolution of the Diavik waste rock, a set of reactive transport simulations were conducted where the proportion of matrix material was varied. The simulations were conducted with the reactive transport code MIN3P [43] and based on the conceptual model presented by Wilson et al. [31]. The simulation mineralogy was consistent with Diavik waste rock including the sulfide minerals pyrrhotite, chalcopyrite, sphalerite, and pentlandite, and the host minerals calcite, dolomite, biotite, muscovite, and albite. Secondary minerals included Fe(III) (oxy)hydroxide, jarosite, gibbsite, gypsum, amorphous silica, and siderite. The simulations were conducted with matrix proportions of 0.1, 0.15, 0.2, 0.25, and 0.3. The simulations were conducted as a 2-D domain using a homogeneous matrix proportion representing the geochemical evolution of waste rock within the core of the pile.

3. Results

In the earth sciences, geostatistical analyses are traditionally applied to materials that occur as a result of natural geological processes. Here these techniques were applied to constructed waste rock piles to provide a framework for characterization of the distribution of physical, geochemical, and mineralogical properties throughout the interior of the test piles.

3.1. Calculation of Hydraulic Conductivity

Based on the K_S values calculated using the methods of Hazen [34], Schlichter [35], Terzaghi [36], and Chapuis [37] (

Table 1. Summary of K_S (m s⁻¹) estimates.

Parameter	Measured	Hazen	Schlichter	Terzaghi	Chapuis
geometric mean	9 × 10−6	5 × 10−5	5 × 10−6	3 × 10−4	4 × 10−5
standard deviation	1 × 10−5	1 × 10−5	1 × 10−6	8 × 10−5	8 × 10−6
minimum	2 × 10−6	2 × 10−5	3 × 10−6	1 × 10−4	2 × 10−5
maximum	3 × 10−5	9 × 10−5	1 × 10−5	6 × 10−4	7 × 10−5

Note: Measured values reprinted from [25], 2013, with permission from Elsevier.

), the Chapuis [37] equation was selected as the preferred method for calculation of K_S for this geostatistical analysis due to comparable geometric mean and standard deviation to the laboratory measured values of Neuner et al. [25]. The values calculated using the Schlichter [35] method were also comparable to the measured values but were rejected due to the reliance of the method on empirical parameters and temperature and a significantly lower standard deviation. The void ratio used in the Chapuis [37] equation was based on the mean porosity of 0.25 reported by Neuner et al. [25].

3.2. Statistical Distribution

Visual inspection of the S, C, and K_S profiles and calculation of sample statistics, including mean and standard deviation for each of the construction and deconstruction S, C, and K_S data sets, did not exhibit any significant trends in the distribution of the parameters (Figure 2). It is important to emphasize that the spatial plot of K_S reflects only the properties of the matrix, and not an effective bulk K_S at the scale of the test pile.

Frequency density histograms of the distribution of the Type I construction, Type I deconstruction, and Type III construction data sets for S, C, and K_S were observed to have skewness to the right, indicating a non-normal distribution (Figure 3). Chi-square goodness-of-fit tests and quantile–quantile plots for each of the log-normally transformed data sets (e.g., Y = Ln(K_S)) indicated acceptance of the normal distribution hypothesis in each case (Figure 3). The assumption of log-normal distribution was further analyzed using 10 subsets of 50 randomly selected values from the Type I deconstruction Ln(K_S). Chi-square goodness-of-fit tests conducted on these subsets indicated 90% of the subsets passed the additional chi-square testing, suggesting the assumption of log-normal distribution for the data sets was reasonable. The construction data set for Type I C contained less than 10 values; therefore, the distribution analysis was not conducted for this data set.

The results of the statistical characterization indicate that the distribution of S, C, and K_S in the waste rock at the Diavik test piles can be approximated using a log-normal distribution with mean and standard deviation calculated from samples collected during construction. Waste rock piles that are segregated according to S content may require individual distributions for S content associated with each section.

The correlation of S, C, and K_S at co-located sampling points was calculated. The S and C data exhibited low correlation with a correlation coefficient of −0.15 based on 120 co-located samples. The S and K_S and C and K_S data sets exhibited low correlation with correlation coefficients of −0.11 and −0.21, respectively, based on 85 co-located samples for each set. The correlation calculations indicated that the S content at a given location is not related to the C content or K_S at that same location and vice versa.

3.3. Stationarity

To characterize the spatial dependence of the Type I deconstruction data set, the assumption of stationarity in the first two moments (mean and variance) for Ln(S), Ln(C), and Ln(K_S) was considered. Stationarity is defined in this context as the absence of significant trends in mean and variance of the analyzed parameters with depth in the test-pile experiment. The stationarity of the Ln(S), Ln(C), and Ln(K_S) was tested by calculating the mean and variance at each bench for the Type I deconstruction data and compared to the overall mean and variance of the deconstruction data set (Figure 4). There was some fluctuation in the mean and variance about the average in each case.; however, it appears that stationarity in the Ln(S), Ln(C), and Ln(K_S) was a reasonable assumption for the Type I test pile. The spatial information for the Type I and Type III construction data sets was limited; however, analysis of the mean and variance of samples by instrumentation face compared to the overall mean and variance indicated that the assumption of stationarity appears reasonable for Ln(S) and Ln(K_S) of the construction samples. The stationarity of Ln(C) was not investigated using construction samples due to the low number of samples analyzed for this parameter during construction of the test piles.

3.4. Experimental Semi-Variogram Estimation

Horizontal and vertical semi-variograms were estimated for Ln(S), Ln(C), and Ln(K_S) using each method (Figure 5). Horizontal semi-variograms were estimated from parameter pairs that were separated by benches. The maximum vertical separation of horizontal Ln(S) and Ln(C) pairs was 2.6 m; the maximum vertical separation of horizontal Ln(K_S) pairs was 3.7 m. Vertical semi-variograms were estimated from parameter pairs that were horizontally separated by a maximum of 3 m. Lag intervals were selected based on number of pairs available to provide a reasonable estimate. The horizontal lag interval selected for Ln(S) and Ln(C) was 4 m; the vertical lag interval selected for Ln(S) and Ln(C) was 1.5 m. The lag intervals selected for Ln(K_S) were 3 m and 1 m for horizontal and vertical lags, respectively. Visual inspection of the experimental semi-variograms indicated that the estimates calculated by the classical method provided the most stable trends in each case, and as such the classical method was selected as the estimator for fitting of the theoretical variograms.

3.5. Theoretical Semi-Variogram Fitting

Visual inspection of the experimental semi-variogram estimates indicated little evidence of spatial dependence in Ln(S) and Ln(C) in the horizontal or vertical directions; limited evidence of spatial dependence is observed for Ln(K_S). To test this observation, three common theoretical semi-variogram models were applied to the estimates. Models were evaluated for fit by calculating the standard error of the regression (sum of squares). The pure nugget model (Equation (5)) was selected for analysis because visual inspection of the semi-variograms for Ln(S) and Ln(C) indicated no spatial dependence at the selected lag intervals. The exponential model (Equation (6)) was selected based on the success of other researchers in fitting the model to semi-variograms of Ln(K_S) for fine-grained sands (e.g., [44,45,46]). The Gaussian model (Equation (7)) was selected based on visual inspection of the experimental Ln(K_S) semi-variograms.

Fitting of the estimated semi-variogram for Ln(S) and Ln(C) indicated that in all cases there was very little to no evidence of spatial dependence at scales of 4 m (horizontal) and 1.5 m (vertical) with negligible differences in standard errors for the three theoretical models. As a result, the pure nugget model was selected to represent the spatial dependence of Ln(S) and Ln(C) (Figure 6). The pure nugget model fitting parameter (σ²_Y) was 0.28 for Ln(S) (horizontal and vertical) and 0.05 for Ln(C) (horizontal and vertical). The horizontal and vertical semi-variograms for Ln(K_S) indicated limited spatial dependence with the Gaussian model providing the best fit with σ²_Y of 0.043, λ of 10 m, and σ²₀ of 0.025 (Figure 6).

3.6. Reactive Transport Simulations

The effluent geochemistry for selected parameters from the simulations (Figure 7) indicated that concentrations were influenced in proportion to the increase or decrease of matrix content for some parameters (e.g., Fe, Ni, Co, Cu). Sulfate and Ca concentrations decreased proportionally with the matrix material content but did not increase with increased matrix material content; the simulation results indicated that as solute concentrations increased, the precipitation of gypsum limited the maximum concentrations of Ca and SO₄ in the higher matrix proportion scenarios. For some parameters (e.g., Cu, Al) the matrix proportion significantly influenced the effluent concentrations, whereas the impact of an increase or decrease in the matrix proportion had little influence on the pH, which was constrained by the dissolution of carbonate minerals and aluminum-bearing phases. The proportion of matrix material within a pile is an important factor in assessing the geochemical evolution of a waste rock pile.

4. Discussion

The geostatistical characterization of S, C, and K_S was used to investigate the spatial distribution of these parameters that are fundamental to the assessment of the geochemical evolution of waste rock. The statistical distributions of S, C, and K_S were determined to be log-normal, which is common in natural systems [46,47,48].

4.1. Statistical Comparison of Construction Samples

Comparison of Ln(S), Ln(C), and Ln(K_S) for the Type I and Type III construction sample results indicated that the assumption of log-normally distributed data and stationarity in the first two moments was reasonable at each test pile. The variances of Ln(S) and Ln(C) were consistent between the Type I and Type III test piles. The means of Ln(S) and Ln(C) were not consistent between the test piles, reflecting the slightly differing S and C contents of the two piles. The mean and variance of Ln(K_S) were consistent between the two test piles. The similarity in parameter distribution, mean variance, and stationarity suggested that the K_S of the two test piles are statistically similar, and that PSD of the matrix material collected at the Type I pile could be used to describe the heterogeneity in K_S of the Type III pile despite difference in S content between the two piles. This consistency is expected given the similar blast design and transport methods used for the source material, the similar construction methods used for the test piles, and the similar mineralogy of host material in the two piles. This result could likely be applied to other waste rock piles where host rock mineralogy is relatively homogeneous.

4.2. Theoretical Semi-Variogram Fitting

The pure nugget theoretical semi-variogram model (Equation (5)) that provided the best fit to the Ln(S) and Ln(C) semi-variogram estimates was somewhat expected in view of the sampling density of the Type I deconstruction program, and the construction techniques used to build the test piles. The use of the pure nugget model indicates that any discernable trends in the spatial distribution of S and C were at a scale smaller than the sampling density. The results of this study indicate that following blasting and transport to the test-pile construction site, and deposition at the test pile, any spatial dependence in S and C in the rock prior to blasting at a scale of >2 m (approximate spatial distribution of samples collected in this study) would likely be lost.

The absence of spatial dependence in S and C could be influenced by factors including (i) mining and pile construction methods, (ii) distribution of the lithology, (iii) sampling density used to determine the spatial dependence of S and/or C, and/or (iv) estimated semi-variograms for S and C considering the spatial dependence in the pile as a whole, which may have eliminated smaller scale detail. The results of the spatial analyses highlighted the importance of considering sampling density that is based on the physical characteristics of the material being investigated. This geostatistical analysis indicated that higher sampling density would be required to investigate small-scale variability of mineralogy and physical parameters within waste rock.

The Gaussian model provided the best fit to the Ln(K_S) horizontal and vertical semi-variograms, with the spatial dependence being slightly more prominent in the horizontal direction. The poorer fit at larger lag distances in the vertical direction is likely the result of relatively low numbers of pairs (i.e., <10 for lag >7.5 m). Both directions were fit using a relatively large nugget value (59% of the underlying variance) and, coupled with the relatively small range of K_S values (within one order of magnitude), suggests limited spatial dependence (i.e., the K_S of the matrix material was not depth dependent).

Research of the internal structure of waste rock piles has indicated in some cases the presence of low K_S lenses associated with traffic surfaces created during pile construction [49,50]. Internal traffic surfaces were not created as part of the construction of the test piles due to the method of construction; however, it is likely that the presence of surface and internal traffic surfaces could significantly alter the hydrogeologic regime within a pile.

Differences in sample collection methods could have also influenced the PSD results (i.e., construction versus deconstruction). Analysis of the stationarity of the mean and variance of Ln(K_S) indicated that K_S did not vary significantly with depth within the core of the Type I test pile experiment. The lack of depth dependence in K_S likely results from using only the PSD of the matrix material within the pile. Using the Chapuis equation (Equation (1)) and the PSD for particle sizes < 75 mm from the Type I test pile deconstruction, Barsi et al. [33] found that the bulk K_S increased with depth in the Type I test pile. The results of this study indicate that the matrix material K_S was not depth dependent.

4.3. Reactive Transport Simulations

The relative uniformity of matrix K_S in the test pile experiments is significant because matrix material is a controlling component in the test pile hydrology; most of the water flow is conveyed through the fine-grained matrix of the test pile. Furthermore, due to the large surface area of the fine fraction, the geochemical evolution of the pore water is dominantly controlled by the matrix material (e.g., [51]). The uniformity of the matrix material within the pile suggests all portions of the test pile that include matrix material could contribute to the solute loading, with the exception of stagnant zones that are isolated from flow by obstructions (e.g., matrix material located below large boulders). This observation is consistent with previous research (e.g., [52]) that indicated relatively uniform solute loading (approximately one order of magnitude) from different parts of a single pile. The results of the reactive transport simulations indicate that most parameters are at least somewhat influenced by the matrix material content within the pile and are consistent with the observation of the importance of matrix material proportion within a waste rock pile.

5. Conclusions

Analysis of the population statistics from samples collected during the construction of the Type I and Type III test piles and deconstruction of the Type I test pile indicated that each data set was log-normally distributed with consistent mean and variance among sub-populations of the same parameter.

Geostatistical analyses of the spatially located deconstruction samples indicated limited spatial dependence of the K_S field within the matrix material of the Type I test pile at a scale of 3 m horizontal and 1 m vertical. The K_S of the matrix material was not depth dependent. It was determined that no spatial dependence in the structures of the S and C distribution was present in the Type I test pile at a scale of 4 m horizontal and 1.5 m vertical. This result was expected given the methods used to construct the test-pile experiments. Sampling density during the deconstruction of the Type I test pile likely contributed to the results of the geostatistical characterization; however, given the extremely large size of typical waste rock piles, a denser characterization of the spatial distribution of S, C, and K_S is not likely to contribute significantly to an overall assessment of the geochemical evolution of waste rock.

The results of the geostatistical characterization suggest that the spatial distribution of S, C, and K_S in the test piles, which can be characterized as randomly heterogeneous, is likely due to construction methods; the distributions of S and C within the test pile could also be influenced by their distribution within the original lithology. The heterogeneity of S, C, and K_S in the waste rock at the Diavik test piles can be approximated using a log-normal distribution with mean and standard deviation calculated from samples collected during construction. This is an important conclusion, suggesting that heterogeneity representative of physical and mineralogical conditions of a pile can be implemented in a spatially random manner according to measured log-normal distribution statistics. Waste rock piles that are portioned according to S content may require individual distributions for S content associated with each segregated component. The lack of spatial dependence in matrix K_S is significant because of the importance of the matrix material in controlling water flow through the test-pile experiments. Due to the large surface area of the matrix material, the geochemical evolution of pore water is dominantly controlled by this finer fraction. Reactive transport simulations conducted using differing proportions of matrix material indicate that most parameters are at least somewhat influenced by the matrix material content within the pile.

The uniformity of the matrix material throughout the pile suggests that contributions to the solute loading from the pile would be dependent on the spatial distribution of matrix material in a pile. The proportion of matrix material within a pile is an important factor in assessing the geochemical evolution of a waste rock pile.