Influence of the amount of available K data on uncertainty about contaminant transport prediction

The effect of increasing the amount of hydraulic conductivity (K) data on the level of uncertainty about contaminant concentration predictions was examined. Contaminant concentrations were predicted using Monte Carlo numerical simulations with progressively more K information, and were compared to reference data obtained from a physical experiment. Increasing the amount of K data used in generating the aquifer's K field improved the measure of confidence in the predicted contaminant concentration. The normalized concentration interval with 95% confidence decreased from 1.29 to 0.75 when the amount of K data increased from 20 to 71. However, the trade-off is the increase in the costs of aquifer investigation. Conducting Monte Carlo numerical simulations may help design an optimal aquifer investigation scheme, so that it can provide enough input data to predict concentrations with the desired level of confidence and at the same time avoid the excessive costs of aquifer investigation.     

Monte Carlo evaluation of microbial-mediated contaminant reactions in heterogeneous aquifers

Monte Carlo simulations are conducted to evaluate microbial-mediated contaminant reactions in an aquifer comprised of spatially variable microbial biomass concentrations, aquifer hydraulic conductivities, and initial electron donor/acceptor concentrations. A finite element simulation model is used that incorporates advection, dispersion, and Monod kinetic expressions to describe biological processes. Comparisons between Monte Carlo simulations of heterogeneous systems and simulations using homogeneous formulation of the same two-dimensional transport problem are presented. For the assumed set of parameters, physical aquifer heterogeneity is found to have a minor effect on the mass of contaminant biodegraded/transformed when compared to a homogeneous system; however, it noticeably changes the dispersion, skewness, and peakness of contaminant concentration distributions. Similarly, for low microbial growth rate, given favorable microbial growth characteristics, biological heterogeneity has minor effect on the mass of contaminant biodegraded/transformed when compared to a homogeneous system. On the other hand, when higher effective growth rates are assumed, biological heterogeneity and spatial heterogeneities in essential electron donor/acceptors reduce the efficiency of biotic contaminant reactions; consequently, model simulations derived from heterogeneous biomass distributions predict remediation time scales that are longer than those simulated for homogeneous systems. When correlations between physical aquifer and biological heterogeneities are considered, the assumed correlation affects predicted mean and variance of contaminant concentration and biomass distributions. For example, an assumed negative correlation between hydraulic conductivity and the initial biomass distribution produces a plume where less efficient biotic contaminant reactions occur at the leading edge of the plume; this is consistent with less degradation/transformation occurring over regions of higher groundwater velocities. However, the presence and absence of these correlations do not appear to affect the efficiency of microbial-mediated contaminant attenuation.     

Optimal allocation of computational resources in hydrogeological models under uncertainty

Flow and transport models in heterogeneous geological formations are usually large-scale with excessive computational complexity and uncertain characteristics. Uncertainty quantification for predicting subsurface flow and transport often entails utilizing a numerical Monte Carlo framework, which repeatedly simulates the model according to a random field parameter representing hydrogeological characteristics of the aquifer. The physical resolution (e.g. spatial grid resolution) for the simulation is customarily chosen based on recommendations in the literature, independent of the number of Monte Carlo realizations. This practice may lead to either excessive computational burden or inaccurate solutions. We develop an optimization-based methodology that considers the trade-off between the following conflicting objectives: time associated with computational costs, statistical convergence of the model prediction and physical errors corresponding to numerical grid resolution. Computational resources are allocated by considering the overall error based on a joint statistical–numerical analysis and optimizing the error model subject to a given computational constraint. The derived expression for the overall error explicitly takes into account the joint dependence between the discretization error of the physical space and the statistical error associated with Monte Carlo realizations. The performance of the framework is tested against computationally extensive simulations of flow and transport in spatially heterogeneous aquifers. Results show that modelers can achieve optimum physical and statistical resolutions while keeping a minimum error for a given computational time. The physical and statistical resolutions obtained through our analysis yield lower computational costs when compared to the results obtained with prevalent recommendations in the literature. Lastly, we highlight the significance of the geometrical characteristics of the contaminant source zone on the optimum physical and statistical resolutions.     

Time-related capture zones for radial flow in two dimensional randomly heterogeneous media

