Parallel Multiple-Point Statistics Algorithm Based on List and Tree Structures

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>

Building a training image with Digital Outcrop Models

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.

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.

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

Resource estimation of mineral deposits requires spatial modelling of orebody boundaries based on a set of exploration borehole data. Given lateritic bauxite deposits, the spacing between the boreholes is often determined based on the grade continuity. As a result, the selected drill spacing might not capture the underlying (true) lateral variability apparent in the orebody boundaries. The purpose of this study is to investigate and address the limitations imposed by such problems in lateritic metal deposits through multiple-point statistics (MPS) framework. Rather than relying on a semivariogram model, we obtain the required structural information from the footwall topographies exposed after previous mining operations. The investigation utilising the MPS was carried out using the Direct Sampling (DS) MPS algorithm. Two historical mine areas along with their mined-out surfaces and ground penetrating radar surveys were incorporated as a bivariate training image to perform the MPS simulations. In addition, geostatistical simulations using the Turning Bands method were also performed to make the comparison against the MPS results. The performances were assessed using several statistical indicators including higher-order spatial cumulants. The results have shown that the DS can satisfactorily simulate the orebody boundaries by using prior information from the previously mined-out areas.

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.

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.

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)

Quantitative data on geobodies are crucial for reservoir modelling. Although abundant quantitative data are available in the literature for siliciclastic depositional systems, equivalent data for carbonate systems are scarce. In this paper we introduce a new approach to the management of quantitative data on carbonate geobodies which is based on a hierarchical classification scheme. The classes to which a carbonate geobody are assigned are: (1) depo‐time (i.e. geological age); (2) depo‐system (i.e. type of carbonate platform); (3) depo‐zone (i.e. facies belt or zone); (4) depo‐shape (i.e. geometry of the geological body); (5) depo‐element (i.e. architectural elements present); and (6) depo‐facies (litho‐ and biofacies). This hierarchical classification is complemented by a set of rules for modifying depo‐shapes which refer to their spatial distribution and patterns of interaction.Based on this classification, an extensive database has been developed which can be used for 3D reservoir modelling. The database holds more than 600 case studies from outcrop analogues and the subsurface and also from satellite images of modern carbonate settings. The database can be used as the basis for a new workflow for reservoir modelling which uses multiple‐point statistics (MPS). MPS makes use of training images to capture and reproduce facies patterns and geometries during stochastic simulations.The application of this new approach is demonstrated by modelling a Cretaceous outcrop reservoir analogue from southern France. The use of MPS allows the generation of geologically realistic and complex facies distributions in the model based on the simplified training images.

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

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