Coupled groundwater flow and transport: 2. Thermohaline and 3D convection systems

This work continues the analysis of variable density flow in groundwater systems. It focuses on both thermohaline (double-diffusive) and three-dimensional (3D) buoyancy-driven convection processes. The finite-element method is utilized to tackle these complex non-linear problems in two and three dimensions. The preferred numerical approaches are discussed regarding appropriate basic formulations, balance-consistent discretization techniques for derivative quantitites, extension of the Boussinesq approximation, proper constraint conditions, time marching schemes, and computational strategies for solving large systems. Applications are presented for the thermohaline Elder and salt dome problem as well as for the 3D extension of the Elder problem with and without thermohaline effects and a 3D Bénard convection process. The simulations are performed by using the package FEFLOW. Conclusions are drawn with respect to numerical efforts and the appropriateness for practical needs.     

Coupled inverse modelling of groundwater flow and mass transport and the worth of concentration data

This paper presents the extension of the self-calibrating method to the coupled inverse modelling of groundwater flow and mass transport. The method generates equally likely solutions to the inverse problem that display the variability as observed in the field and are not affected by a linearisation of the state equations. Conditioning to the state variables is measured by an objective function including, among others, the mismatch between the simulated and measured concentrations. Conditioning is achieved by minimising the objective function by gradient-based methods. The gradient contains the partial derivatives of the objective function with respect to: log conductivities, log storativities, prescribed heads at boundaries, retardation coefficients and mass sources. The derivatives of the objective function with respect to log conductivity are the most cumbersome and need the most CPU-time to be evaluated. For this reason, to compute this derivative only advective transport is considered. The gradient is calculated by the adjoint-state method. The method is demonstrated in a controlled, synthetic study, in which the worth of concentration data is analysed. It is shown that concentration data are essential to improve transport predictions and also help to improve aquifer characterisation and flow predictions, especially in the upstream part of the aquifer, even in the case that a considerable amount of other experimental data like conductivities and heads are available. Besides, conditioning to concentration data reduces the ensemble variances of estimated transmissivity, hydraulic head and concentration.     

Data-based modelling approach for variable density flow and solute transport simulation in a coastal aquifer

Data-based models, namely artificial neural network (ANN), support vector machine (SVM), genetic programming (GP) and extreme learning machine (ELM), were developed to approximate three-dimensional, density-dependent flow and transport processes in a coastal aquifer. A simulation model, SEAWAT, was used to generate data required for the training and testing of the data-based models. Statistical analysis of the simulation results obtained by the four models show that the data-based models could simulate the complex salt water intrusion process successfully. The selected models were also compared based on their computational ability, and the results show that the ELM is the fastest technique, taking just 0.5 s to simulate the dataset; however, the SVM is the most accurate, with a Nash-Sutcliffe efficiency (NSE) ≥ 0.95 and correlation coefficient R ≥ 0.92 for all the wells. The root mean square error (RMSE) for the SVM is also significantly less, ranging from 12.28 to 77.61 mg/L.     

Simulation of groundwater flow and mass transport under uncertainty

An analytical tool which facilitates the analysis of stochastic groundwater problems without resorting to Monte Carlo analysis can be developed using an approach analogous to that used in the field of applied mathematics to investigate the propagation of waves in random media. The fundamental assumption of this approach is that the random part of the porous medium property is small compared to the non-random part. This assumption allows us to treat the problem of solving a partial differential equation with stochastic coefficients using perturbation theory.     

Multi-resolution adaptive modeling of groundwater flow and transport problems

Many groundwater flow and transport problems, especially those with sharp fronts, narrow transition zones, layers and fingers, require extensive computational resources. In this paper, we present a novel multi-resolution adaptive Fup approach to solve the above mentioned problems. Our numerical procedure is the Adaptive Fup Collocation Method (AFCM), based on Fup basis functions and designed through a method of lines (MOL). Fup basis functions are localized and infinitely differentiable functions with compact support and are related to more standard choices such as splines or wavelets. This method enables the adaptive multi-resolution approach to solve problems with different spatial and temporal scales with a desired level of accuracy using the entire family of Fup basis functions. In addition, the utilized collocation algorithm enables the mesh free approach with consistent velocity approximation and flux continuity due to properties of the Fup basis functions. The introduced numerical procedure was tested and verified by a few characteristic groundwater flow and transport problems, the Buckley–Leverett multiphase flow problem, the 1-D vertical density driven problem and the standard 2-D seawater intrusion benchmark–Henry problem. The results demonstrate that the method is robust and efficient particularly when describing sharp fronts and narrow transition zones changing in space and time.     

