Multi-resolution adaptive modeling of groundwater flow and transport problems

doi:10.1016/j.advwatres.2006.10.007

Advances in Water Resources

Volume 30, Issue 5, May 2007, Pages 1105-1126

https://doi.org/10.1016/j.advwatres.2006.10.007 Get rights and content

Abstract

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.

Introduction

Many groundwater flow and transport problems exhibit a wide range of space and/or temporal scales characterized by sharp gradients resulting in fingering and layering with the existence of sharp interface and narrow transition zones. These characteristics are commonly present in problems of unsaturated and multiphase flow [5], [13], [23], [34], density driven flow and transport [6], [17], [43] and reactive transport [3], [39].

The numerical modeling of such processes usually presents significant difficulties in resolving numerical oscillations and dispersion. In order to overcome these difficulties demanding computational resources with a very fine grid and small time steps are needed. In recent numerical approaches adaptive methods with low computational costs are being developed. The first attempt to apply them has been focused upon using classical finite difference and finite element methods [2], [17]. The main difficulty in applying these methods is finding a stable solution at the transition between zones having different discretization. Significant improvements have been obtained by the adaptive discontinuous finite element method, e.g. [5], [34]. Recently, there have been many attempts to develop new adaptive procedures which, among others, are focused upon using the adaptive wavelet Galerkin methods [9], [12] and collocation methods [7], [8], [15], [25], [26], [46], [47]. The adaptive wavelet Galerkin methods have three potential difficulties: treatment of general boundary conditions, treatment of nonlinearities and solving problems with complex domains. The first two difficulties can be successfully solved using the collocation procedure, while the third is still an open research topic.

The wavelet’s main feature is to facilitate the grid space adaptation and multi-resolution data compression. It enables solving problems with a sharp interface and narrow transition zone by changing their location and steepness in time and space. Wavelets utilize non-uniform grids dynamically adaptive according to the solution development. Any function, signal or data can be represented by a linear combination of basis functions (discrete wavelet transform) in multi-resolution fashion with different scales or frequencies and locations. This means that specific frequencies are associated with a particular spatial location that is not possible in classical Fourier transform (coefficients of linear combinations are wavelet coefficients which are associated with a specific resolution level [scale or frequency] and [collocation] point in space/time domain). This procedure is also known as multi-resolution analysis (MRA). Besides the solution variables (pressure, concentration, velocity) any other variables can be represented in the same multi-resolution fashion such as electrical or hydraulic conductivity, porosity, natural recharge or pumping. Furthermore, it is possible to use other basis functions with compact support (non-zero values only in one part of the domain) within the collocation method [46]. The spline adaptive collocation methods are described in [10], [36], [49].

Apart from wavelets and splines, there is a relatively lesser known class of atomic or R_bf basis functions (Rvachev’s basis functions) [37], [38]. Atomic functions are classified between classic polynomials and spline functions. However, in practice, their application as basis functions is closer to splines or wavelets. In this paper we use Fup basis functions which are one type of atomic basis functions. Gotovac and Kozulić [20] systemized the existing knowledge on atomic functions and presented its usage and calculation into a numerically applicable form. The application of Fup basis functions has been demonstrated in signal processing [32], [51], in initial value problems [21], and in the non-adaptive collocation method for boundary value problems [22], [31].

The main objective of this paper is to develop the Adaptive Fup Collocation Method (AFCM) and demonstrate its application to groundwater flow and transport problems. Presented is a novel adaptive Fup collocation method which is well suited to dealing with strong nonlinear groundwater problems with sharp fronts and narrow transition zones. A numerical procedure is implemented through a method of lines. Spatial discretization and grid adaptation are obtained by Fup collocation transform, while time integration is obtained by solving the system of Differential-Algebraic Equations (DAE). Furthermore, this method enables the adaptive multi-resolution evolution of a solution with resolved spatial and temporal scales and a desired level of accuracy. The numerical method has been tested and verified with a few characteristic groundwater flow and transport examples.

The following section presents a brief review of Fup basis functions. Section 3 describes the proposed collocation method; numerical examples are in Section 4, followed by conclusions in Section 5. Finally, two appendixes are added in order to show construction of the finite difference operator and convergence properties of the AFCM.