We consider the effect of randomly heterogeneous hydraulic conductivity on the spatial location of time-related capture zones (isochrones) for a non-reactive tracer in the steady-state radial flow field due to a pumping well in a confined aquifer. A Monte Carlo (MC) procedure is used in conjunction with FFT-based spectral methods. The log hydraulic conductivity field is assumed to be Gaussian and stationary, with isotropic exponential correlation. Various degrees of domain heterogeneity are considered and stability and accuracy of the MC procedure is examined. The location of an isochrone becomes uncertain due to heterogeneity, and it is strongly influenced by hydraulic conductivity variance. The probability that a particle released at a point in the aquifer is pumped by the well within a given time is identified. We propose a new expression for the probabilistic spatial distribution of isochrones, which is formally similar to the analytical solution for a uniform medium and takes into account the effects of heterogeneity.     

Combining pump-and-treat and physical barriers for contaminant plume control

A detailed analysis is presented of the hydraulic efficiency of plume management alternatives that combine a conventional pump-and-treat system with vertical, physical hydraulic barriers such as slurry walls or sheet piles. Various design settings are examined for their potential to reduce the pumping rate needed to obtain a complete capture of a given contaminated area. Using established modeling techniques for flow and transport, those barrier configurations (specified by location, shape, and length) that yield a maximum reduction of the pumping rate are identified assuming homogeneous aquifer conditions. Selected configurations are further analyzed concerning their hydraulic performance under heterogeneous aquifer conditions by means of a stochastic approach (Monte Carlo simulations) with aquifer transmissivity as a random space function. The results show that physical barriers are an appropriate means to decrease expected (mean) pumping rates, as well as the variance of the corresponding pumping rate distribution at any given degree of heterogeneity. The methodology presented can be transferred easily to other aquifer scenarios, provided some basic premises are fulfilled, and may serve as a basis for reducing the pumping rate in existing pump-and-treat systems.     

Probabilistic risk assessment of contaminant transport in groundwater and vapour intrusion following remediation of a contaminant source

Community-scale simulations were performed to investigate the risk to groundwater and indoor air receptors downgradient of a contaminated site following the remediation of a long-term source. Six suites of Monte Carlo simulations were performed using a numerical model that accounted for groundwater flow, reactive solute transport, soil gas flow, and vapour intrusion in buildings. The model was applied to a three-dimensional, community-scale (250 m × 1000 m × 14 m) domain containing heterogeneous, spatially correlated distributions of the hydraulic conductivity, fraction of organic carbon, and biodegradation rate constant, which were varied between realizations. Analysis considered results from both individual realizations as well as the suite of Monte Carlo simulations expressed through several novel, integrated parameters, such as the probability of exceeding a regulatory standard in either groundwater or indoor air. Results showed that exceedance probabilities varied considerably with the consideration of biodegradation in the saturated zone, and were less sensitive to changes in the variance of hydraulic conductivity or the incorporation of heterogeneous distributions of organic carbon at this spatial scale. A sharp gradient in exceedance probability existed at the lateral edges of the plumes due to variability in lateral dispersion, which defined a narrow region of exceedance uncertainty. Differences in exceedance probability between realizations (i.e., due to heterogeneity uncertainty) were similar to differences attributed to changes in the variance of hydraulic conductivity or fraction of organic carbon. Simulated clean-up times, defined by reaching an acceptable exceedance probability, were found to be on the order of decades to centuries in these community-scale domains. Results also showed that the choice of the acceptable exceedance probability level (e.g., 1 vs. 5 %) would likely affect clean up times on the order of decades. Moreover, in the scenarios examined here, the risk of exceeding indoor air standards was greater than that of exceeding groundwater standards at all times and places. Overall, simulations of coupled transport processes combined with novel spatial and temporal quantification metrics for Monte Carlo analyses, provide practical tools for assessing risk in wider communities when considering site remediation.     

An evaluation of methodologies for the generation of stochastic hydraulic conductivity fields in highly heterogeneous aquifers

 Stochastic techniques, such as Monte Carlo experiments, are more and more frequently used for the study of flow and transport in heterogeneous aquifers. When the aquifer is composed of distinct hydrofacies, a common way to model heterogeneity is to first generate equally-possible hydrofacies fields, and then convert these hydrofacies fields into hydraulic conductivity (K) fields by assigning a single K value to each facies. This technique assumes relative homogeneity of K within each facies but may not be appropriate for the most conductive facies that often exhibits substantial variability. In this paper, we assessed the impacts of assigning multiple random K, rather than a uniform K value, to the highly conductive facies on the results of a flow and transport model. A set of fifty stochastic hydrofacies maps depicting an environment similar to the Snake River Plain aquifer (SRPA) in south-east Idaho were generated. Simulations demonstrated that a uniform K value, if carefully chosen, can reasonably reproduce the specific discharges and early particle arrival times produced by multiple K values. Yet, the results obtained with a uniform K value are dramatically less variable than those obtained with multiple K values. It is therefore concluded that stochastic simulations with uniform K assigned to the most conductive and variable facies do not necessarily portray the entire uncertainty in the analyses.     