Scaling of flow and transport behavior in heterogeneous groundwater systems

Three-dimensional numerical simulations using a detailed synthetic hydraulic conductivity field developed from geological considerations provide insight into the scaling of subsurface flow and transport processes. Flow and advective transport in the highly resolved heterogeneous field were modeled using massively parallel computers, providing a realistic baseline for evaluation of the impacts of parameter scaling. Upscaling of hydraulic conductivity was performed at a variety of scales using a flexible power law averaging technique. A series of tests were performed to determine the effects of varying the scaling exponent on a number of metrics of flow and transport behavior. Flow and transport simulation on high-performance computers and three-dimensional scientific visualization combine to form a powerful tool for gaining insight into the behavior of complex heterogeneous systems.Many quantitative groundwater models utilize upscaled hydraulic conductivity parameters, either implicitly or explicitly. These parameters are designed to reproduce the bulk flow characteristics at the grid or field scale while not requiring detailed quantification of local-scale conductivity variations. An example from applied groundwater modeling is the common practice of calibrating grid-scale model hydraulic conductivity or transmissivity parameters so as to approximate observed hydraulic head and boundary flux values. Such parameterizations, perhaps with a bulk dispersivity imposed, are then sometimes used to predict transport of reactive or non-reactive solutes. However, this work demonstrates that those parameters that lead to the best upscaling for hydraulic conductivity and head do not necessarily correspond to the best upscaling for prediction of a variety of transport behaviors. This result reflects the fact that transport is strongly impacted by the existence and connectedness of extreme-valued hydraulic conductivities, in contrast to bulk flow which depends more strongly on mean values. It provides motivation for continued research into upscaling methods for transport that directly address advection in heterogeneous porous media.An electronic version of this article is available online at the journal's homepage at http://www.elsevier.nl/locate/advwatres or http://www.elsevier.com/locate/advwatres (see “Special section on vizualization”. The online version contains additional supporting information, graphics, and a 3D animation of simulated particle movement.©1998 Elsevier Science Limited. All rights reserved     

Stochastic modelling of groundwater flow and solute transport in aquifers

This paper presents an introductory overview of recently developed stochastic theories for tackling spatial variability problems in predicting groundwater flow and solute transport. Advantages and limitations of the theories are discussed. Lastly, strategies based on the stochastic approaches to predict solute transport in aquifers are recommended.     

A variable source area for groundwater evapotranspiration: impacts on modeling stream flow

Evapotranspiration (ET) plays a crucial role in catchment water budgets, typically accounting for more than 50% of annual precipitation falling within temperate deciduous forests. Groundwater ET is a portion of total ET that occurs where plant roots extend to the capillary fringe above the phreatic surface or induce upward movement of water from the water table by hydraulic redistribution. Groundwater ET is spatially restricted to riparian zones or other areas where the groundwater is accessible to plants. Due to the difficulty in measuring groundwater ET, it is rarely incorporated explicitly into hydrological models. In this study, we calibrated Topographic Model (TOPMODEL) using a 14‐year hydrograph record and added a groundwater ET pathway to derive a new model, Groundwater Evapotranspiration TOPMODEL (GETTOP). We inspected groundwater elevations and stream flow hydrographs for evidence of groundwater ET, examined the relationship between groundwater ET and topography, and delineated the area where groundwater ET is likely to take place. The total groundwater ET flux was estimated using a hydrological model. Groundwater ET was larger where the topography was flat and the groundwater table was shallow, occurring within about 10% of the area in a headwater catchment and accounting for 6 to 18% of total annual ET. The addition of groundwater ET to GETTOP improved the simulation of stream discharge and more closely balanced the watershed water budget. Copyright © 2013 John Wiley & Sons, Ltd.     

Mathematical models of groundwater flow with a horizontal drain

