Quantifying geological uncertainty for flow and transport modeling in multi-modal heterogeneous formations

Multiple point statistics (MPS) provides a systematic approach for pattern-based simulation of geologic objects from a conceptual training image (TI). The TI encodes the higher-order spatial statistics of the expected connectivity structures through stationary patterns representing the underlying geologic features. The pattern-imitating nature of MPS simulation implies that the simulated facies inherit the spatial structure of the general features in the TI. This property makes the MPS approach very sensitive to uncertainty in the prior TI. Since TIs are constructed using uncertain data and imperfect assumptions, multiple TIs may be necessary to account for the uncertainty and full range of structural variability in facies descriptions. We present a Bayesian mixture modeling approach for adaptively sampling conditional facies from multiple uncertain TIs using a probability conditioning method (PCM). Using the PCM, we invert the flow data to obtain a facies probability map for drawing conditional facies realizations from each TI. The number of samples drawn from each TI is proportional to the weight assigned to them. The TI weights are assigned based on the predictive performance of its corresponding conditional facies realizations. We demonstrate the suitability of the proposed method using numerical experiments in fluvial formations.

<p><span>In the beginning of the 2000's [1], multiple-point statistics (MPS) was introduced as a novel geostatistical approach to explore the variability of natural phenomena in a realistic way by observing and simulating data patterns, sensibly improving the preservation of connectivity and shape of the modeled structures.</span></p><p><span>A usual requirement for MPS is the presence of complete and representative training images (TI), showing clear and possibly redundant examples of the studied structures. But in the everyday practice, this information is often partially or scarcely available, strongly limiting the use of MPS.</span></p><p><span>In this presentation we start with an overview of MPS strategies proposed to overcome training data limitations. We consider different examples of multisite rain-gauge networks containing sparse data gaps, with the goal of estimating the missing data, using the same incomplete dataset as TI [2]. Another considered study case regards the use of 2D training images of geological outcrops used to reconstruct a 3D volume of fluvioglacial deposits [3]. </span></p><p><span>We then consider a common problem in hydroclimatological studies: the bias correction of weather radar images with ground rainfall measurements. This is a typical no-TI problem where there is no example of unbiased grid image to train MPS. In this case, we propose a novel pattern-to-point approach, where we create a catalog of local grid patterns, each one associated to a rainfall measurement. This way the MPS algorithm 1) selects ungauged locations, 2) searches similar grid patterns in the catalog, and 3) projects the linked historical ground measurements at the ungauged locations.</span></p><p><span>From early results, this technique seems to recover hidden spatial patterns which correct the highly non-linear bias by extracting information from the pattern-to-point catalog. This is a first step for MPS towards the use of TIs integrating variables of different dimensionality, opening a new methodological path for future research.</span></p><p> </p><p><span>BIBLIOGRAPHY</span></p><p><span>[1] Strebelle, S. "Conditional simulation of complex geological structures using multiple-point statistics." Mathematical geology 34.1 (2002): 1-21.</span></p><p><span>[2] Oriani, F. et al. "Missing data imputation for multisite rainfall networks: a comparison between geostatistical interpolation and pattern-based estimation on different terrain types." Journal of Hydrometeorology 21.10 (2020): 2325-2341.</span></p><p><span>[3] Kessler, T. et al. "Modeling fine</span><span>‐</span><span>scale geological heterogeneity—examples of sand lenses in tills." Groundwater 51.5 (2013): 692-705.</span></p>

Multiple-point statistics are widely used for the simulation of categorical variables because the method allows for integrating a conceptual model via a training image and then simulating complex heterogeneous fields. The multiple-point statistics inferred from the training image can be stored in several ways. The tree structure used in classical implementations has the advantage of being efficient in terms of CPU time, but is very RAM demanding and then implies limitations on the size of the template, which serves to make a proper reproduction of complex structures difficult. Another technique consists in storing the multiple-point statistics in lists. This alternative requires much less memory and allows for a straightforward parallel algorithm. Nevertheless, the list structure does not benefit from the shortcuts given by the branches of the tree for retrieving the multiple-point statistics. Hence, a serial algorithm based on list structure is generally slower than a tree-based algorithm. In this paper, a new approach using both list and tree structures is proposed. The idea is to index the lists by trees of reduced size: the leaves of the tree correspond to distinct sublists that constitute a partition of the entire list. The size of the indexing tree can be controlled, and then the resulting algorithm keeps memory requirements low while efficiency in terms of CPU time is significantly improved. Moreover, this new method benefits from the parallelization of the list approach.

