A quasi three-dimensional model for flow and transport in unsaturated and saturated zones: 1. Implementation of the quasi two-dimensional case

A quasi three-dimensional (QUASI 3-D) model is presented for simulating the subsurface water flow and solute transport in the unsaturated and in the saturated zones of soil. The model is based on the assumptions of vertical flow in the unsaturated zone and essentially horizontal groundwater flow. The 1-D Richards equation for the unsaturated zone is coupled at the phreatic surface with the 2-D flow equation for the saturated zone. The latter was obtained by averaging 3-D flow equation in the saturated zone over the aquifer thickness. Unlike the Boussinesq equation for a leaky-phreatic aquifer, the developed model does not contain a storage term with specific yield and a source term for natural replenishment. Instead it includes a water flux term at the phreatic surface through which the Richards equation is linked with the groundwater flow equation. The vertical water flux in the saturated zone is evaluated on the basis of the fluid mass balance equation while the horizontal fluxes, in that equation, are prescribed by Darcy law. A 3-D transport equation is used to simulate the solute migration. A numerical algorithm to solve the problem for the general quasi 3-D case was developed. The developed methodology was exemplified for the quasi 2-D cross-sectional case (QUASI2D). Simulations for three synthetic problems demonstrate good agreement between the results obtained by QUASI2D and two fully 2-D flow and transport codes (SUTRA and 2DSOIL). Yet, simulations with the QUASI2D code were several times faster than those by the SUTRA and the 2DSOIL codes.     

Underdamped slug tests with unsaturated‐saturated flows by considering effects of wellbore skins

An analytical solution is presented for the slug tests conducted in a partially penetrating well in an unconfined aquifer affected from above by an unsaturated zone. The solution considers the effects of wellbore skin and oscillatory responses on underdamped slug tests. The flow in the saturated zone is described by a two‐dimensional, axially symmetric governing equation, and the flow in the unsaturated zone above the water table by a linearized one‐dimensional Richards' equation. The unsaturated medium properties are represented by the exponential constitutive relationships. A Laplace domain solution is derived using the Laplace and finite Fourier transform and the solution in the real‐time domain is evaluated using the numerical inverse Laplace transform method. The solution derived in this study is more general and reduces to the most commonly used solutions for slug tests in their specified conditions. It is found that the unsaturated flow has a significant impact on the slug test conducted in an unconfined aquifer. The impact of unsaturated flow on such a slug test is enhanced with a larger anisotropy ratio, a shorter well screen length, a shorter distance between the well screen and the water table, or a larger well screen radius. The impact of unsaturated flow on slug tests decreases as the degree of penetration (the length of well screen) increases. For a fixed well screen length, the impact of unsaturated flow on slug tests decreases as the distance between the centre of screen and the water table increases. A large dimensionless well screen radius (>0.01) leads to significant effects of unsaturated flow on slug tests. The unsaturated flow reduces the oscillatory responses to underdamped slug tests. The unsaturated zone has significant impact on slug test under high‐permeability wellbore skin.     

Environmental tracer transport (H and SF6) in the saturated and unsaturated zones and its use in nitrate pollution management

An important quantity in groundwater protection is the residence time of water in an aquifer. It relates to both the travel time of a pollutant to arrive at a well and the time span required for self-purification of a polluted aquifer after removal of pollutant inputs. Time scales for aquifers can be gained from artificial tracer experiments or from environmental tracer data, the latter offering the only realistic alternative if time scales of years or decades have to be taken into account.

Different tracers show different time scales due to their different transport mechanisms especially in the unsaturated zone. While solute tracers are moved advectively with the seepage water, gas tracers pass the unsaturated zone diffusively through the air phase. Depending on the properties of the unsaturated zone (hydraulic properties, thickness) this difference in behavior can be used to separate the subsurface transport process into the unsaturated and the saturated parts.