Section snippets

Fup basis functions preliminaries

Atomic or Rvachev’s basis functions – R_bf have a compact support and they are infinitely differentiable functions [20], [37]. They are classified in between classical polynomials and spline functions, but in practice their application as basis functions is closer to splines and wavelets.

Atomic functions, y(·), are defined as solutions of differential-functional equations of the following type: $Ly (x) = λ \sum_{k = 1}^{M} C_{k} y (ax - b_{k})$ where L is a linear differential operator with constant coefficients, λ is a

Adaptive Fup collocation method

This section deals with a detailed description of the Adaptive Fup Collocation Method (AFCM) including its features in space and time. Presented are all necessary steps for AFCM implementation to be reproducible.

Results and discussion

In this section, we present three examples for illustrating the efficiency and robustness of the Fup adaptive collocation method for modeling groundwater flow and transport problems. The first example is the Buckley–Leverett problem [16], [33] with a moving saturation front, i.e. a strong nonlinear problem of a multiphase flow. The next two examples illustrate the strength of the method in modeling the density driven flow and transport problems, i.e. the 1-D vertical density driven flow in a

Conclusions

In this paper we presented a multi-resolution adaptive Fup approach applicable to groundwater flow and transport problems with sharp gradients, fronts and narrow transition zones. AFCM enables adaptive multi-resolution evolution of the system dynamics resolving different spatial and temporal scales with a desired level of accuracy using the entire family of Fup basis functions. AFCM has been tested and verified on a few groundwater flow and transport problems, three 1-D problems with sharp

Acknowledgements

This research was supported and funded by Croatian Ministry of Science under project No. 4-0083-121; http:/www.gradst.hr/cer. The authors are grateful to Reviewers for their constructive and detailed comments which resulted in the improved version of the paper.

References (51)

M.A. Alves et al.
Adaptive multiresolution approach for solution of hyperbolic PDE’s
Comput Methods Appl Mech Eng
(2002)
M. Bause et al.
Computation of variably saturated subsurface flow by adaptive mixed hybrid finite element methods
Adv Water Resour
(2004)
S. Bertoluzza et al.
A wavelet collocation method for the numerical solution of PDEs
Appl Comput Harmonic Anal
(1996)
G. Beylkin et al.
On the adaptive numerical solution of nonlinear partial differential equations in wavelet bases
J Comput Phys
(1997)
G. Chiavassa et al.
Two adaptive wavelet algorithms for nolinear parabolic differential equations
Comput Fluids
(2002)
H. Class et al.
Numerical simulation of non-isothermal multiphase multicomponent processes in porous media. An efficient solution technique
Adv Water Resour
(2002)
P. Cruz et al.
Using wavelets for solving PDEs: and adaptive collocation method
Chem Eng Sci
(2001)
P. Cruz et al.
Solution of hyperbolic PDEs using a stable adaptive multiresolution method
Chem Eng Sci
(2003)
H.-J.G. Diersch et al.
Variable density flow and transport in porous media: approaches and challenges
Adv Water Resour
(2002)
E. Frind
Simulation of long-term transient density-dependent transport in groundwater
Adv Water Resour
(1982)

H. Gotovac et al.

An improved collocation method for solving the Henry problem

J Contam Hydrol

(2003)

J.S. Hesthaven et al.

A wavelet optimized adaptive multi-domain method

J Comput Phys

(1998)

O. Kolditz et al.

Coupled groundwater flow and transport: 1. Verification of variable-density flow and transport models

Adv Water Resour

(1998)

A. Kurganov et al.

New high resolution central scheme for nonlinear conservation laws and convection–diffusion equations

J Comput Phys

(2000)

C.T. Miller et al.

A spatially and temporally adaptive solution of Richards’ equation

Adv Water Resour

(2006)

P.A. Ramachandran et al.

Solution by triple collocation for periodic operation of heat regenerators

Comp Chem Eng

(1984)

J.C. Santos et al.

Adaptive multiresolution approach for two dimensional PDEs

Comp Methods Appl Mech Eng