Building a training image with Digital Outcrop Models

Two approaches are traditionally used to build numerical models for facies distributions within a reservoir. Pixel-based techniques aim at generating simulated realizations that honor the well data values, and reproduce a given variogram which models two-point spatial correlation. However, because the variogram cannot look at spatial continuity between more than two locations at a time, pixel-based algorithms give poor representations of the actual facies geometries. In contrast, object-based techniques allow reproducing crisp geometries, but the conditioning on well data requires iterative "trial-and-error" corrections, which can be time-consuming, particularly when the data are dense with regard to the average object size. This paper presents a new approach that combines the easy conditioning of pixel-based algorithms with the ability to reproduce "shapes" of object-based techniques, without being too time and memory demanding. In this new approach, the complex geological structures expected to be present in the reservoir are characterized by multiple-point statistics, which express joint variability at many more than two locations at a time. Such multiple-point statistics cannot be inferred from typically sparse well data but could be read from training images depicting the expected subsurface heterogeneities. A training image need not carry any locally accurate information on the reservoir; it need only reflect a prior stationary geological/structural concept. Thus training images can be generated by object-based algorithms freed of the constraint of data conditioning. The multiple-point statistics inferred from the training image(s) are then exported to the reservoir model, where they are anchored to the well data using a pixel-based sequential simulation algorithm. This algorithm is tested for the simulation of a turbidite system where flow is controlled by meandering channels with cross-bedding. The training image reflecting the channel patterns is an unconditional realization generated by an object-based algorithm. The final simulated numerical models reproduce these channel patterns, and honor exactly all well data values at their locations. The methodology proposed appears to be practical, general, and fast.

Reservoirs result from complex and intricate natural formation processes, and are therefore very heterogeneous at multiple scales. As fluid displacements are very sensitive to heterogeneity, it is crucial to understand and characterize how petrophysical properties vary in space to be able to predict oil production. Multiple point statistics has become a common tool to build models representing geological formations comprising objects characterized by complex geometries such as curvilinear shapes. The idea consists in referring to a training image depicting the expected geological structures, and then in generating images reproducing the multiple point statistics inferred from the training image. We investigate a substitute to the commonly referenced multiple-point methodologies that was initially developed in computer graphics: it is called texture synthesis. We modify the way it works so that it can handle reservoir simulation issues. Additionally, we focus on its multiscale variant to make it more efficient. Within this multiscale context, the basic training image is used to create a set of images with decreasing resolutions that are all viewed as training images. Then, a realization is simulated one scale at a time, starting from the level of lowest resolution. When simulating the realization at a given scale, the selection of the value to be assigned to a grid block involves various levels of comparison. We still compare the neighbors of the grid block under consideration with the events extracted from the training image corresponding to the same level of resolution as done for multipoint simulation. In addition, we compare its ancestors with the events extracted from the coarser training images. The ancestors are the values already simulated to populate the coarser grids. In such conditions, the first realization generated at the coarsest scale is able to capture the largest geological objects or heterogeneities. Then, it is used as a constraint that is propagated scale by scale towards the finest resolution grid. This makes it possible to use small size neighborhoods. A few improvements are then suggested to make the algorithm more efficient without deteriorating the simulated images. Last, we present numerical tests performed with complex training images.

We propose a smooth formulation of multiple-point statistics that enables us to solve inverse problems using gradient-based optimization techniques. We introduce a differentiable function that quantifies the mismatch between multiple-point statistics of a training image and of a given model. We show that, by minimizing this function, any continuous image can be gradually transformed into an image that honors the multiple-point statistics of the discrete training image. The solution to an inverse problem is then found by minimizing the sum of two mismatches: the mismatch with data and the mismatch with multiple-point statistics. As a result, in the framework of the Bayesian approach, such a solution belongs to a high posterior region. The methodology, while applicable to any inverse problem with a training-image-based prior, is especially beneficial for problems which require expensive forward simulations, as, for instance, history matching. We demonstrate the applicability of the method on a two-dimensional history matching problem. Starting from different initial models we obtain an ensemble of solutions fitting the data and prior information defined by the training image. At the end we propose a closed form expression for calculating the prior probabilities using the theory of multinomial distributions, that allows us to rank the history-matched models in accordance with their relative posterior probabilities.

