Numerical implementation of a backward probabilistic model of ground water contamination

Backward location and travel time probabilities can be used to characterize known and unknown sources or prior positions of ground water contamination. Backward location probability describes the position of the observed contamination at some time in the past; backward travel time probability describes the amount of time prior to observation that the contamination was released from its source or was at a particular upgradient location. The governing equation for backward probabilities is the adjoint of the governing equation for contaminant transport, but with new load terms. Numerical codes that have been written to solve the forward equations of contaminant transport, e.g., the advection-dispersion equation, can also be used to solve the adjoint equation for location and travel time probabilities; however, the interpretation of the results is different and some new approximations must be made for the load terms. We present the governing equations for backward location and travel time probabilities, and provide appropriate numerical approximations for these load terms using the cell-centered finite difference method, one of the most popular numerical methods in ground water hydrology. We discuss some additional numerical considerations for the backward model including boundary conditions, reversal of the flow field, and interpretation of the results. We illustrate the implementation of the backward probability model using hypothetical examples in one- and two-dimensional domains. We also present a three-dimensional application of a pump-and-treat remediation capture zone delineation at the Massachusetts Military Reservation. The illustrations are performed using MODFLOW-96 for flow simulations and MT3DMS for transport simulations.     

Identification of contaminant point source in surface waters based on backward location probability density function method

A backward location probability density function (BL-PDF) method capable of identifying location of point sources in surface waters is presented in this paper. The relation of forward location probability density function (FL-PDF) and backward location probability density, based on adjoint analysis, is validated using depth-averaged free-surface flow and mass transport models and several surface water test cases. The solutions of the backward location PDF transport equation agreed well to the forward location PDF computed using the pollutant concentration at the monitoring points. Using this relation and the distribution of the concentration detected at the monitoring points, an effective point source identification method is established. The numerical error of the backward location PDF simulation is found to be sensitive to the irregularity of the computational meshes, diffusivity, and velocity gradients. The performance of identification method is evaluated regarding the random error and number of observed values. In addition to hypothetical cases, a real case was studied to identify the source location where a dye tracer was instantaneously injected into a stream. The study indicated the proposed source identification method is effective, robust, and quite efficient in surface waters; the number of advection–diffusion equations needed to solve is equal to the number of observations.     

On the relations between the hydrological dynamical systems of water budget,travel time,response time and tracer concentrations

Using previous results on extended Petri Nets (EPN), we present the relations between various hydrological dynamical systems (HDSys) derived from the water budget. Once the water budget has been implemented, there is a consistent way of getting the equations for backward travel time distributions, for forward response time distributions and for the concentration of a solute or tracer. We show that the water budget has a correspondence of one to many with the backward travel time distributions. In fact, to any one of the water budget equations there correspond as many equations as there are input precipitation events. The backward travel time distributions are related to the response time distributions by Niemi's relationship and, if there are n outputs, by the definition of n − 1 partition functions. These determine what fraction of the water volume injected into the control volume at a specific time t_in goes asymptotically to a specific output. Given the backward travel time distributions, the output concentration of a solute or tracer also depends on the input concentration. The paper clarifies the complicated relations described above by taking [Hydrology and Earth System Sciences, 20, 299–328] as an example from literature. Once the appropriate information is presented, it is shown how these HDSys can be solved simultaneously without duplicating calculations. Then, it is also shown that, under the hypothesis of uniform mixing of water ages within each reservoir, these systems can be solved exactly.     

Modeling Transport in Transient Ground-Water Flow: An Unacknowledged Approximation

Abstract. During unsteady or transient ground-water flow, the fluid mass per unit volume of aquifer changes as the potentiometric head changes, and solute transport is affected by this change in fluid storage. Three widely applied numerical models of two-dimensional transport partially account for the effects of transient flow by removing terms corresponding to the fluid continuity equation from the transport equation, resulting in a simpler governing equation. However, fluid-storage terms remaining in the transport equation that change during transient flow are, in certain cases, held constant in time in these models. For the case of increasing heads, this approximation, which is unacknowledged in these models'documentation, leads to transport velocities that are too high, and increased concentration at fluid and solute sources. If heads are dropping in time, computed transport velocities are too low. Using parameters that somewhat exaggerate the effects of this approximation, an example numerical simulation indicates solute travel time error of about 14 percent but only minor errors due to incorrect dilution volume. For horizontal flow and transport models that assume fluid density is constant, the product of porosity and aquifer thickness changes in time: initial porosity times initial thickness plus the change in head times the storage coefficient. This formula reduces to the saturated thickness in unconfined aquifers if porosity is assumed to be constant and equal to specific yield. The computational cost of this more accurate representation is insignificant and is easily incorporated in numerical models of solute transport.     