A numerical Eulerian method of moment for solute transport in a nonstationary dual-porosity medium

An Eulerian perturbation approach was applied to develop a method of moment for solute transport in a nonstationary, fractured medium. The conceptualized fractured medium is described through a dual-porosity model. Stochastic governing equations for mean concentration and concentration covariance were analytically derived to the first-order accuracy of log-conductivity variance and solved with a numerical method––a finite difference method. The developed method is called a numerical Eulerian method of moment (NEMM). This method was compared with the stationary transport theory [Water Resour. Res. 36(7) (2000) 1665] for predicting mean concentration and its spatial moments. The comparison indicated that the two methods matched very well in predicting first and second spatial moments. NEMM solutions were also compared with Monte Carlo simulations for solute transport in stationary fractured media. The results of the two methods were consistent for calculating small log conductivity variance. The theory was then used to study effects of various parameters and nonstationarity of the medium on flow and transport processes. Results indicated that medium nonstationarity would significantly influence the solute transport process. The nonstationary transport theory relaxes many assumptions adopted in stationary theories and paves the way for applying the NEMM to many environmental projects, especially in analyzing uncertainty of solute transport.     

Modeling groundwater velocity uncertainty in nonstationary composite porous media

Despite the intensive research over the past decades in the field of stochastic subsurface hydrology, our ability to analyze and model heterogeneous groundwater systems remains limited. Most existing theories are either too restrictive to handle practical complexity or too expensive to be applied to realistic problem sizes. In this paper we present approximate, closed-form equations that allow modeling 2D nonstationary flows in statistically inhomogeneous aquifers, including composite aquifers containing multiple zones characterized by different statistical models. The composite representation has the effect of decreasing the variance of deviations from the mean, relaxing the limitation of the small-perturbation assumption. The simple formulas are illustrated with a number of examples and compared with a corresponding first-order nonstationary numerical analysis and Monte Carlo simulation. The results show that, despite the gross simplifications, the closed-form equations are robust and able to capture complex variance dynamics, reproducing surprisingly well the first-order numerical solutions and the Monte Carlo simulation even in highly nonstationary, variable situations.     

Flow in unsaturated random porous media,nonlinear numerical analysis and comparison to analytical stochastic models

This work presents a rigorous numerical validation of analytical stochastic models of steady state unsaturated flow in heterogeneous porous media. It also provides a crucial link between stochastic theory based on simplifying assumptions and empirical field and simulation evidence of variably saturated flow in actual or realistic hypothetical heterogeneous porous media. Statistical properties of unsaturated hydraulic conductivity, soil water tension, and soil water flux in heterogeneous soils are investigated through high resolution Monte Carlo simulations of a wide range of steady state flow problems in a quasi-unbounded domain. In agreement with assumptions in analytical stochastic models of unsaturated flow, hydraulic conductivity and soil water tension are found to be lognormally and normally distributed, respectively. In contrast, simulations indicate that in moderate to strong variable conductivity fields, longitudinal flux is highly skewed. Transverse flux distributions are leptokurtic. the moments of the probability distributions obtained from Monte Carlo simulations are compared to modified first-order analytical models. Under moderate to strong heterogeneous soil flux conditions (σ²_y≥1), analytical solutions overestimate variability in soil water tension by up to 40% as soil heterogeneity increases, and underestimate variability of both flux components by up to a factor 5. Theoretically predicted model (cross-)covariance agree well with the numerical sample (cross-)covarianaces. Statistical moments are shown to be consistent with observed physical characteristics of unsaturated flow in heterogeneous soils.©1998 Elsevier Science Limited. All rights reserved     

Polynomial chaos expansion for global sensitivity analysis applied to a model of radionuclide migration in a randomly heterogeneous aquifer