In a field study in Germany, SF₆ and ³H were used as environmental tracers. Both have a relatively well-known input function. Interpretation of data from observation wells by a box model approach led to spatially and temporally varying residence times. This was an indication that the influence of the unsaturated zone could not be neglected. While the gas tracer SF₆ shows only residence times in the saturated zone, the tracer ³H reflects the whole travel time of water including both the unsaturated and saturated zones. Using a one-dimensional plug-flow model for the unsaturated zone combined with a detailed two-dimensional flow and transport model for the saturated zone leads to a holistic and consistent interpretation of the measured tracer concentrations. The observed pattern of old water under thick loess cover and younger water under areas where the fractured basalt aquifer crops out is reproduced after adjusting only two parameters: the effective porosity of the saturated aquifer and the product of field capacity and thickness of the unsaturated zone. While the effective porosity of the saturated zone is adjusted by means of the SF₆ data, the field capacity of the loess layer is adjusted by means of the ³H observations. The thickness of the unsaturated zone is deduced from geological and pedological maps. All flow data are obtained from a calibrated flow model, which is based on geological data, observed heads and pumping tests only.

The transport model for the saturated zone was calibrated by fitting the porosity by means of gaseous tracer concentrations (SF₆). The combined saturated–unsaturated zone model was then calibrated by fitting the field capacity of the unsaturated zone by means of ³H concentrations. With this model it was possible to verify the observed NO₃ concentrations at the drinking water wells and to develop predictions for their future development under various scenarios of fertilizer input reduction in specific areas.     

On the leakage assumption applied to equations of groundwater flow

The behaviour and the delayed response of flow to a well in an unconfined aquifer is considered on the basis of recognition of the variable vertical movement of flow inherent in unconfined flow. The effect of elastic storage is taken into account.The typical S-shaped form of the characteristic time-drawdown curve of an unconfined aquifer is analysed. A procedure for determining the unconfined parameters is described.It is shown that the vertical components of flow are a significant factor in the delayed response of an unconfined aquifer while the contribution of the unsaturated zone is essentially subordinate. An equation including allowance for the drainage from the unsaturated zone is derived.     

Influence of Groundwater Hydraulic Gradient on Bank Storage Metrics

The hydraulic gradient between aquifers and rivers is one of the most variable properties in a river/aquifer system. Detailed process understanding of bank storage under hydraulic gradients is obtained from a two‐dimensional numerical model of a variably saturated aquifer slice perpendicular to a river. Exchange between the river and the aquifer occurs first at the interface with the unsaturated zone. The proportion of total water exchanged through the river bank compared to the river bed is a function of aquifer hydraulic conductivity, partial penetration, and hydraulic gradient. Total exchange may be estimated to within 50% using existing analytical solutions provided that unsaturated zone processes do not strongly influence exchange. Model‐calculated bank storage is at a maximum when no hydraulic gradient is present and increases as the hydraulic conductivity increases. However, in the presence of a hydraulic gradient, the largest exchange flux or distance of penetration does not necessarily correspond to the highest hydraulic conductivity, as high hydraulic conductivity increases the components of exchange both into and out of an aquifer. Flood wave characteristics do not influence ambient groundwater discharge, and so in large floods, hydraulic gradients must be high to reduce the volume of bank storage. Practical measurement of bank storage metrics is problematic due to the limitations of available measurement technologies and the nested processes of exchange that occur at the river‐aquifer interface. Proxies, such as time series concentration data in rivers and groundwater, require further development to be representative and quantitative.     

Subsurface transfer of chloride after a lake retreat in the central Andes

The area under study covers 3500 km2 in the upstream part of the closed catchment basin of the salt crust of Uyuni. This crust is a remnant of the saline Lake Tauca, which covered the area about 15,000 years ago. In the downstream part of the aquifer, the Cl concentration of ground water and Cl content in the unsaturated zone exceed 20 meq/L and 18 kg/m2, respectively. With the present hydrological conditions under semiarid conditions, the ground water residence time in the study area exceeds 3000 years. Transient simulations over 11,000 years were made using initial conditions as the retreat of Lake Tauca and taking into account a low recharge during the arid mid-Holocene period. The modeling simulates ground water flow, Cl transport, and ground water residence time. It includes the evaporation from the aquifer that leads to the accumulation of chloride in the unsaturated zone. Results of the modeling are consistent with the observations if it is assumed that the Cl previously accumulated in the unsaturated zone was flushed back into the aquifer around 2000 years B.P., contemporaneously with the end of the arid period.     

