Transport In Heterogeneous Porous Media Research Articles

Krylov subspace approximations for the exponential Euler method: error estimates and the harmonic Ritz approximant

We study Krylov subspace methods for approximating the matrix-function vector product $\varphi(tA)b$ where $\varphi(z) = [\exp(z)-1]/z$. This product arises in the numerical integration of large stiff systems of differential equations by the Exponential Euler Method, where $A$~is the Jacobian matrix of the system. Recently, this method has found application in the simulation of transport phenomena in porous media within mathematical models of wood drying and groundwater flow. We develop an a posteriori upper bound on the Krylov subspace approximation error and provide a new interpretation of a previously published error estimate. This leads to an alternative Krylov approximation to $\varphi(tA)b$, the so-called Harmonic Ritz approximant, which we find does not exhibit oscillatory behaviour of the residual error. References E. J. Carr, T. J. Moroney and I. W. Turner. Efficient simulation of unsaturated flow using exponential time integration. Appl. Math. Comput., 217(14):6587--6596, 2011. doi:10.1016/j.amc.2011.01.041 E. J. Carr, I. W. Turner and P. Perre. A Jacobian-free exponential integrator for simulating transport in heterogeneous porous media: application to the drying of softwood, submitted for publication. E. J. Carr, I. W. Turner and P. Perre. A new control-volume finite-element scheme for heterogeneous porous media:application to the drying of softwood. Chem. Eng. Technol., 34(7):1143--1150, 2011. doi:10.1002/ceat.201100060 E. Celledoni and I. Moret. A Krylov projection method for systems of ODEs. Appl. Numer. Math., 24(2--3):365--378, 1997. doi:10.1016/S0168-9274(97)00033-0 N. J. Higham. Functions of matrices: theory and computation. SIAM, Philadelphia, PA, USA, 2008. M. Hochbruck, M. E. Hochstenbach. Subspace extraction for matrix functions, submitted for publication. http://www.win.tue.nl/ hochsten/pdf/funext.pdf M. Hochbruck, C. Lubich and H. Selhofer. Exponential integrators for large systems of differential equations. SIAM J. Sci. Comput., 19(5):1552--1574, 1998. doi:10.1137/S1064827595295337 B. V. Minchev and W. M. Wright. A review of exponential integrators for first order semi-linear problems. Numerics No. 2/05, Norwegian University of Science and Technology, Trondheim, Norway, 2005. http://www.math.ntnu.no/preprint/numerics/2005/N2-2005.ps P. Perre and I. Turner. A heterogeneous wood drying computational model that accounts for material property variation across growth rings. Chem. Eng. J., 86(1--2):117--131, 2002. doi:10.1016/S1385-8947(01)00270-4

Combining numerical simulations with time-domain random walk for pathogen risk assessment in groundwater

We present a methodology that combines numerical simulations of groundwater flow and advective transport in heterogeneous porous media with analytical retention models for computing the infection risk probability from pathogens in aquifers. The methodology is based on the analytical results presented in [1,2] for utilising the colloid filtration theory in a time-domain random walk framework. It is shown that in uniform flow, the results from the numerical simulations of advection yield comparable results as the analytical TDRW model for generating advection segments. It is shown that spatial variability of the attachment rate may be significant, however, it appears to affect risk in a different manner depending on if the flow is uniform or radially converging. In spite of the fact that numerous issues remain open regarding pathogen transport in aquifers on the field scale, the methodology presented here may be useful for screening purposes, and may also serve as a basis for future studies that would include greater complexity.

Direct Multiple-Point Geostatistical Simulation of Edge Properties for Modeling Thin Irregularly Shaped Surfaces

