Hydrogeophysical model calibration and uncertainty analysis via full integration of PEST/PEST++ and COMSOL

https://doi.org/10.1016/j.envsoft.2021.105183Get rights and content

Highlights

  • We present a strategy for hydrogeophysical inversion and uncertainty analysis.

  • PEST/PEST++ and COMSOL Multiphysics are fully integrated.

  • The approach is applied to electrical resistivity monitoring in a coastal aquifer.

  • The coupled inversion achieves better delineation of seawater intrusion.

  • Uncertainty analysis is done using PEST++ Iterative Ensemble Smoother.

Abstract

Calibration of groundwater models is frequently hampered by limited availability of hydrogeological data. Non-intrusive geophysical methods are increasingly used to provide higher spatio-temporal resolution datasets for identification of hydrological processes and estimation of hydraulic properties. In groundwater model calibration performed through joint or coupled hydrogeophysical inversion, the hydrogeological datasets are supplemented with auxiliary geophysical data. In this work we propose a methodological approach to perform coupled inversion by integrating the calibration software PEST/PEST++ with COMSOL Multiphysics using MATLAB. The strategy provides multiple options for calibration and uncertainty analysis relevant for a broad range of environmental models. To illustrate the approach, we show a hydrogeophysical application in which electrical resistivity is used jointly with borehole data for the identification of seawater intrusion in a coastal aquifer.

Introduction

Groundwater management and decision making is frequently supported by outcomes from numerical models (Doherty and Simmons, 2013, Ferré, 2017). To have predictive capabilities numerical groundwater models are traditionally calibrated (e.g., Anderson et al., 2015) and further evaluated to quantify prediction uncertainty (Linde et al., 2017). The uncertainty analysis can be conducted either using a deterministic (error propagation analysis) or stochastic approach, even avoiding the necessity of a pre-calibrated model (e.g., Scheidt et al., 2018; Hermans et al., 2018, 2019). Model calibration and uncertainty analysis in groundwater modelling, as for other subsurface systems, is frequently hampered by the scarcity of direct hydrogeological data, such as head or concentration measurements from wells or boreholes. As drilling is expensive, intrusive and, generally, provides spatially scattered point information, key groundwater processes and spatial distribution of hydrogeological properties are not fully captured, more so in heterogeneous aquifer systems (e.g., Zheng and Gorelick, 2003).

Geophysical methods have been, and are, routinely used in hydrogeological studies to provide indirect information of the subsurface (Kirsch, 2006). Traditionally and most often, geophysical methods are processed and interpreted independently and used qualitatively in the initial stages of the groundwater modelling workflow, typically for geometrical delineation of geological bodies and structures or the mapping of hydrological features (e.g. Mastrocicco et al., 2010) or processes of interest, such as the location of the water table (Buchanan, S., & Triantafilis, J., 2009), the identification of saltwater intrusion (de Franco et al., 2009), soil moisture variations (Chambers et al., 2014) or the tracking of contaminant plumes (Gasperikova et al., 2012).

In the last two decades, there has been a growing interest in a more quantitative use of geophysical information in hydrogeology including their integration in groundwater models, eventually leading to the emergence of a new scientific sub-discipline known as hydrogeophysics (Binley et al., 2015). At the core of the hydrogeophysical approach lies the use or definition of petrophysical relationships (e.g., Archie, 1942 and its variants; Slater, 2007), commonly derived at the laboratory. These are required to link geophysical properties (e.g. electrical resistivities) to hydrogeological parameters (e.g. storage properties; Mezquita-González et al., 2021) or state variables (e.g. water salinity content; Klotzsche et al., 2018).