Monitoring of drawdown pattern during pumping in an unconfined physical aquifer model

In this study, we attempted to analyse a drawdown pattern around a pumping well in an unconfined sandy gravelly aquifer constructed in a laboratory tank by means of both experimental and numerical modelling of groundwater flow. The physical model consisted of recharge, aquifer and discharge zones. Permeability and specific yield of the aquifer material were determined by Dupuit approximation under steady‐state flow and stepwise gravitational drainage of groundwater, respectively. The drawdown of water table in pumping and neighbouring observation wells was monitored to investigate the effect of no‐flow boundary on the drawdown pattern during pumping for three different boundary conditions: (i) no recharge and no discharge with four no‐flow boundaries (Case 1); (ii) no recharge and reservoir with three no‐flow boundaries (Case 2); (iii) recharge and discharge with two no‐flow boundaries (Case 3). Based on the aquifer parameters, numerical modelling was also performed to compare the simulated drawdown with that observed. Results showed that a large difference existed between the simulated drawdown and that observed in wells for all cases. The reason for the difference could be explained by the formation of a curvilinear type water table between wells rather than a linear one due to a delayed response of water table in the capillary fringe. This phenomenon was also investigated from a mass balance study on the pumping volume. The curvilinear type of water table was further evidenced by measurement of water contents at several positions in the aquifer between wells using time domain reflectometry (TDR). This indicates that the existing groundwater flow model applicable to an unconfined aquifer lacks the capacity to describe a slow response of water table in the aquifer and care should be taken in the interpretation of water table formation in the aquifer during pumping. Copyright © 2001 John Wiley & Sons, Ltd.     

Groundwater recharge estimation from ephemeral streams. Case study: Wadi Tabalah,Saudi Arabia

Estimation of groundwater recharge to an unconfined aquifer is studied using analytical and numerical techniques and results are compared with field observations. There is an acute need for such estimation in water balance studies in arid climates, and the case study in this paper is for such a region. The wetting front movement in the unsaturated zone depends on antecedent soil moisture, the ponded water depth and its duration, and on the position of the water table and the hydraulic properties of the unsaturated zone. A hydraulic connection between the recharge basin and the aquifer is not immediately established because the wetting front is unsaturated. A numerical model is applied to estimate recharge in an arid-zone wadi, and its validity is tested by comparing it with an analytical solution of the equations. The calculated recharge values matched the piezometric levels observed at a well site at the edge of the wadi channel. The total recharge depths found by integration in the time domain provided a good estimate of the transmitted volume of water per unit length of wadi channel. The findings were confirmed by runoff volume measurements at gauging stations located in the basin. © 1997 John Wiley & Sons, Ltd.     

Unsaturated Zone Flow Processes and Aquifer Response Time in the Chalk Aquifer,Brighton, South East England