Thin, irregularly shaped surfaces such as clay drapes often have a major control on flow and transport in heterogeneous porous media. Clay drapes are often complex, curvilinear three-dimensional surfaces and display a very complex spatial distribution. Variogram-based stochastic approaches are also often not able to describe the spatial distribution of clay drapes since complex, curvilinear, continuous, and interconnected structures cannot be characterized using only two-point statistics. Multiple-point geostatistics aims to overcome the limitations of the variogram. The premise of multiple-point geostatistics is to move beyond two-point correlations between variables and to obtain (cross) correlation moments at three or more locations at a time using training images to characterize the patterns of geological heterogeneity. Multiple-point geostatistics can reproduce thin irregularly shaped surfaces such as clay drapes, but this is often computationally very intensive. This paper describes and applies a methodology to simulate thin, irregularly shaped surfaces with a smaller CPU and RAM demand than the conventional multiple-point statistical methods. The proposed method uses edge properties for indicating the presence of thin irregularly shaped surfaces. Instead of pixel values, edge properties indicating the presence of irregularly shaped surfaces are simulated using snesim. This method allows direct simulation of edge properties instead of pixel properties to make it possible to perform multiple-point geostatistical simulations with a larger cell size and thus a smaller computation time and memory demand. This method is particularly valuable for three-dimensional applications of multiple-point geostatistics.

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

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

PDF equations for advective–reactive transport in heterogeneous porous media with uncertain properties

We consider advective–reactive solute transport in porous media whose hydraulic and transport properties are uncertain. These properties are treated as random fields, which renders nonlinear advection–reaction transport equations stochastic. We derive a deterministic equation for the probability density function (PDF) of the concentration of a solute that undergoes heterogeneous reactions, e.g., precipitation or dissolution. The derivation treats exactly (without linearization) a reactive term in the transport equation which accounts for uncertainty (randomness) in both flow velocity and kinetic rate constants but requires a closure, such as a Large-Eddy-Diffusivity (LED) approximation used in the present analysis. No closure is required when reaction rates are the only source of uncertainty. We use exact concentration PDFs obtained for this setting to analyze the accuracy of our general, LED-based PDF equations.

Macroscale lattice‐Boltzmann methods for low Peclet number solute and heat transport in heterogeneous porous media

This paper introduces new methods for simulating subsurface solute and heat transport in heterogeneous media using large‐scale lattice‐Boltzmann models capable of representing both macroscopically averaged porous media and open channel flows. Previous examples of macroscopically averaged lattice‐Boltzmann models for solute and heat transport are only applicable to homogeneous media. Here, we extend these models to properly account for heterogeneous pore‐space distributions. For simplicity, in the majority of this paper we assume low Peclet number flows with an isotropic dispersion tensor. Nevertheless, this approach may also be extended to include anisotropic‐dispersion by using multiple relaxation time lattice‐Boltzmann methods. We describe two methods for introducing heterogeneity into macroscopically averaged lattice‐Boltzmann models. The first model delivers the desired behavior by introducing an additional time‐derivative term to the collision rule; the second model by separately weighting symmetric and anti‐symmetric components of the fluid packet densities. Chapman‐Enskog expansions are conducted on the governing equations of the two models, demonstrating that the correct constitutive behavior is obtained in both cases. In addition, methods for improving model stability at low porosities are also discussed: (1) an implicit formulation of the model; and (2) a local transformation that normalizes the lattice‐Boltzmann model by the local porosity. The model performances are evaluated through comparisons of simulated results with analytical solutions for one‐ and two‐dimensional flows, and by comparing model predictions to finite element simulations of advection isotropic‐dispersion in heterogeneous porous media. We conclude by presenting an example application, demonstrating the ability of the new models to couple with simulations of reactive flow and changing flow geometry: a simulation of groundwater flow through a carbonate system.

Analytical solution for one-dimensional advection–dispersion transport equation with distance-dependent coefficients