Watershed structural influences on the distributions of stream network water and solute travel times under baseflow conditions

Watershed structure influences the timing, magnitude, and spatial location of water and solute entry to stream networks. In turn, stream reach transport velocities and stream network geometry (travel distances) further influence the timing of export from watersheds. Here, we examine how watershed and stream network organization can affect travel times of water from delivery to the stream network to arrival at the watershed outlet. We analysed watershed structure and network geometry and quantified the relationship between stream discharge and solute velocity across six study watersheds (11.4 to 62.8 km²) located in the Sawtooth Mountains of central Idaho, USA. Based on these analyses, we developed stream network travel time functions for each watershed. We found that watershed structure, stream network geometry, and the variable magnitude of inputs across the network can have a pronounced affect on water travel distances and velocities within a stream network. Accordingly, a sample taken at the watershed outlet is composed of water and solutes sourced from across the watershed that experienced a range of travel times in the stream network. We suggest that understanding and quantifying stream network travel time distributions are valuable for deconvolving signals observed at watershed outlets into their spatial and temporal sources, and separating terrestrial and in‐channel hydrological, biogeochemical, and ecological influences on in‐stream observations. Copyright © 2016 John Wiley & Sons, Ltd.     

An analytical solution for ground water transit time through unconfined aquifers

An exact, closed-form analytical solution is developed for calculating ground water transit times within Dupuit-type flow systems. The solution applies to steady-state, saturated flow through an unconfined, horizontal aquifer recharged by surface infiltration and discharging to a downgradient fixed-head boundary. The upgradient boundary can represent, using the same equation, a no-flow boundary or a fixed head. The approach is unique for calculating travel times because it makes no a priori assumptions regarding the limit of the water table rise with respect to the minimum saturated aquifer thickness. The computed travel times are verified against a numerical model, and examples are provided, which show that the predicted travel times can be on the order of nine times longer relative to existing analytical solutions.     

Laboratory sandbox validation of pollutant source location methods

Inverse methods can be used to recover the pollutant source location from concentration data. In this paper, the relative effectiveness of two proposed methods, simultaneous release function and source location identification (SRSI) and backward probability model based on adjoint state method (BPM-ASM) are evaluated using real data collected by using experimental equipment. The device is a sandbox that reproduces an unconfined aquifer in which all the variables are controlled. A numerical model was calibrated using experimental observations. The SRSI is a stochastic procedure which finds the source location and the release history by means of a Bayesian geostatistical approach (GA). The BPM-ASM provides the backward probability location of the pollutant detected at a monitoring point by means of a reverse transport simulation. The results show that both methods perform well. While the simultaneous release function and SRSI method requires a preliminary delineation of a probable source area and some weak hypotheses about the statistical structure of the unknown release function, the backward probability model requires some hypothesis about the contaminant release time. A case study was performed using two observation points only, and despite the scarcity of data, both methodologies were able to accurately reconstruct the true source location. The GA has the advantage to recover the release history function too, whilst the backward probability model works well with fewer data. If there are many observations, both methodologies may be computationally heavy. A transfer function approach has been adopted for the numerical definition of the sensitivity matrix in the SRSI method. The reliability of the experimental equipment was tested in previous laboratory works, conducted under several different conditions.     

Analytical solutions to the transient,unsaturated transport of water and contaminants through horizontal porous media

We present a range of analytical solutions to the combined transient water and solute transport for horizontal flow. We adopt the concept of a scale and time dependent dispersivity used for contaminant transport in aquifers and apply it to transient, unsaturated horizontal flow to develop similarity solutions for both constant solute concentration and solute flux boundary conditions. Through the use of a specific form of the water profile as used by Brutsaert [Water Resour Res 1968:4;785], the solute profiles can be reduced to a simple quadrature. We also derive a solution for the instantaneous injection of water and solute into a horizontal media for an arbitrary dispersivity. It is found that the solute concentration remains constant in both space and time as the water redistributes, suggesting that the solute does not disperse relative to the water.     

Single-Well Injection-Drift Test to Estimate Groundwater Velocity