The Chalk aquifer is one of the main sources of water in South East England. The unsaturated zone in the aquifer plays an important role controlling the time and magnitude of recharge and is major pathway for contaminant transport to the water table. A range of previous work has addressed flow processes in the Chalk unsaturated zone, but physical understanding is still incomplete. Here we present the results of a study on flow mechanism in the Chalk unsaturated zone using a combination of statistical analysis and novel laboratory methods. The study was undertaken at three sites (North Heath Barn [NHB], Pyecombe East [PE], and Preston Park [PP]) on the Chalk of the Brighton block, South East England. Daily and hourly time series data of groundwater level and rainfall were correlated. The results show that a slower groundwater level response to rainfall occurs during dry seasons (summer and autumn) when the amount of effective rainfall is less than 4 mm/day, with a thicker and drier unsaturated zone. A faster response occurs during wet seasons (winter and spring) when the daily effective rainfall exceeds 4 mm/day with a thinner and wetter unsaturated zone. Periods of very rapid response (within 15 h) were observed during wet seasons at NHB and PE sites, with unsaturated hydraulic conductivity (K_u) inferred to reach 839 mm/day. A slower response was observed at an urbanized site (PP) as a result of reduction in direct recharge due to reduced infiltration, due to presences of impermeable infrastructure covering the area around PP borehole. Laboratory measurements of K_u of the Chalk matrix using a geotechnical centrifuge show variation from 4.27 to 0.07 mm/day, according to the level of saturation. Thus, the rapid responses cannot be linked to matrix flow only but indicate the contribution of fracture and karstic flow processes in conducting water.     

Analytical Solutions for Water Flow and Solute Transport in the Unsaturated Zone

Analytical solutions for the water flow and solute transport equations in the unsaturated zone are presented. We use the Broadbridge and White nonlinear model to solve the Richards’ equation for vertical flow under a constant infiltration rate. Then we extend the water flow solution and develop an exact parametric solution for the advection-dispersion equation. The method of characteristics is adopted to determine the location of a solute front in the unsaturated zone. The dispersion component is incorporated into the final solution using a singular perturbation method. The formulation of the analytical solutions is simple, and a complete solution is generated without resorting to computationally demanding numerical schemes. Indeed, the simple analytical solutions can be used as tools to verify the accuracy of numerical models of water flow and solute transport. Comparison with a finite-element numerical solution indicates that a good match for the predicted water content is achieved when the mesh grid is one-fourth the capillary length scale of the porous medium. However, when numerically solving the solute transport equation at this level of discretization, numerical dispersion and spatial oscillations were significant.     

Seepage face height, water table position, and well efficiency at steady state

When a fully penetrating well pumps an ideal unconfined aquifer at steady state, the water table usually does not join the water level in the well. There is a seepage face inside the well, which is a key element in evaluating the well performance. This problem is analyzed using the finite-element method, solving the complete equations for saturated and unsaturated flow. The seepage face position is found to be almost independent of the unsaturated zone properties. The numerical results are used to test the validity of several analytic approximations. Equations are proposed to predict the seepage face position at the pumping well for any well drawdown, and the water table position at any distance from the pumping well for any in-well drawdown. Practical hints are provided for installing monitoring wells and evaluating well efficiency.     

Analysis of the thickness of a fresh water lens and of the transition zone between this lens and upwelling saline water

In regions with saline groundwater, fresh water lenses may develop due to rainwater infiltration. The amount of fresh water that is available for e.g. agricultural crops depends on the thickness of the lens and the extent of mixing between fresh and saline water. In this paper, we consider the mixing of fresh water and upward moving saline ground water in low-lying deltaic areas. The parameters that dominate the flow and transport problem are investigated using dimensionless groups and scaled sensitivities. We characterize the numerically simulated thicknesses of the lens and of the mixing zone by spatial moments. Rayleigh number and mass flux ratio, which is the ratio of the salt water seepage and the precipitation, determine the thickness of the fresh water lens. The local thickness of the mixing zone is mainly influenced by the dispersive/diffusive groups and the mass flux ratio. In addition, convergence of streamlines towards an outflow boundary affects the thickness, particularly in the vicinity of this boundary.Analytical and numerical steady state solutions for lens thickness are compared, taking into account upward seepage, for the two cases with and without a density difference between lens and underlying groundwater. Agreement between the numerical and analytical solutions for the lens thickness is good except when the mass flux ratio becomes small. For zero mass flux ratio, it is implicitly assumed in the analytical solution that salt water is stagnant, and that is unrealistic. Relative contributions of longitudinal and transversal hydrodynamic dispersion and diffusion to the thickness of the mixing zone are quantified numerically for different phases of lens formation. Longitudinal dispersion dominates in the early stages of lens formation, while diffusion and transversal dispersion dominate at steady state.     