Previous authors have shown limitations of the capability of geophysical inversion to accurately resolve subsurface features (Day-Lewis et al., 2005; Singha and Moysey, 2006), typically providing a blurry, smoother, and sometimes distorted representation of reality, as well as significant uncertainties in the values of the parameter measured. Although alternative regularization approaches may improve the inverted model (Hermans et al., 2012; Bouchedda et al., 2017; Thibaut et al., 2021), recent studies (e.g. Revil et al., 2017; Brunetti and Linde, 2018; González-Quirós and Comte, 2020) have also shown the importance of conceptual and structural errors in the hydrogeophysical workflow, such as the assumption of homogeneity in heterogeneous systems when parameterizing the petrophysical model or the incorrect selection of the petrophysical model itself. As a result, the quantitative use of information derived from geophysical inversion in groundwater model parametrization (i.e. direct mapping) or calibration could lead, if incorrectly applied, to unrealistic or erroneous property distributions or estimations, and consequently, inaccurate hydrological interpretations. Awareness of these limitations and incorporation of their associated uncertainty in the workflow (e.g. Beaujean et al., 2014; Hermans and Irving, 2017) have been shown to reduce some of these errors associated with the use quantitative geophysical information in hydrogeophysical applications.

To overcome these limitations, it was proposed (e.g. Lebbe, 1999) that a more appropriate strategy for an efficient integration of geophysics in groundwater modelling calibration is through a coupled hydrogeophysical modelling. In it, both the geophysical and hydrogeological observations, rather than inversion results, are simultaneously used as model verification datasets (Comte and Banton, 2007), or in an automated way through coupled hydrogeophysical inversion in which geophysical and hydrogeological measurements are used as input datasets for automatic groundwater model calibration (Hinnell et al., 2010;Linde and Doestch, 2016). While the approach is still subject to practical and structural errors, their main advantage is that the significantly extended multiphysical observation datasets provide additional constraints in the estimation of the ranges and distributions of hydrological parameters and variables (or quantities), and in the resolution of groundwater processes, which are the main objectives in hydrogeological studies (Linde et al., 2015). Coupled modelling approaches start with the simulation of the hydrological model for some defined distribution of hydrological properties and boundary conditions. A petrophysical relationship, which can be uncertain (e.g. Irving and Singha, 2010), is used to obtain a spatial distribution of geophysical properties from which the forward geophysical model is simulated to obtain the geophysical response at locations of interest. Model predicted geophysical and hydrologic observations are compared with field measurements through an iterative procedure where the hydrological-geophysical model properties are sequentially adjusted until an acceptable fit is obtained. In the coupled inversion approaches the iterative comparison and minimisation of differences between modelled and observed values is automated.

A key difficulty for development and application of coupled hydrogeophysical inversion for real world applications is the requirement to code or program the necessary governing equations of both, the groundwater and —at least one— geophysical problems linked with the appropriate petrophysical relationship. Previous authors have used specific software –e.g. MODFLOW (Harbaugh et al., 2005), SEAWAT (Langevin et al., 2008), SUTRA (Voss and Provost, 2002) or FEFLOW (Diersch, 2013) for the groundwater problem, and RES2DMOD/RES2DINV (Loke, 2018), BERT (Rücker et al., 2006) or R2 (Binley and Kemna, 2005) for the geophysical problem– to solve the hydrogeophysical model sequentially (e.g. Herckenrath et al., 2013; Kang et al., 2019) or have developed solutions for hydrogeophysical inversion that usually require the development or programming of in-house codes for specific applications of (Pollock and Cirpka, 2012; Steklova and Haber, 2017). On this sense, open-source libraries under continuous development, such as SimPEG (Cockett et al., 2015) or PyGIMLI (Rücker et al., 2017), provide a wide range of solutions to perform coupled hydrogeophysical inversion or be integrated with a groundwater model.