(2004)

C.T. Simmons et al.

Variable-density groundwater flow and solute transport in heterogeneous media: approaches, resolutions and future challenges

J Contam Hydrol

(2001)

O.V. Vasilyev et al.

A fast adaptive wavelet collocation algorithm for multidimensional PDEs

J Comput Phys

(1997)

O.V. Vasilyev et al.

Second Generation wavelet collocation method for the solution of the partial differential equations

J Comput Phys

(2000)

R. Wang et al.

A High order global spatially adaptive collocation method for 1-D parabolic PDEs

Appl Numer Math

(2004)

P. Ackerer et al.

Modeling variable density flow and solute transport in porous medium: 1. Numerical model and verification

Transport Porous Media

(1999)

R. Andričević et al.

Relative dispersion for solute flux in aquifers

J Fluid Mech

(1998)

U. Ascher et al.

Computer methods for ordinary differential equations and differential-algebraic equations

(1998)

J. Bear

Conceptual and mathematical modeling

Cited by (32)

2-D local hp adaptive isogeometric analysis based on hierarchical Fup basis functions
2022, Computer Methods in Applied Mechanics and Engineering
Citation Excerpt :
For example, for function approximation, the function derivatives can be an ideal option in some cases. Furthermore, for BVP, the solution error between two resolution levels can be defined as criteria [45]. Moreover, satisfaction of the Peclet number can be very valuable for ADE problems (see [46]).
In this paper, 2-D hp local adaptive procedure is developed based on Control Volume Isogeometric Analysis (CV-IGA) and Hierarchical Fup (HF) basis functions. Contrary to the most common truncated hierarchical splines, HF enables local hp adaptation because higher resolution levels do not include only basis with smaller compact support or higher frequencies, but also with higher order. Consequence of this property is spectral convergence of the proposed adaptive algorithm which is presented on classical benchmarks such as L-shape benchmark and advection dominated problems. Even in non-smooth problems, spectral convergence is achieved contrary to the application of uniform grid. CV-IGA ensures local and global mass conservation which is potentially very important for fluid mechanics problems. 2-D proposed algorithm chooses regular control volumes in parametric space at all resolution levels closely related to the Greville points (vertices) of basis functions. Therefore, methodology is very simple requiring only overlapping of control volumes in the areas where different levels are connected, while its computational cost lies between Galerkin and collocation formulations.
Control Volume Isogeometric Analysis for groundwater flow modeling in heterogeneous porous media
2021, Advances in Water Resources
Citation Excerpt :
The choice of basis functions depends on the problem and required continuity. Future research can include the development of fully space and time adaptive (Kamber et al., 2020) IGA methods for efficient modeling of flow, as well as for contaminant transport or saltwater intrusion, which needs to involve a velocity dependent dispersion tensor with continuous velocity field and its derivatives (e.g. Gotovac et al., 2007). Finally, resolving of singularities due to non-smooth description of, for example, geometry, pumping/recharge or heterogeneity needs to be solved by using adaptive approaches with higher order spline basis functions.
Many important groundwater transport applications require solving the Darcy flow in heterogeneous porous media. Flow simulations, especially in large, highly heterogeneous aquifers, require extensive computational resources, a multiresolution (multiscale) approach to resolve the different heterogeneity scales and an accurate calculation of the velocity field. Common methods, such as finite volumes and elements, assume a discontinuous conductivity field introducing velocity discontinuities along the cell or element interfaces due to using classic discrete operators or Lagrangian basis functions. Over the last decade, the development of isogeometric analysis (IGA) eliminates many of the aforementioned limitations bridging the gap between CAD and numerical analysis. Since classic IGA uses the Galerkin and collocation approach, in this paper, we present a third concept in the form of Control Volume IsoGeometric Analysis (CV-IGA) enabling local and global mass conservation as well as multiresolution description of all heterogeneity scales. Due to the approximation properties of spline basis functions, the velocity field and its derivatives are continuous and are obtained by an optimal convergence rate. The CV-IGA methodology is verified with 2-D numerical and stochastic benchmark flow simulations, including comparisons with classic methods and two other IGA formulations as well as the convergence analysis of the head and velocity fields for different orders of Fup and B-spline basis functions.
Numerical benchmark study for flow in highly heterogeneous aquifers
2020, Advances in Water Resources
Citation Excerpt :
Such order relations also served as guidelines for implementing solutions with finite difference methods (FDM) (de Dreuzy et al., 2007; Kurbanmuradov and Sabelfeld, 2010; Trefry et al., 2003) or finite element methods (FEM) (Bellin et al., 1992; Radu et al., 2011; Salandin and Fiorotto, 1998) used in Monte Carlo simulations. In spectral approaches (Gotovac et al., 2007; Li and McLaughlin, 1991) or probabilistic collocation and chaos expansion approaches (Li and Zhang, 2007; Li et al., 2009) the spatial resolution is also related to the number of terms in the series representation of the coefficients and of the solution of the flow problem. Despite the advances in numerical methods, computing accurate flow solutions for highly heterogeneous formations faces computational challenges in terms of code efficiency and computational resources (see e.g., Gotovac et al., 2009; Li and Zhang, 2007).
This article presents numerical investigations on accuracy and convergence properties of several numerical approaches for simulating steady state flows in heterogeneous aquifers. Finite difference, finite element, discontinuous Galerkin, spectral, and random walk methods are tested on two-dimensional benchmark flow problems. Realizations of log-normal hydraulic conductivity fields are generated by Kraichnan algorithms in closed form as finite sums of random periodic modes, which allow direct code verification by comparisons with manufactured reference solutions. The quality of the methods is assessed for increasing number of random modes and for increasing variance of the log-hydraulic conductivity fields with Gaussian and exponential correlation. Experimental orders of convergence are calculated from successive refinements of the grid. The numerical methods are further validated by comparisons between statistical inferences obtained from Monte Carlo ensembles of numerical solutions and theoretical first-order perturbation results. It is found that while for Gaussian correlation of the log-conductivity field all the methods perform well, in the exponential case their accuracy deteriorates and, for large variance and number of modes, the benchmark problems are practically not solvable with reasonably large computing resources, for all the methods considered in this study.
Impact of hydrogeological data on measures of uncertainty, site characterization and environmental performance metrics
2012, Advances in Water Resources
The significance of conditioning predictions of environmental performance metrics (EPMs) on hydrogeological data in heterogeneous porous media is addressed. Conditioning EPMs on available data reduces uncertainty and increases the reliability of model predictions. We present a rational and concise approach to investigate the impact of conditioning EPMs on data as a function of the location of the environmentally sensitive target receptor, data types and spacing between measurements. We illustrate how the concept of comparative information yield curves introduced in de Barros et al. [de Barros FPJ, Rubin Y, Maxwell R. The concept of comparative information yield curves and its application to risk-based site characterization. Water Resour Res 2009;45:W06401. doi:10.1029/2008WR007324] could be used to assess site characterization needs as a function of flow and transport dimensionality and EPMs. For a given EPM, we show how alternative uncertainty reduction metrics yield distinct gains of information from a variety of sampling schemes. Our results show that uncertainty reduction is EPM dependent (e.g., travel times) and does not necessarily indicate uncertainty reduction in an alternative EPM (e.g., human health risk). The results show how the position of the environmental target, flow dimensionality and the choice of the uncertainty reduction metric can be used to assist in field sampling campaigns.
Risk characterization for toxic chemicals transported in aquifers
2012, Advances in Water Resources
Citation Excerpt :
The AFMCM methodology is based on Fup basis functions with compact support (related to the other localized basis functions such as splines or wavelets) and the Fup collocation transform (FCT), which is closely related to the discrete Fourier transform needed for the multi-resolution representation of all flow and transport variables. Fup basis functions and the FCT are presented in [21]. Other improved Monte Carlo (MC) methodology steps were also utilized: (a) the Fup regularized transform (FRT) for data or function (e.g., log-conductivity) approximations in the same multi-resolution fashion as FCT, but computationally more efficient, (b) the adaptive Fup collocation method (AFCM) for approximation of the flow differential equation, (c) the random walk particle tracking algorithm based on the Runge–Kutta-Verner explicit time integration scheme and FRT and (d) MC statistics represented by Fup basis functions.
The risk characterization resulting from the introduction of toxic chemicals in a subsurface flow field is presented. The concept of concentration threshold is used to quantify the risk associated with non-carcinogenic chemicals introduced to the population by ingestion of groundwater as the exposure pathway. The risk assessment methodology presented uses the expected mass fraction (EMF) concept with exposure duration to identify the distribution of dosage over different concentrations during the plume migration over the well location. The numerical simulation of the subsurface transport by advection and local diffusion is used to produce the concentration plume that passes different locations of interest. The EMF obtained presents the probability of the expected mass above some concentration threshold found at the location of interest. The risk formulation is defined with the risk reliability (safety) and its complement, the risk exceedance (failure) value. The risk characterization is obtained as a probability for exceeding the human reference dose which is considered uncertain due to the necessary extrapolation between concentration used in toxicological studies and the concentration to which humans could be exposed in nature. The final risk assessment expression is derived in a closed form by coupling the expected mass fraction with the safe human threshold concentration probability density function (pdf) inferred from the toxicological studies. The results indicate the importance of estimating the probability of a concentration mass found at locations of interest together with its exposure duration. The exposure duration was revealed to be an important parameter that needs to be estimated depending on the human concentration threshold selected and the distance from the source. The results in terms of the risk safety and risk failure also indicate the dilution effect on the passing concentration plume in the subsurface as a function of the distance and orientation from the source. Inclusion of uncertainty in the selection of the human concentration threshold shows the important effect on the risk quantification.
Dense, viscous brine behavior in heterogeneous porous medium systems
2010, Journal of Contaminant Hydrology
The behavior of dense, viscous calcium bromide brine solutions used to remediate systems contaminated with dense nonaqueous phase liquids (DNAPLs) is considered in laboratory and field porous medium systems. The density and viscosity of brine solutions are experimentally investigated and functional forms fit over a wide range of mass fractions. A density of 1.7 times, and a corresponding viscosity of 6.3 times, that of water is obtained at a calcium bromide mass fraction of 0.53. A three-dimensional laboratory cell is used to investigate the establishment, persistence, and rate of removal of a stratified dense brine layer in a controlled system. Results from a field-scale experiment performed at the Dover National Test Site are used to investigate the ability to establish and maintain a dense brine layer as a component of a DNAPL recovery strategy, and to recover the brine at sufficiently high mass fractions to support the economical reuse of the brine. The results of both laboratory and field experiments show that a dense brine layer can be established, maintained, and recovered to a significant extent. Regions of unstable density profiles are shown to develop and persist in the field-scale experiment, which we attribute to regions of low hydraulic conductivity. The saturated–unsaturated, variable-density groundwater flow simulation code SUTRA is modified to describe the system of interest, and used to compare simulations to experimental observations and to investigate certain unobserved aspects of these complex systems. The model results show that the standard model formulation is not appropriate for capturing the behavior of sharp density gradients observed during the dense brine experiments.

View all citing articles on Scopus

View full text

Multi-resolution adaptive modeling of groundwater flow and transport problems

Abstract

Introduction

Section snippets

Fup basis functions preliminaries

Adaptive Fup collocation method

Results and discussion

Conclusions

Acknowledgements

Comput Methods Appl Mech Eng

Adv Water Resour

Appl Comput Harmonic Anal

J Comput Phys

Comput Fluids

Adv Water Resour

Chem Eng Sci

Chem Eng Sci

Adv Water Resour

Adv Water Resour

J Contam Hydrol

J Comput Phys

Adv Water Resour

J Comput Phys

Adv Water Resour

Comp Chem Eng

Comp Methods Appl Mech Eng

J Contam Hydrol

J Comput Phys

J Comput Phys

Appl Numer Math

Modeling variable density flow and solute transport in porous medium: 1. Numerical model and verification

Transport Porous Media

Relative dispersion for solute flux in aquifers

J Fluid Mech

Computer methods for ordinary differential equations and differential-algebraic equations

Conceptual and mathematical modeling