Hydraulics of a partially penetrating well with skin zone in a confined aquifer

A steady/quasi-steady model is developed for predicting flow into a partially penetrating well with skin zone in a confined aquifer overlying an impervious layer. The model takes into account flow through the bottom of the wellbore, finite skin thickness and finite horizontal and vertical extent of the aquifer. Moreover, the solution can be easily extended to include the mixed-type boundary condition at the well face, where a Dirichlet in the form of a specified hydraulic head and a Neumann in the form of zero flux coexist at the same time at different portions of the well face. The validity of the proposed solution is tested by comparing a few results obtained from the developed model with corresponding results obtained by analytical and numerical means. The study shows that, among other factors remaining constant, both the horizontal and vertical extent of an artesian aquifer, thickness of the skin zone, bottom flow and conductivity contrast of the skin and formation zones, play an important part in deciding flow to a well dug in the aquifer, and hence these factors must be considered while analyzing the problem. The model proposed here can be used to estimate skin thickness as well as hydraulic conductivities of the skin and formation zones of a well with skin zone in an artesian aquifer underlain by an impervious layer by utilizing pumping test data falling in the steady or quasi-steady state of a typical pumping test. As the proposed solution is of a general nature in the sense that it can handle, apart from partial penetration and bottom flow, the finite size skin zone and finite horizontal and vertical extent of an artesian aquifer together with the mixed-type boundary condition at the well face, it is hoped that the predictions coming out of the model will be more realistic than those obtained using solutions developed with more stringent assumptions.     

A Stable and Efficient Numerical Algorithm for Unconfined Aquifer Analysis

The nonlinearity of equations governing flow in unconfined aquifers poses challenges for numerical models, particularly in field-scale applications. Existing methods are often unstable, do not converge, or require extremely fine grids and small time steps. Standard modeling procedures such as automated model calibration and Monte Carlo uncertainty analysis typically require thousands of model runs. Stable and efficient model performance is essential to these analyses. We propose a new method that offers improvements in stability and efficiency and is relatively tolerant of coarse grids. It applies a strategy similar to that in the MODFLOW code to the solution of Richard's equation with a grid-dependent pressure/saturation relationship. The method imposes a contrast between horizontal and vertical permeability in gridblocks containing the water table, does not require "dry" cells to convert to inactive cells, and allows recharge to flow through relatively dry cells to the water table. We establish the accuracy of the method by comparison to an analytical solution for radial flow to a well in an unconfined aquifer with delayed yield. Using a suite of test problems, we demonstrate the efficiencies gained in speed and accuracy over two-phase simulations, and improved stability when compared to MODFLOW. The advantages for applications to transient unconfined aquifer analysis are clearly demonstrated by our examples. We also demonstrate applicability to mixed vadose zone/saturated zone applications, including transport, and find that the method shows great promise for these types of problem as well.     

Factors Governing Sustainable Groundwater Pumping near a River

The objective of this paper was to provide new insights into processes affecting riverbank filtration (RBF). We consider a system with an inflatable dam installed for enhancing water production from downstream collector wells. Using a numerical model, we investigate the impact of groundwater pumping and dam operation on the hydrodynamics in the aquifer and water production. We focus our study on two processes that potentially limit water production of an RBF system: the development of an unsaturated zone and riverbed clogging. We quantify river clogging by calibrating a time‐dependent riverbed permeability function based on knowledge of pumping rate, river stage, and temperature. The dynamics of the estimated riverbed permeability reflects clogging and scouring mechanisms. Our results indicate that (1) riverbed permeability is the dominant factor affecting infiltration needed for sustainable RBF production; (2) dam operation can influence pumping efficiency and prevent the development of an unsaturated zone beneath the riverbed only under conditions of sufficient riverbed permeability; (3) slow river velocity, caused by dam raising during summer months, may lead to sedimentation and deposition of fine‐grained material within the riverbed, which may clog the riverbed, limiting recharge to the collector wells and contributing to the development of an unsaturated zone beneath the riverbed; and (4) higher river flow velocities, caused by dam lowering during winter storms, scour the riverbed and thus increase its permeability. These insights can be used as the basis for developing sustainable water management of a RBF system.     