When using coarse-resolution remote sensing images, super-resolution reconstruction is widely desired, and can be realized by reproducing the intrinsic features from a set of coarse-resolution fraction data to fine-resolution remote sensing images that are consistent with the coarse fraction information. Prior models of spatial structures that encode the expected features at the fine (target) resolution are helpful to constrain the spatial patterns of remote sensing images to be generated at that resolution. These prior models can be used properly by multiple-point statistics (MPS), capable of extracting the intrinsic features of patterns from prior models such as training images, and copying them to the simulated regions using hard and soft conditional data, or even without any conditional data. However, because traditional MPS methods based on linear dimensionality reduction are not suitable to deal with nonlinear data, and isometric mapping (ISOMAP) can reduce the dimensionality of nonlinear data effectively, this paper presents a sequential simulation framework for generating super-resolution remote sensing images using ISOMAP and MPS. Using four different examples, it is demonstrated that the structural characteristics of super-resolution reconstruction of remote sensing images using this method, are similar to those of training images.

A nearest neighbor multiple-point statistics method for fast geological modeling

We seek to develop a methodology enabling fast geostatistical simulations honoring both geophysical data and a complex prior model. Particularly, we consider a multiple-point statistics (MPS) framework in which a training image (TI) describes the available prior knowledge. Accurate posterior sampling is then possible by using a so-called extended Metropolis algorithm in which proposals are drawn from the prior using sequential geostatistical resampling. Such a Markov chain Monte Carlo (MCMC) algorithm will eventually locate and sample proportionally to the posterior distribution, however, it is often exceedingly slow and typically demands millions of MCMC iterations before the posterior is sampled sufficiently. We are developing a methodology in which the MPS simulation is built up iteratively pixel-by-pixel starting from an empty grid. At each pixel, multiple proposals are generated using an MPS algorithm and the proposals are accepted proportionally to the likelihood considering conditioning data in terms of linear averages (for instance geophysical data). The likelihood function is generally intractable as it depends on the pixels that have not yet been sampled. We approximate the likelihood function using a Gaussian model in which the posterior mean and covariance are updated sequentially as the simulation builds up. The posterior statistics are approximated by running the algorithm multiple times (sequentially or in parallel). Considering crosshole first-arrival ground-penetrating radar data, we assess the accuracy of our methodology both for multi-Gaussian priors for which analytical posteriors are available and for more complex training images against the extended Metropolis method. Our approach is inherently approximate due to the use of a finite training image, a finite number of candidates for each pixel and the need to approximate intractable likelihood functions. Nevertheless, preliminary results are promising as this method allows directly obtaining a reasonable estimation at a reduced computational cost compared to MCMC.

I. Principles, Basic Equations and Analytical Solutions.- Review of stochastic theory of transport in groundwater flow.- On the identification of the parameters of groundwater mass transport.- Heat and mass transfer in unsaturated porous media with application to thermal energy storage.- Basin-scale transport of dissolved species in groundwater.- Variable density fluid flow in the brackish transition zone between fresh and saline groundwater.- Stationary principles for flow and transport in aquifers.- II. Modelling Flow and Transport in Porous-Like Media.- An overwiew of groundwater modelling.- Incorporating assurance into groundwater quality management models.- The eigenvalues approach for solving linear groundwater flow problems.- New method for diffusive transport.- Comparison of fast equation solvers for groundwater flow problems.- Modelling flow and transport through porous media in vector computers.- Eulerian-Lagrangian method for solving transport in aquifers.- Hydrodynamic dispersion in model porous media.- The random walk method in pollutant transport simulation.- Modeling of solute transport with the random walk method.- III. Modelling Flow and Transport in Fractured Media.- Modelling of flow through fractured rocks: geostatistical generation of fractures networks. Stress-flow relationship in fractures.- Combined seismic and hydraulic method of modeling flow in fractured low permeability rocks.- Stochastic continuum representation of fractured rock permeability as an alternative to the REV and fracture network concepts.- Flow in three-dimensional fracture networks using a discrete approach.- MINC: an approach for analyzing transport in strongly heterogeneous systems.- A stochastic particle transport model based on directional statistics of flow through fracture networks.- IV. Consideration of Fluid-Solid Phase Interactions and Heterogeneities in Modelling.- Advances in modelling water-rock interaction in aquifers.- Parameters for modelling the transport of cadmium as influenced by the chemical properties of groundwater and aquifer material.- Micro-scale modelling in the study of plume evolution in heterogeneous media.- Modelling the increasing dispersivity with FE transport models using the multilayer concept.- Random-walk method to simulate pollutant transport in alluvial aquifers or fractured rocks.- V. Multiphase Flow and Transport Modelling.- Advances in modeling of water in the unsaturated zone.- FLuid-mechanical aspects of the migration of chemicals in fractured media.- About the numerical analysis of dynamics in multi-component-continua.- VI. Aquifer Parameter Identification by Models.- State-of-the-art of the inverse problem applied to the flow and solute transport equations.- Estimation of spatial covariance structures with application to hydrological, hydro-chemical and isotopic data from aquifers: state-of-the-art and adjoint state maximum likelihood cross-validation methods.- Bayesian identification of steady-state (anisotropic) groundwater flow models.- VII. Data Gathering and Utilization of Models.- Present limitations and perspectives on modelling pollution problems in aquifers.- Advances in the assessment of data worth for engineering decision analysis in groundwater contamination problems.- A comparative analysis of mathematical mass transport codes for groundwater pollution studies.- New approaches and applications in subsurface flow modeling: 3-D finite element analysis of dewatering for an electro-nuclear plant.- Large scale 3-D groundwater flow modeling in highly heterogeneous geologic medium.- Gathering of data for modelling.- Present status of coastal aquifer modelling: short review.- Conclusions.- Transport modeling.- Solution methods in groundwater flow and transport.- Transport processes in groundwater: chemical and biological aspects.- Application of models.- Author Index.- List of Participants.