We perform global sensitivity analysis (GSA) through polynomial chaos expansion (PCE) on a contaminant transport model for the assessment of radionuclide concentration at a given control location in a heterogeneous aquifer, following a release from a near surface repository of radioactive waste. The aquifer hydraulic conductivity is modeled as a stationary stochastic process in space. We examine the uncertainty in the first two (ensemble) moments of the peak concentration, as a consequence of incomplete knowledge of (a) the parameters characterizing the variogram of hydraulic conductivity, (b) the partition coefficient associated with the migrating radionuclide, and (c) dispersivity parameters at the scale of interest. These quantities are treated as random variables and a variance-based GSA is performed in a numerical Monte Carlo framework. This entails solving groundwater flow and transport processes within an ensemble of hydraulic conductivity realizations generated upon sampling the space of the considered random variables. The Sobol indices are adopted as sensitivity measures to provide an estimate of the role of uncertain parameters on the (ensemble) target moments. Calculation of the indices is performed by employing PCE as a surrogate model of the migration process to reduce the computational burden. We show that the proposed methodology (a) allows identifying the influence of uncertain parameters on key statistical moments of the peak concentration (b) enables extending the number of Monte Carlo iterations to attain convergence of the (ensemble) target moments, and (c) leads to considerable saving of computational time while keeping acceptable accuracy.     

An efficient probabilistic finite element method for stochastic groundwater flow

We present an efficient numerical method for solving stochastic porous media flow problems. Single-phase flow with a random conductivity field is considered in a standard first-order perturbation expansion framework. The numerical scheme, based on finite element techniques, is computationally more efficient than traditional approaches because one can work with a much coarser finite element mesh. This is achieved by avoiding the common finite element representation of the conductivity field. Computations with the random conductivity field only arise in integrals of the log conductivity covariance function. The method is demonstrated in several two- and three-dimensional flow situations and compared to analytical solutions and Monte Carlo simulations. Provided that the integrals involving the covariance of the log conductivity are computed by higher-order Gaussian quadrature rules, excellent results can be obtained with characteristic element sizes equal to about five correlation lengths of the log conductivity field. Investigations of the validity of the proposed first-order method are performed by comparing nonlinear Monte Carlo results with linear solutions. In box-shaped domains the log conductivity standard deviation σ_Y may be as large as 1.5, while the head variance is considerably influenced by nonlinear effects as σ_Y approaches unity in more general domains.     

Comparison of approaches for predicting solute transport: sandbox experiments

The main purpose of this paper was to compare three approaches for predicting solute transport. The approaches include: (1) an effective parameter/macrodispersion approach (Gelhar and Axness 1983); (2) a heterogeneous approach using ordinary kriging based on core samples; and (3) a heterogeneous approach based on hydraulic tomography. We conducted our comparison in a heterogeneous sandbox aquifer. The aquifer was first characterized by taking 48 core samples to obtain local-scale hydraulic conductivity (K). The spatial statistics of these K values were then used to calculate the effective parameters. These K values and their statistics were also used for kriging to obtain a heterogeneous K field. In parallel, we performed a hydraulic tomography survey using hydraulic tests conducted in a dipole fashion with the drawdown data analyzed using the sequential successive linear estimator code (Yeh and Liu 2000) to obtain a K distribution (or K tomogram). The effective parameters and the heterogeneous K fields from kriging and hydraulic tomography were used in forward simulations of a dipole conservative tracer test. The simulated and observed breakthrough curves and their temporal moments were compared. Results show an improvement in predictions of drawdown behavior and tracer transport when the K tomogram from hydraulic tomography was used. This suggests that the high-resolution prediction of solute transport is possible without collecting a large number of small-scale samples to estimate flow and transport properties that are costly to obtain at the field scale.     

Accuracy of numerical simulations of contaminant transport in heterogeneous aquifers: A comparative study

This work deals with a comparison of different numerical schemes for the simulation of contaminant transport in heterogeneous porous media. The numerical methods under consideration are Galerkin finite element (GFE), finite volume (FV), and mixed hybrid finite element (MHFE). Concerning the GFE we use linear and quadratic finite elements with and without upwind stabilization. Besides the classical MHFE a new and an upwind scheme are tested. We consider higher order finite volume schemes as well as two time discretization methods: backward Euler (BE) and the second order backward differentiation formula BDF (2). It is well known that numerical (or artificial) diffusion may cause large errors. Moreover, when the Péclet number is large, a numerical code without some stabilising techniques produces oscillating solutions. Upwind schemes increase the stability but show more numerical diffusion. In this paper we quantify the numerical diffusion for the different discretization schemes and its dependency on the Péclet number. We consider an academic example and a realistic simulation of solute transport in heterogeneous aquifer. In the latter case, the stochastic estimates used as reference were obtained with global random walk (GRW) simulations, free of numerical diffusion. The results presented can be used by researchers to test their numerical schemes and stabilization techniques for simulation of contaminant transport in groundwater.     

Conventional and combined pump-and-treat systems under nonuniform background flow