Small‐Scale ASR Between Flow Barriers in a Saline Aquifer

Regular aquifer storage recovery, ASR, is often not feasible for small‐scale storage in brackish or saline aquifers because fresh water floats to the top of the aquifer where it is unrecoverable. Flow barriers that partially penetrate a brackish or saline aquifer prevent a stored volume of fresh water from expanding sideways, thus increasing the recovery efficiency. In this paper, the groundwater flow and mixing is studied during injection, storage, and recovery of fresh water in a brackish or saline aquifer in a flow‐tank experiment and by numerical modeling to investigate the effect of density difference, hydraulic conductivity, pumping rate, cyclic operation, and flow barrier settings. Two injection and recovery methods are investigated: constant flux and constant head. Fresh water recovery rates on the order of 65% in the first cycle climbing to as much as 90% in the following cycles were achievable for the studied configurations with constant flux whereas the recovery efficiency was somewhat lower for constant head. The spatial variation in flow velocity over the width of the storage zone influences the recovery efficiency, because it induces leakage of fresh water underneath the barriers during injection and upconing of salt water during recovery.     

Aquifer response to stream-stage and recharge variations. II. Convolution method and applications

In this second of two papers, analytical step-response functions, developed in the companion paper for several cases of transient hydraulic interaction between a fully penetrating stream and a confined, leaky, or water-table aquifer, are used in the convolution integral to calculate aquifer heads, streambank seepage rates, and bank storage that occur in response to stream-stage fluctuations and basinwide recharge or evapotranspiration. Two computer programs developed on the basis of these step-response functions and the convolution integral are applied to the analysis of hydraulic interaction of two alluvial stream–aquifer systems in the northeastern and central United States. These applications demonstrate the utility of the analytical functions and computer programs for estimating aquifer and streambank hydraulic properties, recharge rates, streambank seepage rates, and bank storage. Analysis of the water-table aquifer adjacent to the Blackstone River in Massachusetts suggests that the very shallow depth of water table and associated thin unsaturated zone at the site cause the aquifer to behave like a confined aquifer (negligible specific yield). This finding is consistent with previous studies that have shown that the effective specific yield of an unconfined aquifer approaches zero when the capillary fringe, where sediment pores are saturated by tension, extends to land surface. Under this condition, the aquifer's response is determined by elastic storage only. Estimates of horizontal and vertical hydraulic conductivity, specific yield, specific storage, and recharge for a water-table aquifer adjacent to the Cedar River in eastern Iowa, determined by the use of analytical methods, are in close agreement with those estimated by use of a more complex, multilayer numerical model of the aquifer. Streambank leakance of the semipervious streambank materials also was estimated for the site. The streambank-leakance parameter may be considered to be a general (or lumped) parameter that accounts not only for the resistance of flow at the river–aquifer boundary, but also for the effects of partial penetration of the river and other near-stream flow phenomena not included in the theoretical development of the step-response functions.     

Numerical analysis of one-dimensional unsaturated flow in layered soils

The paper deals with numerical solutions to the Richards equation to simulate one-dimensional flow processes in the unsaturated zone of layered soil profiles. The equation is expressed in the pressure-based form and a finite-difference algorithm is developed for accurately estimating the values of the hydraulic conductivity between two neighboring nodes positioned in different soil layers, often referred to as the interlayer hydraulic conductivity. The algorithm is based upon flux conservation and continuity of pressure potential at the interface between two consecutive layers, and does not add significantly to simulation run time. The validity of the model is established for a number of test problems by comparing numerical results with the analytical solutions developed by Srivastava and Yeh²⁹ which hold for vertical infiltration towards the water table in a two-layer soil profile. The results show a significant reduction in relative mass balance errors when using the proposed model. Some specific insights into its numerical performance are also gained by comparisons with a numerical model in which the more common geometric averaging operator acts on the interlayer conductivities.     