A simple algebraic equation is presented here to estimate the magnitude of groundwater velocity based on data from a single-well injection-drift test thereby eliminating the time-consuming and costly extraction phase. A volume of tracer-amended water was injected by forced-gradient into a single well followed by monitoring of the conservative solute tracers under natural-gradient conditions as their upgradient portions drifted back through the well. The breakthrough curve data from the single well during the drift phase was analyzed to determine the mean travel times of the tracers. The estimated mean upgradient travel distance back through the single well and the mean travel times of the tracers were used in a simple algebraic equation to estimate groundwater velocity. The groundwater velocity based on the single-well injection-drift test was estimated to be approximately 0.64 ft per day. Two transects of observation wells were used to monitor the natural-gradient tracer transport downgradient of the injection well. The one-dimensional, or dual-well, transport of the tracer from the injection well to the nearest downgradient observation well indicated that the groundwater velocity was 0.55 ft per day. The two-dimensional, or multi-well, transport of the center of mass of the tracers indicated that the groundwater velocity was 0.60 ft per day; the dual- and multi-well results were in excellent agreement with those from the single-well and validated the simple algebraic equation. The new single-well method presented here is relatively simple, rapid, and does not require an extraction phase.     

Direct simulation of solute recycling in irrigated areas

Solute recycling from irrigation can be described as the process that occurs when the salt load that is extracted from irrigation wells and distributed on the fields is returned to the groundwater below irrigated surfaces by deep percolation. Unless the salt load leaves the system by means of drains or surface runoff, transfer to the groundwater will take place, sooner or later. This can lead to solute accumulation and thus to groundwater degradation, particularly in areas where extraction rates exceed infiltration rates (semi-arid and arid regions). Thus, considerable errors can occur in a predictive solute mass budget if the recycling process is not accounted for in the calculation. A method is proposed which allows direct simulation of solute recycling. The transient solute response at an extraction well is shown to be a superposition of solute mass flux contributions from n recycling cycles and is described as a function of the travel time distribution between a recycling point and a well. This leads to an expression for a transient ‘recycling source’ term in the advection–dispersion equation, which generates the effect of solute recycling. At long times, the ‘recycling source’ is a function of the local capture probability of the irrigation well and the solute mass flux captured by the well from the boundaries. The predicted concentration distribution at steady state reflects the maximum spatial concentration distribution in response to solute recycling and can thus be considered as the solute recycling potential or vulnerability of the entire domain for a given hydraulic setting and exploitation scheme. Simulation of the solute recycling potential is computationally undemanding and can therefore, for instance, be used for optimisation purposes. Also, the proposed method allows transient simulation of solute recycling with any standard flow and transport code.     

Well vulnerability: a quantitative approach for source water protection

The concept of vulnerability of drinking water sources is reviewed, and a quantitative approach for assessing well vulnerability for complex three-dimensional ground water systems is developed. The approach focuses on the relative expected impact of potential contaminant sources at unknown locations within a well capture zone, providing relative measures of intrinsic well vulnerability, including the expected times of arrival of a contaminant, the dispersion-related reduction in concentration, the time taken to breach a certain quality objective, and the corresponding exposure times. Thus, the result of the analysis includes the usual advective travel time information used in conventional wellhead protection analysis, plus a set of selected quantitative measures expressing the expected impact. The technique is based on adjoint theory and combines forward- and backward-in-time transport modeling using a standard numerical flow and transport code. The methodology is demonstrated using the case study of a complex glacial multiaquifer system in Ontario. The new approach will be useful in helping water managers develop more physically based and quantitative wellhead protection strategies.     

The response of playa and sabkha hydraulics and mineralogy to climate forcing

Dry playa lakes and sabkhat often represent the terminus of large ground water flow systems and act as integrators of both upgradient (recharge) and downgradient discharge (evaporation). Ground water levels beneath playa/sabkha systems show a variety of surprising responses driven by large evaporation demands and chemical processes not typically encountered in more humid regions. When the water table is very close to the land surface, almost instantaneous rises can be observed with little observed change in either upgradient ground water recharge or potential evaporation. Conversely, when water tables are several meters below the playa surface, water table responses to interannual variability of recharge can be damped and lag significantly behind such changes. This review of the dynamics of shallow water tables in playa lakes and sabkhat discusses the pertinent hydraulic and solute processes and extracts a simple but comprehensive model based on soil physics for predicting the water table response to either upstream recharge changes or changes in potential evaporation at the playa/sabkha. Solutes and associated authigenic minerals are also shown to be important in discriminating both the causes and effects of water level fluctuations.     

Identification of Multiple Sources of Groundwater Contamination by Dual Isotopes