Summary Mathematical models describing contaminant transport in heterogeneous porous media are often formulated as an advection–dispersion transport equation with distance-dependent transport coefficients. In this work, a general analytical solution is presented for the linear, one-dimensional advection–dispersion equation with distance-dependent coefficients. An integrating factor is employed to obtain a transport equation that has a self-adjoint differential operator, and a solution is found using the generalized integral transform technique (GITT). It is demonstrated that an analytical expression for the integrating factor exists for several transport equation formulations of practical importance in groundwater transport modeling. Unlike nearly all solutions available in the literature, the current solution is developed for a finite spatial domain. As an illustration, solutions for the particular case of a linearly increasing dispersivity are developed in detail and results are compared with solutions from the literature. Among other applications, the current analytical solution will be particularly useful for testing or benchmarking numerical transport codes because of the incorporation of a finite spatial domain.

Significance of higher moments for complete characterization of the travel time probability density function in heterogeneous porous media using the maximum entropy principle

The travel time formulation of advective transport in heterogeneous porous media is of interest both conceptually, e.g., for incorporating retention processes, and in applications where typically the travel time peak, early, and late arrivals of contaminants are of major concern in a regulatory or remediation context. Furthermore, the travel time moments are of interest for quantifying uncertainty in advective transport of tracers released from point sources in heterogeneous aquifers. In view of this interest, the travel time distribution has been studied in the literature; however, the link to the hydraulic conductivity statistics has been typically restricted to the first two moments. Here we investigate the influence of higher travel time moments on the travel time probability density function (pdf) in heterogeneous porous media combining Monte Carlo simulations with the maximum entropy principle. The Monte Carlo experimental pdf is obtained by the adaptive Fup Monte Carlo method (AFMCM) for advective transport characterized by a multi‐Gaussian structure with exponential covariance considering two injection modes (in‐flux and resident) and lnK variance up to 8. A maximum entropy (MaxEnt) algorithm based on Fup basis functions is used for the complete characterization of the travel time pdf. All travel time moments become linear with distance. Initial nonlinearity is found mainly for the resident injection mode, which exhibits a strong nonlinearity within first 5IY for high heterogeneity. For the resident injection mode, the form of variance and all higher moments changes from the familiar concave form predicted by the first‐order theory to a convex form; for the in‐flux mode, linearity is preserved even for high heterogeneity. The number of moments sufficient for a complete characterization of the travel time pdf mainly depends on the heterogeneity level. Mean and variance completely describe travel time pdf for low and mild heterogeneity, skewness is dominant for lnK variance around 4, while kurtosis and fifth moment are required for lnK variance higher than 4. Including skewness seems sufficient for describing the peak and late arrivals. Linearity of travel time moments enables the prediction of asymptotic behavior of the travel time pdf which in the limit converges to a symmetric distribution and Fickian transport. However, higher‐order travel time moments may be important for most practical purposes and in particular for advective transport in highly heterogeneous porous media for a long distance from the source.

Uncertainty quantification in modeling flow and transport in porous media