Extensive outcrop data or photographs of present day depositions or even simple drawings from expert geologists contain precious structural information about spatial continuity that is beyond the present tools of geostatistics essentially limited to two-point statistics (histograms and covariances). A neural net can be learned to collect multiple point statistics from various training images, these statistics are then used to generate stochastic models conditioned to actual data. In petroleum applications, the methodology developed can be a substitute for objects based-algorithms when facies geometry and reservoir continuity are too complex to be modeled by simple object such as channels. The performance of the neural net approach is illustrated using training images of increasing complexity, and attempts at explaining that performance is provided in each case. The neural net builds local probability distributions for the facies types (categorical case) or for petrophysical properties (continuous case). These local probabilities include multiple point statistics learned from training images and are sampled with a Metropolis-Hastings sampler which also ensures reproduction of statistics coming from the actual subsurface data such as locally varying facies proportions. P. 321

The construction of training images (TIs) depicting the geological prior is one of the most critical step in multiple-point statistics. Geophysical techniques may be used to reduce the uncertainty in the understanding of prior geological scenarios. We developed a methodology to verify the consistency of geophysical data with independently-built TIs. Synthetic geophysical models built from TI scenarios are compared, using multidimensional scaling, with inverted models from field surveys to check if TIs are consistent with geophysical models. Then, the probability of each TI scenario is computed. A cluster analysis enables to determine which parameters used in building the TIs are most impacting the geophysical response. The methodology is tested using ERT to analyze TI scenarios in the Meuse River alluvial aquifer (Belgium)

The reconstruction of porous media is of great importance in predicting fluid transport properties, which are widely used in various fields such as catalysis, oil recovery, medicine and aging of building materials. The real three-dimensional structural data of porous media are helpful to describe the irregular topologic structures of porous media. By using multiple-point statistics (MPS) to extract the characteristics of real porous media acquired from micro computed tomography (micro-CT) scanning, the probabilities of structural characteristics of pore spaces are obtained first, and then reproduced in the reconstructed regions. One solution to overcome the anisotropy of training images is to use real 3D volume data as a training image (TI). The CPU cost and memory burden brought up by 3D simulations can be reduced greatly by selecting the optimal multiple-grid template size that is determined by the entropy of a TI. Moreover, both soft data and hard data are integrated in MPS simulation to improve the accuracy of reconstructed images. The variograms and permeabilities, computed by lattice Boltzmann method, of the reconstructed images and the target image obtained from real volume data are compared, showing that the structural characteristics of reconstructed porous media using our method are similar to those of real volume data.

Multiple-point geostatistical simulation has recently become popular in stochastic hydrogeology, primarily because of its capability to derive geologically reasonable patterns and multivariate distributions from a training image and conditioning training image data to multiple hard and soft data sources. This article resents, evaluates and contrasts the results of using multiplepoint geostatistical simulation for producing geologically realistic models of a Quaternary inner-alpine aquifer. Borehole data, expert-designed geological profiles and training image data were subject to conditional simulation with the Single Normal Equation Simulation (SNESIM) algorithm, one of the most widely used Matrix Product State (MPS) algorithms. The sensitivity of model predictions to the training image and hard as well as soft data input was evaluated. Modeling results indicate that soft conditioning in MPS is a convenient and efficient way for integrating secondary data such as geological expert drawings. Keywords: Multiple point statistics; Secondary data integration; Aquifer heterogeneity; Quaternary over-deepened alpine basin