However, although coupled modelling is usually achievable for models with simple geometries, representing scenarios that can be straightforwardly conceptualized, and for small multiphysical calibration datasets, the implementation of real-world aquifers targeting complex hydrological processes and/or large datasets (e.g. Comte et al., 2017) can be time consuming and would benefit from use of automatic calibration and a much higher degree of flexibility for geometry, meshing, parametrization and implementation of the governing equations and boundary conditions. Additionally, multiphysical model calibration and uncertainty analysis demands a robust framework —aside from the mathematical engine itself— to handle the multiple types of observation datasets and program the specific relationships between parameter types with geologically realistic distributions of properties (Linde et al., 2015).

In this work, we propose a methodology to implement a fully coupled multiphysical inversion for the calibration and uncertainty analysis in hydrogeophysical environmental models. The strategy is achieved by integrating a fully coupled hydrogeophysical forward model, developed with the commercial finite-element software COMSOL Multiphysics®, and the model-independent calibration software PEST (Doherty, 2020) and PEST++ (White et al., 2020), widely used by the groundwater community. The coupling of PEST and COMSOL was performed previously by Halloran et al. (2019), who developed a Java interface for communication. Here we follow a different strategy for communication using MATLAB instead of Java, and we include some novelties that were not implemented in COMPEST, such as the use of pilot points for spatial parameterization, Tikhonov regularization, application of the novel iterative Ensemble Smoother from PEST++, parallelization and a first application of the strategy to a hydrogeophysical problem.

The paper is organised as follows. In section 2.1 we present some basic concepts and theory behind PEST and PEST++. In section 2.2 we explain the methodological workflow for integration of PEST/PEST++ and COMSOL with model parallelization. In section 2.3 we describe the hydrogeophysical modelling and inversion and present the example model. In section 3 we show the hydrogeophysical example application results: in section 3.1 the results of a Monte Carlo forward coupled hydrogeophysical model, in section 3.2 the results of model calibration and in section 3.3 the results of the uncertainty analysis obtained with application of the iterative ensemble smoother. Finally, in section 4, we discuss the capabilities of the multiphysical inversion for calibration of environmental models.

Section snippets

PEST and PEST++

The model-independent calibration software PEST (Doherty, 2020) has been widely used since the 90s by the groundwater community for model calibration. PEST relies on the use of input and output files to interact with any numerical model. It includes a suite of functionalities for model parameterization, data conversion or setup of the calibration workflow, including parallelization. PEST++ (White et al., 2020) development began in 2009 within the USGS and contained much (but not all) of the

Stochastic groundwater modelling results

Each single forward coupled model was computed in less than 80 s in an Intel® Core™ i5-8500 with 40 GB RAM. For the forward stochastic modelling exercise, from the total of 500 model realizations of the random field, 8 failed to converge (i.e., less than 2%), mainly because of a particular arrange of low hydraulic conductivities in the density-dependent flow model.

Representative computed results are compiled in Fig. 4. Fig. 4a shows the positions of the drinking water limit contour

Comparison with standalone ERT inversion

It is well known that stand-alone inversion of geophysical data is a blurry, usually a smoother, representation of reality in which small scale features of the subsurface are not well recovered because of resolution limitations (Day-Lewis et al., 2005; Singha and Moysey, 2006). Additionally, there is a loss of resolution with depth that, when studying saltwater intrusion in coastal aquifers, has resulted in discrepancies when comparing tomograms obtained from surface ERI with borehole salinity

Conclusions

The full integration of geophysical data in the groundwater modelling workflow has been identified as a step forward to provide more reliable models used for decision support. Fully coupled hydrogeophysical inversion has been proposed as one of the solutions to improve the integration of geophysics in hydrogeological studies. However, its use in real-world settings has been limited because of the difficulty of simulating complex field conditions, the unavailability of flexible strategies and

Software availability

PEST (Doherty, 2020) is freely available in its dedicated webpage https://pesthomepage.org/.

PEST++ (White et al., 2020) is accessible and available in the US Geological Survey webpage https://www.usgs.gov/software/pest-parameter-estimation-code-optimized-large-environmental-models.