This special issue of Stochastic Environmental Research and Risk Assessment in honor of Shlomo P. Neuman presents recent developments in uncertainty quantification in predictions of flow and transport in heterogeneous porous media. It covers subjects as diverse as stochastic hydrogeology, aquifer characterization, effective (average) descriptions of subsurface processes, inverse modeling, and scaling theory, all of which benefited from the insightful contributions of Shlomo P. Neuman. We hope that this volume will become a reference point for students learning these subject, scientists pursuing innovative research, and professionals and water resources planners looking to quantify predicting uncertainty arising from the intrinsic complexity of the subsurface environment. S. P. Neuman is one of the pioneers in the use of stochastic modeling as a means of dealing with uncertainty in hydraulic and transport parameters of heterogeneous subsurface environments. This special issue contains a number of contributions in this field. Riva et al. introduce a stochastic approach to interpret field-scale tracer tests conducted in heterogeneous aquifers. By using Monte Carlo simulations, the authors demonstrate that directly linking porosity and permeability distribution to the spatial variability of soil particle sizes renders the best prediction of the observed heavy tailing of measured breakthrough curves. Zhang et al. compare the performance of Monte Carlo and stochastic collocation methods to evaluate the probability of the concentration of a conservative solute attaining a given threshold. Their analysis shows that Latin Hypercube sampling and Quasi Monte Carlo outperform standard Monte Carlo simulations, and that a sparse-grid collocation and a probabilistic collocation method provide an accurate and efficient alternative to Monte Carlo simulations. Broyda et al. develop partial differential equations satisfied by the probability density function (PDF) of the concentration of solutes undergoing heterogeneous reactions in porous media with uncertain chemical properties. The authors obtain a semi-analytical solution for the concentration PDF and demonstrate that its shape and evolution are significantly affected by both spatial variability of advective velocity and the Damkohler number. Tartakovsky presents a novel Lagrangian particle model based on smoothed particle hydrodynamics (SPH). His analysis reveals the advantages of Lagrangian methods for the stochastic analysis of miscible density-driven fluid flows in the presence of RayleighTaylor instability. When used in the context of Monte Carlo simulations, the method can provide reliable estimated of the probability of occurrence of rare events. Two papers are concerned with the estimation of hydrogeological parameters governing groundwater flow. Liu et al. present an approach to efficient and accurate parameter estimation for nonlinear flow and transport problems in porous media. They discuss the relative merits and weaknesses of different schemes, including Maximum a posteriori (MAP) and Monte Carlo based methods, used to parameterize a DNAPL dissolution/transport model. Methods to diagnose Markov Chain Monte Carlo (MCMC) samples are introduced, compared and discussed in the context of computationally demanding parameter estimation A. Guadagnini (&) Dipartimento di Ingeneria Idraulca, Ambientale, Infrastrutture Viarie e Rilevamento, Politecnico di Milano, Piazza L. Da Vinci, 32, 20133 Milan, Italy e-mail: alberto.guadagnini@polimi.it

A modified Lagrangian-volumes method to simulate nonlinearly and kinetically sorbing solute transport in heterogeneous porous media

Transport in subsurface environments is conditioned by physical and chemical processes in interaction, with advection and dispersion being the most common physical processes and sorption the most common chemical reaction. Existing numerical approaches become time-consuming in highly-heterogeneous porous media. In this paper, we discuss a new efficient Lagrangian method for advection-dominated transport conditions. Modified from the active-walker approach, this method comprises dividing the aqueous phase into elementary volumes moving with the flow and interacting with the solid phase. Avoiding numerical diffusion, the method remains efficient whatever the velocity field by adapting the elementary volume transit times to the local velocity so that mesh cells are crossed in a single numerical time step. The method is flexible since a decoupling of the physical and chemical processes at the elementary volume scale, i.e. at the lowest scale considered, is achieved. We implement and validate the approach to the specific case of the nonlinear Freundlich kinetic sorption. The method is relevant as long as the kinetic sorption-induced spreading remains much larger than the dispersion-induced spreading. The variability of the surface-to-volume ratio, a key parameter in sorption reactions, is explicitly accounted for by deforming the shape of the elementary volumes.

A multidimensional streamline-based method to simulate reactive solute transport in heterogeneous porous media

We present a new streamline-based numerical method for simulating reactive solute transport in porous media. The key innovation of the method is that both longitudinal and transverse dispersion are incorporated accurately without numerical dispersion. Dispersion is approximated in a flow-oriented grid using a combination of a one-dimensional finite difference scheme and a meshless approximation. In contrast to previous hybrid alternatives to incorporate dispersion in streamline-based simulations, the proposed scheme does not require a grid and, hence, it does not introduce numerical dispersion. In addition, the proposed scheme eliminates numerical oscillations and negative concentration values even when the dispersion tensor includes the off-diagonal coefficients and the flow field is non-uniform. We demonstrate that for a set of two- and three-dimensional benchmark problems, the new proposed streamline-based formulation compares favorably to two state of the art finite volume and hybrid Eulerian–Lagrangian solvers.

An exponential integrator for advection-dominated reactive transport in heterogeneous porous media