Chlorinated solvents are one of the most commonly detected groundwater contaminants in industrial areas. Identification of polluters and allocation of contaminant sources are important concerns in the evaluation of complex subsurface contamination with multiple sources. In recent years, compound‐specific isotope analyses (CSIA) have been employed to discriminate among different contaminant sources and to better understand the fate of contaminants in field‐site studies. In this study, the usefulness of dual isotopes (carbon and chlorine) was shown in assessments of groundwater contamination at an industrial complex in Wonju, Korea, where groundwater contamination with chlorinated solvents such as trichloroethene (TCE) and carbon tetrachloride (CT) was observed. In November 2009, the detected TCE concentrations at the study site ranged between nondetected and 10,066 µg/L, and the CT concentrations ranged between nondetected and 985 µg/L. In the upgradient area, TCE and CT metabolites were detected, whereas only TCE metabolites were detected in the downgradient area. The study revealed the presence of separate small but concentrated TCE pockets in the downgradient area, suggesting the possibility of multiple contaminant sources that created multiple comingling plumes. Furthermore, the variation of the isotopic (δ¹³C and δ³⁷Cl) TCE values between the upgradient and downgradient areas lends support to the idea of multiple contamination sources even in the presence of detectable biodegradation. This case study found it useful to apply a spatial distribution of contaminants coupled with their dual isotopic values for evaluation of the contaminated sites and identification of the presence of multiple sources in the study area.     

Implementation of Solute Transport in the Vadose Zone into the “HYDRUS Package for MODFLOW”

The “HYDRUS package for MODFLOW” is an existing MODFLOW package that allows MODFLOW to simultaneously evaluate transient water flow in both unsaturated and saturated zones. The package is based on incorporating parts of the HYDRUS-1D model (to simulate unsaturated water flow in the vadose zone) into MODFLOW (to simulate saturated groundwater flow). The coupled model is effective in addressing spatially variable saturated-unsaturated hydrological processes at the regional scale. However, one of the major limitations of this coupled model is that it does not have the capability to simulate solute transport along with water flow and therefore, the model cannot be employed for evaluating groundwater contamination. In this work, a modified unsaturated flow and transport package (modified HYDRUS package for MODFLOW and MT3DMS) has been developed and linked to the three-dimensional (3D) groundwater flow model MODFLOW and the 3D groundwater solute transport model MT3DMS. The new package can simulate, in addition to water flow in the vadose zone, also solute transport involving many biogeochemical processes and reactions, including first-order degradation, volatilization, linear or nonlinear sorption, one-site kinetic sorption, two-site sorption, and two-kinetic sites sorption. Due to complex interactions at the groundwater table, certain modifications of the pressure head (compared to the original coupling) and solute concentration profiles were incorporated into the modified HYDRUS package. The performance of the newly developed model is evaluated using HYDRUS (2D/3D), and the results indicate that the new model is effective in simulating the movement of water and contaminants in the saturated-unsaturated flow domains.     

Mineral Precipitation Upgradient from a Zero-Valent Iron Permeable Reactive Barrier

Core samples taken from a zero-valent iron permeable reactive barrier (ZVI PRB) at Cornhusker Army Ammunition Plant, Nebraska, were analyzed for physical and chemical characteristics. Precipitates containing iron and sulfide were present at much higher concentrations in native aquifer materials just upgradient of the PRB than in the PRB itself. Sulfur mass balance on core solids coupled with trends in ground water sulfate concentrations indicates that the average ground water flow after 20 months of PRB operation was approximately twenty fold less than the regional ground water velocity. Transport and reaction modeling of the aquifer PRB interface suggests that, at the calculated velocity, both iron and hydrogen could diffuse upgradient against ground water flow and thereby contribute to precipitation in the native aquifer materials. The initial hydraulic conductivity ( K ) of the native materials is less than that of the PRB and, given the observed precipitation in the upgradient native materials, it is likely that K reduction occurred upgradient to rather than within the PRB. Although not directly implicated, guar gum used during installation of the PRB is believed to have played a role in the precipitation and flow reduction processes by enhancing microbial activity.     

Characterizing heterogeneous soil water flow and solute transport using information measures

Heterogeneous water flow and solute transport in soils are an important phenomenon and difficult to be characterized. The objectives of this study were to investigate the heterogeneity of solute transport related to heterogeneous soil water flow using dye infiltration experiments, and to characterize heterogeneous water flow and solute transport in soils using the information theory. Field experiments of dye infiltration were performed in four plots. Various information measures were applied to characterize information content and complexity of water flow and solute transport in soils. Information contents and complexities of the maximum and apparent infiltration depths, and the mean and standard deviation of concentrations in the vertical direction of the plots were calculated. More heterogeneous processes of soil water flow and transport result in higher information/complexity values. The probability distributions of mean concentration were similar to those of the corresponding apparent infiltration depths for the plots, indicating that heterogeneity of dye concentrations was closely related to that of soil water flow. However, the range of information entropy and complexity of the water flow sequences was much narrower than that of the sequences of the concentrations. The results suggested that the transport processes were more heterogeneous than the water flow processes. Compared with the probability distributions of flow parameters, the information measures appeared to be a more versatile tool to describe flow and transport heterogeneities in soils.     