An analytical review is given to drainage theory, the foundation of which was laid down by researchers in our country in the first half of the 20th century. Later the author of this paper developed and thoroughly studied a complex of mathematical models describing two-dimensional groundwater flow in drained soils based on the boundary problems of the theory of analytical functions. These models were used to provide a hydrodynamic substantiation for the regularities in subsoil water dynamics that were established by ameliorators in relation to the leaching of salts from soils; to identify the specific features of groundwater flow with drainage in the presence in the flow of one or two unknown free boundaries; and to develop approaches to the analysis of such flows in a direct formulation.     

Stochastic perturbation analysis of groundwater flow. Spatially variable soils,semi-infinite domains and large fluctuations

As is well known, a complete stochastic solution of the stochastic differential equation governing saturated groundwater flow leads to an infinite hierarchy of equations in terms of higher-order moments. Perturbation techniques are commonly used to close this hierarchy, using power-series expansions. These methods are applied by truncating the series after a finite number of terms, and products of random gradients of conductivity and head potential are neglected. Uncertainty regarding the number or terms required to yield a sufficiently accurate result is a significant drawback with the application of power series-based perturbation methods for such problems. Low-order series truncation may be incapable of representing fundamental characteristics of flow and can lead to physically unreasonable and inaccurate solutions of the stochastic flow equation. To support this argument, one-dimensional, steady-state, saturated groundwater flow is examined, for the case of a spatially distributed hydraulic conductivity field. An ordinary power-series perturbation method is used to approximate the mean head, using second-order statistics to characterize the conductivity field. Then an interactive perturbation approach is introduced, which yields improved results compared to low-order, power-series perturbation methods for situations where strong interactions exist between terms in such approximations. The interactive perturbation concept is further developed using Feynman-type diagrams and graph theory, which reduce the original stochastic flow problem to a closed set of equations for the mean and the covariance functions. Both theoretical and practical advantages of diagrammatic solutions are discussed; these include the study of bounded domains and large fluctuations.     

Impact of log-transmissivity variogram structure on groundwater flow and transport predictions

We analyze the impact of the choice of the variogram model adopted to characterize the spatial variability of natural log-transmissivity on the evaluation of leading (statistical) moments of hydraulic heads and contaminant travel times and trajectories within mildly (randomly) heterogeneous two-dimensional porous systems. The study is motivated by the fact that in several practical situations the differences between various variogram types and a typical noisy sample variogram are small enough to suggest that one would often have a hard time deciding which of the tested models provides the best fit. Likewise, choosing amongst a set of seemingly likely variogram models estimated by means of geostatistical inverse models of flow equations can be difficult due to lack of sensitivity of available model discrimination criteria. We tackle the problem within the framework of numerical Monte Carlo simulations for mean uniform and radial flow scenarios. The effect of three commonly used isotropic variogram models, i.e., Gaussian, Exponential and Spherical, is analyzed. Our analysis clearly shows that (ensemble) mean values of the quantities of interest are not considerably influenced by the variogram shape for the range of parameters examined. Contrariwise, prediction variances of the quantities examined are significantly affected by the choice of the variogram model of the log-transmissivity field. The spatial distribution of the largest/lowest values of the relative differences observed amongst the tested models depends on a combination of variogram shape and parameters and relative distance from internal sources and the outer domain boundary. Our findings suggest the need of developing robust techniques to discriminate amongst a set of seemingly equally likely alternative variogram models in order to provide reliable uncertainty estimates of state variables.     

Soil pipe flow tracer experiments: 1. Connectivity and transport characteristics

Much debate has occurred in catchment hydrology regarding the connectivity of flow paths from upslope areas to catchment outlets. This study was conducted in two catchments, one with three upper branches, in a loess soil with a fragipan that fosters lateral flow and exhibits an extensive distribution of soil pipe collapse features. The study aimed to determine the connectivity of multiple soil pipe networks as well as determine pipe flow velocities during storm events. Fluorescein dye was injected directly into soil pipes at the upper most pipe collapse feature of four different hillslopes. Breakthrough curves (BTC) were determined by sampling multiple pipe collapse features downslope. The BTCs were used to determine the ‘average’ (centre of mass) and ‘maximum’ (first arrival) flow velocities. This study confirmed that these catchments contain individual continuous soil pipe networks that extend over 190 m and connect the upper most hillslopes areas with the catchment outlet. While the flow paths are continuous, the individual pipe networks consist of alternating reaches of subsurface flow through soil pipes and reaches of surface flow through gullies formed by pipe collapses. In addition, flow can be occurring both through the subsurface soil pipes simultaneous with surface flow generated by artesian flow from the soil pipes. The pipe flow velocities were as high as 0.3 m/s, which was in the range of streamflow velocities. These pipe flow velocities were also in the range of velocities observed in pinhole erosion tests suggesting that these large, mature soil pipes are still actively eroding. Copyright © 2015 John Wiley & Sons, Ltd.     