We present an exponential time integrator in conjunction with a finite volume discretisation in space for simulating transport by advection and diffusion including chemical reactions in highly heterogeneous porous media representative of geological reservoirs. These numerical integrators are based on the variation of constants solution and solving the linear system exactly. This is at the expense of computing the exponential of the stiff matrix comprising the finite volume discretisation. Using real Léja points or a Krylov subspace technique compared to standard finite difference-based time integrators. We observe for a variety of example applications that numerical solutions with exponential methods are generally more accurate and require less computational cost. They hence comprise an efficient and accurate method for simulating non-linear advection-dominated transport in geological formations.

Lattice Boltzmann Models for Flow and Transport in Saturated Karst

Flow and transport simulation in karst aquifers remains a significant challenge for the ground water modeling community. Darcy's law-based models cannot simulate the inertial flows characteristic of many karst aquifers. Eddies in these flows can strongly affect solute transport. The simple two-region conduit/matrix paradigm is inadequate for many purposes because it considers only a capacitance rather than a physical domain. Relatively new lattice Boltzmann methods (LBMs) are capable of solving inertial flows and associated solute transport in geometrically complex domains involving karst conduits and heterogeneous matrix rock. LBMs for flow and transport in heterogeneous porous media, which are needed to make the models applicable to large-scale problems, are still under development. Here we explore aspects of these future LBMs, present simple examples illustrating some of the processes that can be simulated, and compare the results with available analytical solutions. Simulations are contrived to mimic simple capacitance-based two-region models involving conduit (mobile) and matrix (immobile) regions and are compared against the analytical solution. There is a high correlation between LBM simulations and the analytical solution for two different mobile region fractions. In more realistic conduit/matrix simulation, the breakthrough curve showed classic features and the two-region model fit slightly better than the advection-dispersion equation (ADE). An LBM-based anisotropic dispersion solver is applied to simulate breakthrough curves from a heterogeneous porous medium, which fit the ADE solution. Finally, breakthrough from a karst-like system consisting of a conduit with inertial regime flow in a heterogeneous aquifer is compared with the advection-dispersion and two-region analytical solutions.

A meshless method to simulate solute transport in heterogeneous porous media

We derive a meshless numerical method based on smoothed particle hydrodynamics (SPH) for the simulation of conservative solute transport in heterogeneous geological formations. We demonstrate that the new proposed scheme is stable, accurate, and conserves global mass. We evaluate the performance of the proposed method versus other popular numerical methods for the simulation of one- and two-dimensional dispersion and two-dimensional advective–dispersive solute transport in heterogeneous porous media under different Pèclet numbers. The results of those benchmarks demonstrate that the proposed scheme has important advantages over other standard methods because of its natural ability to control numerical dispersion and other numerical artifacts. More importantly, while the numerical dispersion affecting traditional numerical methods creates artificial mixing and dilution, the new scheme provides numerical solutions that are “physically correct”, greatly reducing these artifacts.

Spatial Markov processes for modeling Lagrangian particle dynamics in heterogeneous porous media

We investigate the representation of Lagrangian velocities in heterogeneous porous media as Markov processes. We use numerical simulations to show that classical descriptions of particle velocities using Markov processes in time fail because low velocities are much more strongly correlated in time than high velocities. We demonstrate that Lagrangian velocities describe a Markov process at fixed distances along the particle trajectories (i.e., a spatial Markov process). This remarkable property has significant implications for modeling effective transport in heterogeneous velocity fields: (i) the spatial Lagrangian velocity transition densities are sufficient to fully characterize these complex velocity field organizations, (ii) classical effective transport descriptions that rely on Markov processes in time for the particle velocities are not suited for describing transport in heterogeneous porous media, and (iii) an alternative effective transport description derives from the Markovian nature of the spatial velocity transitions. It expresses particle movements as a random walk in space time characterized by a correlated random temporal increment and thus generalizes the continuous time random walk model to transport in correlated velocity fields.

Numerical modeling of the flow and transport of radionuclides in heterogeneous porous media