Aquifer vulnerability to pesticide migration through till aquitards

This study investigates the influence of key factors-mainly recharge rate and degradation half-life--on downward migration of the widely used pesticide mecoprop (MCPP) through a typical clayey till aquitard. The study uses the numerical model FRAC3Dvs, which is a three-dimensional discrete fracture/matrix diffusion (DFMD) numerical transport model. The model was calibrated with laboratory and field data from a site near Havdrup, Denmark, but the overall findings are expected to be relevant to many other sites in similar settings. Fracture flow and MCPP transport parameters for the model were obtained through calibration using well-characterized laboratory experiments with large (0.5 m diameter by 0.5 m high) undisturbed columns of the fractured till and a field experiment. A second level of upscaling and sensitivity analysis was then carried out using data on hydraulic head, fracture spacing, and water budget from the field site. The simulations of downward migration of MCPP show that MCPP concentration and mass flux into the underlying aquifer, and hence the aquifer vulnerability to this pesticide compound, is mainly dependent on the degradation rate of the pesticide, the overall aquitard water budget, and the ground water recharge rate into the aquifer. The influence of flow rate, matrix diffusion, and degradation rate are intertwined. This results in one to four orders of magnitude higher MCPP flux into the aquifer from aquifer recharge rates of 20 and 120 mm/yr, respectively, for no degradation and MCPP half-life of 0.5 yr. From a sensitivity analysis it was found that the range of MCPP flux into the aquifer varied less than one order of magnitude due to (1) changing fracture spacing from 1 to 10 m, or (2) preferential flow along inclined thin sand layers, which represent common conditions for the current and other settings of clayey till in Denmark and other glaciated areas in Europe and North America. The results indicate that for aquifers overlain by fractured clayey tills, the vulnerability to contamination with pesticides (pesticide flux into the aquifer) and other widespread agricultural contaminants is going to vary strongly in the watershed as a function of the distribution of aquitard water budget (flow rate) and aquitard redox environment (controlling contaminant degradation rates), even if the thickness of the till is relatively constant. DFMD modeling of cause-effect relationships within such systems has great potential to support decisions in planning, regulation, and contaminant remediation.     

Numerical modelling of stream–aquifer interaction: Quantifying the impact of transient streambed permeability and aquifer heterogeneity

Stream–aquifer interaction plays a vital role in the water cycle, and a proper study of this interaction is needed for understanding groundwater recharge, contaminants migration, and for managing surface water and groundwater resources. A model‐based investigation of a field experiment in a riparian zone of the Schwarzbach river, a tributary of the Rhine River in Germany, was conducted to understand stream–aquifer interaction under alternative gaining and losing streamflow conditions. An equivalent streambed permeability, estimated by inverting aquifer responses to flood waves, shows that streambed permeability increased during infiltration of stream water to aquifer and decreased during exfiltration. Aquifer permeability realizations generated by multiple‐point geostatistics exhibit a high degree of heterogeneity and anisotropy. A coupled surface water groundwater flow model was developed incorporating the time‐varying streambed permeability and heterogeneous aquifer permeability realizations. The model was able to reproduce varying pressure heads at two observation wells near the stream over a period of 55 days. A Monte Carlo analysis was also carried out to simulate groundwater flow, its age distribution, and the release of a hypothetical wastewater plume into the aquifer from the stream. Results of this uncertainty analysis suggest (a) stream–aquifer exchange flux during the infiltration periods was constrained by aquifer permeability; (b) during exfiltration, this flux was constrained by the reduced streambed permeability; (c) the effect of temporally variable streambed permeability and aquifer heterogeneity were found important to improve the accurate capture of the uncertainty; and (d) probabilistic infiltration paths in the aquifer reveal that such pathways and the associated prediction of the extent of the contaminant plume are highly dependent on aquifer heterogeneity.