Parameter-independent model reduction of transient groundwater flow models: Application to inverse problems

A new methodology is proposed for the development of parameter-independent reduced models for transient groundwater flow models. The model reduction technique is based on Galerkin projection of a highly discretized model onto a subspace spanned by a small number of optimally chosen basis functions. We propose two greedy algorithms that iteratively select optimal parameter sets and snapshot times between the parameter space and the time domain in order to generate snapshots. The snapshots are used to build the Galerkin projection matrix, which covers the entire parameter space in the full model. We then apply the reduced subspace model to solve two inverse problems: a deterministic inverse problem and a Bayesian inverse problem with a Markov Chain Monte Carlo (MCMC) method. The proposed methodology is validated with a conceptual one-dimensional groundwater flow model. We then apply the methodology to a basin-scale, conceptual aquifer in the Oristano plain of Sardinia, Italy. Using the methodology, the full model governed by 29,197 ordinary differential equations is reduced by two to three orders of magnitude, resulting in a drastic reduction in computational requirements.     

Alternative split-operator approach for solving chemical reaction/groundwater transport models

Various schemes are available to solve coupled transport/reaction mathematical models, one of the most efficient and easy to apply being the two-step split-operator method in which the transport and reaction steps are performed separately. Operator splitting, however, does not solve exactly the fully coupled numerical model derived from the governing partial differential and algebraic equations describing the transport and reaction processes. An error, proportional to Δt (the time step used in the numerical solution) is introduced. Thus, small time steps must be used to ensure that accurate solutions result. An alternative scheme is presented, which iterates to the exact solution of the fully coupled numerical model. The new scheme enables accurate solutions to be calculated more efficiently than the two-step method, while maintaining separation of the transport and reaction steps in the calculations. As in the two-step method, the reaction calculations are performed node-wise throughout the computation grid. However, because the scheme relies on LU factorisation of the coefficient matrix in the transport equation solution, the reaction calculations must be performed in sequence, the sequence order being determined by the ordering of the nodes in the grid. Also, because LU factorisation is used, the scheme is limited to solute transport problems for which LU factorisation is a practical solution method.     

Cluster-based fuzzy models for groundwater flow in the unsaturated zone

In this paper fuzzy models are used as an alternative to describe groundwater flow in the unsaturated zone. The core of these models consists of a fuzzy rule-based model of the Takagi–Sugeno type. Various fuzzy clustering algorithms are compared in the data-driven identification of these Takagi–Sugeno models. The performance of the resulting fuzzy models is evaluated on the training surface on which they were identified, and on time series measurements of water content values obtained through an experiment carried out by the non-vegetated terrain (NVT) workgroup of the European Microwave Signature Laboratory (EMSL) (see [Mancini M, Hoeben R, Troch PA. Multifrequency radar observations of bare surface soil moisture content: a laboratory experiment. Water Resour Res 1999;35(6):1827–38] and [Hoeben R, Troch PA. Assimilation of active microwave observation data for soil moisture profile estimation. Water Resour Res 2000;36(10):2805–19]). Despite higher errors at the borders of high water content values in the training surface, good results are obtained on the simulation of the time series.     

Bias and uncertainty in regression-calibrated models of groundwater flow in heterogeneous media