This paper is concerned with numerical methods for the modeling of flow and transport of contaminant in porous media. The numerical methods feature the mixed finite element method over triangles as a solver to the Darcy flow equation and a conservative finite volume scheme for the concentration equation. The convective term is approximated with a Godunov scheme over the dual finite volume mesh, whereas the diffusion–dispersion term is discretized by piecewise linear conforming triangular finite elements. It is shown that the scheme satisfies a discrete maximum principle. Numerical examples demonstrate the effectiveness of the methodology for a coupled system that includes an elliptic equation and a diffusion–convection–reaction equation arising when modeling flow and transport in heterogeneous porous media. The proposed scheme is robust, conservative, efficient, and stable, as confirmed by numerical simulations.

Special issue on multiscale methods for flow and transport in heterogeneous porous media

Special issue on multiscale methods for flow and transport in heterogeneous porous media

Upscaling models of solute transport in porous media through genetic programming

Due to the considerable computational demands of modeling solute transport in heterogeneous porous media, there is a need for upscaled models that do not require explicit resolution of the small-scale heterogeneity. This study investigates the development of upscaled solute transport models using genetic programming (GP), a domain-independent modeling tool that searches the space of mathematical equations for one or more equations that describe a set of training data. An upscaling methodology is developed that facilitates both the GP search and the implementation of the resulting models. A case study is performed that demonstrates this methodology by developing vertically averaged equations of solute transport in perfectly stratified aquifers. The solute flux models developed for the case study were analyzed for parsimony and physical meaning, resulting in an upscaled model of the enhanced spreading of the solute plume, due to aquifer heterogeneity, as a process that changes from predominantly advective to Fickian. This case study not only demonstrates the use and efficacy of GP as a tool for developing upscaled solute transport models, but it also provides insight into how to approach more realistic multi-dimensional problems with this methodology.

Stochastic Methods in Subsurface Contaminant Hydrology

Stochastic methods for flow and solute transport in heterogeneous porous media have become well received by the hydrology community in the last two to three decades. While several books on stochastic analysis for flow problems have been published in the last few years, this is the first book that

Two-domain description of solute transport in heterogeneous porous media: Comparison between theoretical predictions and numerical experiments

This paper deals with two-equation models describing solute transport in highly heterogeneous porous systems and more particularly dual permeability structures composed of high- and low-permeability regions. A macroscopic two-equation model has been previously proposed in the literature based on the volume averaging technique [Ahmadi A, Quintard M, Whitaker S. Transport in chemically and mechanically heterogeneous porous media V: two-equation model for solute transport with adsorption, Adv Water Resour 1998;22:59–86; Cherblanc F, Ahmadi A, Quintard M. Two-medium description of dispersion in heterogeneous porous media: calculation of macroscopic properties. Water Resour Res 2003;39(6):1154–73]. Through this theoretical upscaling method, both convection and dispersion mechanisms are taken into account in both regions, allowing one to deal with a large range of heterogeneous systems. In this paper, the numerical tools associated with this model are developed in order to test the theory by comparing macroscopic concentration fields to those obtained by Darcy-scale numerical experiments. The heterogeneous structures considered are made up of low-permeability nodules embedded in a continuous high-permeability region. Several permeability ratios are used, leading to very different macroscopic behaviours. Taking advantage of the Darcy-scale simulations, the role of convection and dispersion in the mass exchange between the two regions is investigated. Large-scale averaged concentration fields and elution curves are extracted from the Darcy-scale numerical experiments and compared to the theoretical predictions given by the two-equation model. Very good agreement is found between experimental and theoretical results. A permeability ratio around 100 presents a behaviour characteristic of “mobile–mobile” systems emphasizing the relevance of this two-equation description. Eventually, the theory is used to set-up a criterion for the existence of local equilibrium conditions. The potential importance of local-scale dispersion in reducing large-scale dispersion is highlighted. The results also confirm that a non-equilibrium description may be necessary in such systems, even if local-equilibrium behaviour could be observed.