We investigate the performance of vertical hydraulic barriers in combination with extraction wells for the partial hydraulic isolation of contaminated aquifer areas. The potential advantage of such combinations compared to a conventional pump-and-treat system has already been demonstrated in a previous study. Here we extend the scope of the performance analysis to the impact of uncertainty in the regional flow direction as well as to highly heterogeneous aquifer transmissivity distributions. In addition, two new well-barrier scenarios are proposed and analyzed. The hydraulic efficiency of the scenarios is rated based on the expected (mean) reduction of the pumping rate that is required to achieve downgradient contaminant capture. The uncertain spatial distribution of aquifer transmissivity is considered by means of unconditioned Monte Carlo simulations. The significance of uncertain background flow conditions is incorporated by computing minimized pumping rates for deviations of the regional flow direction up to 30 degrees from a normative base case. The results give an answer on how pumping rates have to be changed for each barrier-well combination in order to achieve robust systems. It is exposed that in comparison to installing exclusively wells, the barrier-supported approach generally yields savings in the (average) pumping rate. The particular efficiency is shown to be highly dependent on the interaction of variance and integral scale of transmissivity distribution, well and barrier position, as well as direction of background flow.     

Analytical solutions for the coefficient of variation of the volume-averaged solute concentration in heterogeneous aquifers

Under the assumption that local solute dispersion is negligible, a new general formula (in the form of a convolution integral) is found for the arbitrary k-point ensemble moment of the local concentration of a solute convected in arbitrary m spatial dimensions with general sure initial conditions. From this general formula new closed-form solutions in m=2 spatial dimensions are derived for 2-point ensemble moments of the local solute concentration for the impulse (Dirac delta) and Gaussian initial conditions. When integrated over an averaging window, these solutions lead to new closed-form expressions for the first two ensemble moments of thevolume-averaged solute concentration and to the corresponding concentration coefficients of variation (CV). Also, for the impulse (Dirac delta) solute concentration initial condition, the second ensemble moment of thesolute point concentration in two spatial dimensions and the corresponding CV are demonstrated to be unbound. For impulse initial conditions the CVs for volume-averaged concentrations axe compared with each other for a tracer from the Borden aquifer experiment. The point-concentration CV is unacceptably large in the whole domain, implying that the ensemble mean concentration is inappropriate for predicting the actual concentration values. The volume-averaged concentration CV decreases significantly with an increasing averaging volume. Since local dispersion is neglected, the new solutions should be interpreted as upper limits for the yet to be derived solutions that account for local dispersion; and so should the presented CVs for Borden tracers. The new analytical solutions may be used to test the accuracy of Monte Carlo simulations or other numerical algorithms that deal with the stochastic solute transport. They may also be used to determine the size of the averaging volume needed to make a quasi-sure statement about the solute mass contained in it.     

Analytical solutions to steady state unsaturated flow in layered,randomly heterogeneous soils via Kirchhoff transformation

In this study, we derive analytical solutions of the first two moments (mean and variance) of pressure head for one-dimensional steady state unsaturated flow in a randomly heterogeneous layered soil column under random boundary conditions. We first linearize the steady state unsaturated flow equations by Kirchhoff transformation and solve the moments of the transformed variable up to second order in terms of σ_Y and σ_β, the standard deviations of log hydraulic conductivity Y=ln(K_s) and of the log pore size distribution parameter β=ln(α). In addition, we also give solutions for the mean and variance of the unsaturated hydraulic conductivity. The analytical solutions of moment equations are validated via Monte Carlo simulations.     

Assessment of uncertainty associated with the estimation of well catchments by moment equations

Non-local stochastic moment equations are used successfully to analyze groundwater flow in randomly heterogeneous media. Here we present a moment equations-based approach to quantify the uncertainty associated with the estimation of well catchments. Our approach is based on the development of a complete second order formalism which allows obtaining the first statistical moments of the trajectories of conservative solute particles advected in a generally non-uniform groundwater flow. Approximate equations of moments of particles’ trajectories are then derived on the basis of a second order expansion in terms of the standard deviation of the aquifer log hydraulic conductivity. Analytical expressions are then obtained for the predictors of locations of mean stagnation points, together with their associated uncertainties. We implement our approach on heterogeneous media in bounded two-dimensional domains, with and without including the effect of conditioning on hydraulic conductivity information. The impact of domain size, boundary conditions, heterogeneity and non-stationarity of hydraulic conductivity on the prediction of a well catchment is explored. The results are compared against Monte Carlo simulations and semi-analytical solutions available in the literature. The methodology is applicable to both infinite and bounded domains and is free of distributional assumptions (and so applies to both Gaussian and non-Gaussian log hydraulic conductivity fields) and formally includes the effect of conditioning on available information.