Groundwater models need to account for detailed but generally unknown spatial variability (heterogeneity) of the hydrogeologic model inputs. To address this problem we replace the large, m-dimensional stochastic vector β that reflects both small and large scales of heterogeneity in the inputs by a lumped or smoothed m-dimensional approximation γθ_∗, where γ is an interpolation matrix and θ_∗ is a stochastic vector of parameters. Vector θ_∗ has small enough dimension to allow its estimation with the available data. The consequence of the replacement is that model function f(γθ_∗) written in terms of the approximate inputs is in error with respect to the same model function written in terms of β, f(β), which is assumed to be nearly exact. The difference f(β) − f(γθ_∗), termed model error, is spatially correlated, generates prediction biases, and causes standard confidence and prediction intervals to be too small. Model error is accounted for in the weighted nonlinear regression methodology developed to estimate θ_∗ and assess model uncertainties by incorporating the second-moment matrix of the model errors into the weight matrix. Techniques developed by statisticians to analyze classical nonlinear regression methods are extended to analyze the revised method. The analysis develops analytical expressions for bias terms reflecting the interaction of model nonlinearity and model error, for correction factors needed to adjust the sizes of confidence and prediction intervals for this interaction, and for correction factors needed to adjust the sizes of confidence and prediction intervals for possible use of a diagonal weight matrix in place of the correct one. If terms expressing the degree of intrinsic nonlinearity for f(β) and f(γθ_∗) are small, then most of the biases are small and the correction factors are reduced in magnitude. Biases, correction factors, and confidence and prediction intervals were obtained for a test problem for which model error is large to test robustness of the methodology. Numerical results conform with the theoretical analysis.     

Coupled groundwater flow and heat transport simulation for estimating transient aquifer–stream exchange at the lowland River Spree (Germany)

Subsurface flow and heat transport near Freienbrink, NE Germany, was simulated in order to study groundwater–surface water exchange between a floodplains aquifer and a section of the lowland River Spree and an adjacent oxbow. Groundwater exfiltration was the dominant process, and only fast surface water level rises resulted in temporary infiltration into the aquifer. The main groundwater flow paths are identified based on a 3D groundwater flow model. To estimate mass fluxes across the aquifer–surface water interfaces, a 2D flow and heat transport modelling approach along a transect of 12 piezometers was performed. Results of steady‐state and transient water level simulations show an overall high accuracy with a Spearman coefficient ρ = 0.9996 and root mean square error (RMSE) = 0.008 m. Based on small groundwater flow velocities of about 10^?7 to 10^?6 ms^?1, mean groundwater exfiltration rates of 233 l m^?2 d^?1 are calculated. Short periods of surface water infiltration into the aquifer do not exceed 10 days, and the infiltration rates are in the same range. The heat transport was modelled with slightly less accuracy (ρ = 0.8359 and RMSE = 0.34 °C). In contrast to the predominant groundwater exfiltration, surface water temperatures determine the calculated temperatures in the upper aquifer below both surface water bodies down to 10 m during the whole simulation period. These findings emphasize prevailing of heat conduction over advection in the upper aquifer zones, which seems to be typical for lowland streams with sandy aquifer materials and low hydraulic gradients. Moreover, this study shows the potential of coupled numerical flow and heat transport modelling to understand groundwater–surface water exchange processes in detail. Copyright © 2013 John Wiley & Sons, Ltd.     

Variable-density groundwater flow and solute transport in irregular 2D fracture networks

Numerical simulations of variable-density flow and solute transport have been conducted to investigate dense plume migration for various configurations of 2D fracture networks. For orthogonal fractures, simulations demonstrate that dispersive mixing in fractures with small aperture does not stabilize vertical plume migration in fractures with large aperture. Simulations in non-orthogonal 2D fracture networks indicate that convection cells form and that they overlap both the porous matrix and fractures. Thus, transport rates in convection cells depend on matrix and fracture flow properties. A series of simulations in statistically equivalent networks of fractures with irregular orientation show that the migration of a dense plume is highly sensitive to the geometry of the network. If fractures in a random network are connected equidistantly to the solute source, few equidistantly distributed fractures favor density-driven transport. On the other hand, numerous fractures have a stabilizing effect, especially if diffusive transport rates are high. A sensitivity analysis for a network with few equidistantly distributed fractures shows that low fracture aperture, low matrix permeability and high matrix porosity impede density-driven transport because these parameters reduce groundwater flow velocities in both the matrix and the fractures. Enhanced molecular diffusion slows down density-driven transport because it favors solute diffusion from the fractures into the low-permeability porous matrix where groundwater velocities are smaller. For the configurations tested, variable-density flow and solute transport are most sensitive to the permeability and porosity of the matrix, which are properties that can be determined more accurately than the geometry and hydraulic properties of the fracture network, which have a smaller impact on density-driven transport.