Evaluating the reliability of time series analysis to estimate variable riparian travel times by numerical groundwater modelling

The transition zones between rivers and adjacent riparian aquifers are locations of high biogeochemical activities that contribute to a removal of potentially hazardous substances in the aquatic system. The potential of the removal processes depends highly on subsurface water travel times, which can be determined by using the propagation of electrical conductivity (EC) signal from the river into the riparian aquifer. Although this method has been applied and verified in many studies, we observe possible limitations for the usage of EC fluctuation analysis. Our findings are based on EC time series analyses during storm events and artificial hydropeaks induced by watermill operations. Travel times derived by cross‐correlation analysis were compared with travel times calculated based on backward particle tracking of a calibrated transient numerical groundwater flow model. The cross‐correlation method produced only reasonable travel times for the artificial hydropeaks. In contrast, cross‐correlation analysis of the EC data during natural storm events resulted in implausibly negative or unrealistically low travel times for the bulk of the data sets. We conclude that the reason for this behaviour is, first, the low EC contrast between river and groundwater in connection with a strong damping of the infiltrating river EC signal into the subsurface during storm events. Second, the existence of old and less‐mineralized riparian water between the river and the monitoring well resulted in bank‐storage‐driven EC breakthrough curves with earlier arrival times and the subsequent estimation of implausible riparian travel times.     

Specification of Matrix Cleanup Goals in Fractured Porous Media

Semianalytical transient solutions have been developed to evaluate what level of fractured porous media (e.g., bedrock or clay) matrix cleanup must be achieved in order to achieve compliance of fracture pore water concentrations within a specified time at specified locations of interest. The developed mathematical solutions account for forward and backward diffusion in a fractured porous medium where the initial condition comprises a spatially uniform, nonzero matrix concentration throughout the domain. Illustrative simulations incorporating the properties of mudstone fractured bedrock demonstrate that the time required to reach a desired fracture pore water concentration is a function of the distance between the point of compliance and the upgradient face of the domain where clean groundwater is inflowing. Shorter distances correspond to reduced times required to reach compliance, implying that shorter treatment zones will respond more favorably to remediation than longer treatment zones in which back‐diffusion dominates the fracture pore water response. For a specified matrix cleanup goal, compliance of fracture pore water concentrations will be reached sooner for decreased fracture spacing, increased fracture aperture, higher matrix fraction organic carbon, lower matrix porosity, shorter aqueous phase decay half‐life, and a higher hydraulic gradient. The parameters dominating the response of the system can be measured using standard field and laboratory techniques.     

One-dimensional uniform and time varying solute dispersion along transient groundwater flow in a semi-infinite aquifer

An analytical solution for the space-time variation of contaminant concentration in one-dimensional transient groundwater flow in a homogenous semi-infinite aquifer, subjected to time-dependent source contamination, is derived. The uniform and time varying dispersion along transient groundwater flow is investigated under two conditions. First, the flow velocity distribution in the aquifer is considered as a sinusoidally varying function, and second, the flow velocity distribution is treated as an exponentially increasing function of time. It is assumed that initially the aquifer is not solute free, so the initial background concentration is considered as an exponentially decreasing function of the space variable which is tending to zero at infinity. It is assumed that dispersion is directly proportional to the square of the velocity, noting that experimental observations indicate that dispersion is directly proportional to the velocity with a power ranging from 1 to 2. The analytical solution is illustrated using an example and may help benchmark numerical codes and solutions.     

Transport with spatially variable kinetic sorption: recursion formulation

A recursion formulation for the transport of linearly sorbing solutes undergoing non-equilibrium sorption is developed. Constant or spatially varying sorption kinetics can be modeled using the recursion approach. The sorption and desorption rates are modeled as two independent random processes with a prescribed mean and covariance structure with spatial variability in the rate parameters included as well. The recursion solution, in terms of the probability density function for solute travel times, is derived by specifying transition probabilities for moving between the aqueous and sorbed phases. A few simple examples are used to illustrate the approach. The computer implementation leads to a very rapid algorithm that is easily extended to cover cases beyond the basic model presented here.