MATLAB is a commercial numerical computing platform (https://uk.mathworks.com/, or the respective webpage for each county). Version 2018a was used in this work.

COMSOL Multiphysics is a commercial finite element

Declaration of competing interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Acknowledgments

Andrés González Quirós is the recipient of a Royal Society – Newton International Fellowship (NIF\R1\182210), hosted at the University of Aberdeen. We thank the Scottish Funding Council/Scottish Alliance for Geoscience, Environment and Society for seed funding the development of the project. We would also like to thank the Associate Editor Tim Green, Thomas Hermans, Landon Halloran and one anonymous reviewer for their comments and suggestions in the revision process.

References (97)

  • L.J. Halloran et al.

    COMPEST, a PEST-COMSOL interface for inverse multiphysics modelling: development and application to isotopic fractionation of groundwater contaminants

    Comput. Geosci.

    (2019)
  • L.J. Halloran et al.

    Sorption-and diffusion-induced isotopic fractionation in chloroethenes

    Sci. Total Environ.

    (2021)
  • T. Hermans et al.

    Imaging artificial salt water infiltration using electrical resistivity tomography constrained by geostatistical data

    J. Hydrol.

    (2012)
  • X. Kang et al.

    Coupled hydrogeophysical inversion to identify non-Gaussian hydraulic conductivity field by jointly assimilating geochemical and time-lapse geophysical data

    J. Hydrol.

    (2019)
  • L. Lebbe

    Parameter identification in fresh-saltwater flow based on borehole resistivities and freshwater head data

    Adv. Water Resour.

    (1999)
  • N. Linde et al.

    Geological realism in hydrogeological and geophysical inverse modeling: a review

    Adv. Water Resour.

    (2015)
  • N. Linde et al.

    On uncertainty quantification in hydrogeology and hydrogeophysics

    Adv. Water Resour.

    (2017)
  • J. Lopez-Alvis et al.

    A cross-validation framework to extract data features for reducing structural uncertainty in subsurface heterogeneity

    Adv. Water Resour.

    (2019)
  • J.A. Mezquita-González et al.

    Quantification of groundwater storage heterogeneity in weathered/fractured basement rock aquifers using electrical resistivity tomography: sensitivity and uncertainty associated with petrophysical modelling

    J. Hydrol.

    (2021)
  • J.C. Refsgaard et al.

    Review of strategies for handling geological uncertainty in groundwater flow and transport modeling

    Adv. Water Resour.

    (2012)
  • C. Rücker et al.

    pyGIMLi: an open-source library for modelling and inversion in geophysics

    Comput. Geosci.

    (2017)
  • R. Thibaut et al.

    A new workflow to incorporate prior information in minimum gradient support (MGS) inversion of electrical resistivity and induced polarization data

    J. Appl. Geophys.

    (2021)
  • C.H.M. Tso et al.

    Integrated hydrogeophysical modelling and data assimilation for geoelectrical leak detection

    J. Contam. Hydrol.

    (2020)
  • J.T. White

    A model-independent iterative ensemble smoother for efficient history-matching and uncertainty quantification in very high dimensions

    Environ. Model. Software

    (2018)
  • H. Zhou et al.

    Inverse methods in hydrogeology: evolution and recent trends

    Adv. Water Resour.

    (2014)
  • M.P. Anderson et al.

    Applied Groundwater Modeling: Simulation of Flow and Advective Transport

    (2015)
  • G.E. Archie

    The electrical resistivity log as an aid in determining some reservoir characteristics

    Trans. Am. Inst. Min. Metall. Pet. Eng.

    (1942)
  • J. Beaujean et al.

    Calibration of seawater intrusion models: inverse parameter estimation using surface electrical resistivity tomography and borehole data

    Water Resour. Res.

    (2014)
  • C.R. Berg

    A simple effective-medium model for water saturation in porous rocks

    Geophysics

    (1995)
  • A. Binley et al.

    DC resistivity and induced polarization methods

  • A. Binley et al.

    The emergence of hydrogeophysics for improved understanding of subsurface processes over multiple scales

    Water Resour. Res.

    (2015)
  • A. Bouchedda et al.

    Constrained ERT Bayesian inversion using inverse Matérn covariance matrix

  • S. Buchanan et al.

    Mapping water table depth using geophysical and environmental variables

    Ground Water

    (2009)
  • M. Camporese et al.

    Coupled and uncoupled hydrogeophysical inversions using ensemble Kalman filter assimilation of ERT‐monitored tracer test data

    Water Resour. Res.

    (2015)
  • J. Carrera et al.

    Computational and conceptual issues in the calibration of seawater intrusion models

    Hydrogeol. J.

    (2010)
  • J.E. Chambers et al.

    4D electrical resistivity tomography monitoring of soil moisture dynamics in an operational railway embankment

    Near Surf. Geophys.

    (2014)
  • Y. Chen et al.

    Levenberg–Marquardt forms of the iterative ensemble smoother for efficient history-matching and uncertainty quantification

    Comput. Geosci.

    (2013)
  • J.C. Comte et al.

    Cross‐validation of geo‐electrical and hydrogeological models to evaluate seawater intrusion in coastal aquifers

    Geophys. Res. Lett.

    (2007)
  • J.C. Comte et al.

    Effect of volcanic dykes on coastal groundwater flow and saltwater intrusion: a field-scale multiphysics approach and parameter evaluation

    Water Resour. Res.

    (2017)
  • S.C. Constable et al.

    Occam's inversion: a practical algorithm for generating smooth models from electromagnetic sounding data

    Geophysics

    (1987)
  • A. Costall et al.

    Electrical resistivity imaging and the saline water interface in high-quality coastal aquifers

    Surv. Geophys.

    (2018)
  • A.R. Costall et al.

    Groundwater throughflow and seawater intrusion in high quality coastal aquifers

    Sci. Rep.

    (2020)
  • F.D. Day-Lewis et al.

    The application of petrophysical models to radar and electrical resistivity tomograms: resolution dependent limitations

    J. Geophys. Res.

    (2005)
  • G. de Marsily et al.

    Interpretation of interference tests in a well field using geostatistical techniques to fit the permeability distribution in a reservoir model

  • H.J.G. Diersch

    FEFLOW: Finite Element Modeling of Flow, Mass and Heat Transport in Porous and Fractured Media

    (2013)
  • J. Doherty

    Calibration and uncertainty analysis for complex environmental models

    Watermark Numerical Comp

    (2015)
  • J. Doherty

    PEST, Model-independent Parameter Estimation

    (2020)
  • J.E. Doherty et al.

    Approaches to highly parameterized inversion: pilot-point theory, guidelines, and research directions

    US Geological Survey scientific investigations report

    (2010)
  • Cited by (12)

    • Joint estimation of groundwater salinity and hydrogeological parameters using variable-density groundwater flow, salt transport modelling and airborne electromagnetic surveys

      2022, Advances in Water Resources
      Citation Excerpt :

      As the 3D-VDG model computation step cost ∼5 h per iteration, one could simply run this step externally on large computational clusters using iMOD-SEAWAT, which utilizes distributed memory parallelization for faster computation times (Verkaik et al., 2021). In a practical setting however, we suggest fully parallelizing the optimization itself, using methods such as evolutionary algorithms (Brauer et al., 2002; Mühlenbein et al., 1991) and parallel Bayesian optimization (e.g. González et al., 2015; Kandasamy et al., 2017) where function evaluations can be done in parallel rather than sequentially. The latter approach has been successfully used in hydrocarbon reservoir modelling (e.g. Abdollahzadeh et al., 2011), where similar to the optimization used in this study, inverse problems are solved that require multiple flow-simulations, which are more computationally expensive compared to our problem.

    View all citing articles on Scopus
    View full text