The use of connectivity clusters in stochastic modelling of highly heterogeneous porous media

Stochastic geostatistical techniques are essential tools for groundwater flow and transport modelling in highly heterogeneous media. Typically, these techniques require massive numbers of realizations to accurately simulate the high variability and account for the uncertainty. These massive numbers of realizations imposed several constraints on the stochastic techniques (e.g. increasing the computational effort, limiting the domain size, grid resolution, time step and convergence issues). Understanding the connectivity of the subsurface layers gives an opportunity to overcome these constraints. This research presents a sampling framework to reduce the number of the required Monte Carlo realizations utilizing the connectivity properties of the hydraulic conductivity distributions in a three-dimensional domain. Different geostatistical distributions were tested in this study including exponential distribution with the Turning Bands (TBM) algorithm and spherical distribution using Sequential Gaussian Simulation (SGSIM). It is found that the total connected fraction of the largest clusters and its tortuosity are highly correlated with the percentage of mass arrival and the first arrival quantiles at different control planes. Applying different sampling techniques together with several indicators suggested that a compact sample representing only 10% of the total number of realizations can be used to produce results that are close to the results of the full set of realizations. Also, the proposed sampling techniques specially utilizing the low conductivity clustering show very promising results in terms of matching the full range of realizations. Finally, the size of selected clusters relative to domain size significantly affects transport characteristics and the connectivity indicators.     

Assessing the risk of soil vulnerability to wind erosion through conditional simulation of soil water content in Sistan plain, Iran

Surface soil water content (SWC) is one of the key factors controlling wind erosion in Sistan plain, southeast of Iran. Knowledge of the spatial variability of surface SWC is then important to identify high-risk areas over the region. Sequential Gaussian simulation (SGSIM) is used to produce a series of equiprobable models of SWC spatial distribution across the study area. The simulated realizations are used to model the uncertainty attached to the surface SWC estimates through producing a probability map of not exceeding a specified critical threshold when soil becomes vulnerable to wind erosion. The results show that SGSIM is a suitable approach for modelling SWC uncertainty, generating realistic representations of the spatial distribution of SWC that honour the sample data and reproduce the sample semivariogram model. The uncertainty model obtained using SGSIM is compared with the model achieved through sequential indicator simulation (SISIM). According to accuracy plots, goodness statistics and probability interval width plots, SGSIM performs better for modelling local uncertainty than SISIM. Sequential simulation provided a probabilistic approach to assess the risk that SWC does not exceed a critical threshold that might cause soil vulnerability to wind erosion. The resulted risk map can be used in decision-making to delineate “vulnerable” areas where a treatment is needed.     

Gradual Deformation and Iterative Calibration of Gaussian-Related Stochastic Models

This paper describes a new method for gradually deforming realizations of Gaussian-related stochastic models while preserving their spatial variability. This method consists in building a stochastic process whose state space is the ensemble of the realizations of a spatial stochastic model. In particular, a stochastic process, built by combining independent Gaussian random functions, is proposed to perform the gradual deformation of realizations. Then, the gradual deformation algorithm is coupled with an optimization algorithm to calibrate realizations of stochastic models to nonlinear data. The method is applied to calibrate a continuous and a discrete synthetic permeability fields to well-test pressure data. The examples illustrate the efficiency of the proposed method. Furthermore, we present some extensions of this method (multidimensional gradual deformation, gradual deformation with respect to structural parameters, and local gradual deformation) that are useful in practice. Although the method described in this paper is operational only in the Gaussian framework (e.g., lognormal model, truncated Gaussian model, etc.), the idea of gradually deforming realizations through a stochastic process remains general and therefore promising even for calibrating non-Gaussian models.     

应用随机反演改善储层表征的过程(英文)

为了适当地完成储层表征的过程,一个有效的方法就是把现场所有可以利用的信息融合成一个一致性的模型。在实际生产中实现这种融合并非简单的任务,所以有必要运用如地震反演等特殊方法。应用地震反演可以使测井数据和地震数据的有效结合成为可能,并且可以得到一个模型,该模型在预测过程中可通过流体数字模拟来验证。地震反演可以通过多种方法进行,主要分为两大类:一类是确定性方法(其代表是回归反演和约束稀疏脉冲反演),另一类是随机方法(其代表是地质统计学反演)。在本次研究中,通过随机反演结果和确定性反演结果的对比展示了随机反演是如何改进储层表征过程的。事实上,随机反演,可以运用较高的采样率(和储层模型的网格大小相接近),来产生一个更可靠的模型。随机反演的另一个好处就是随机方法可产生一些基本的统计测量值来改进解释精度,并且在储层表征过程中能生成大量的实现,从而使储层模型的不确定性研究成为可能。     

Integration of Markov mesh models and data assimilation techniques in complex reservoirs

We present a methodology that allows conditioning the spatial distribution of geological and petrophysical properties of reservoir model realizations on available production data. The approach is fully consistent with modern concepts depicting natural reservoirs as composite media where the distribution of both lithological units (or facies) and associated attributes are modeled as stochastic processes of space. We represent the uncertain spatial distribution of the facies through a Markov mesh (MM) model, which allows describing complex and detailed facies geometries in a rigorous Bayesian framework. The latter is then embedded within a history matching workflow based on an iterative form of the ensemble Kalman filter (EnKF). We test the proposed methodology by way of a synthetic study characterized by the presence of two distinct facies. We analyze the accuracy and computational efficiency of our algorithm and its ability with respect to the standard EnKF to properly estimate model parameters and assess future reservoir production. We show the feasibility of integrating MM in a data assimilation scheme. Our methodology is conducive to a set of updated model realizations characterized by a realistic spatial distribution of facies and their log permeabilities. Model realizations updated through our proposed algorithm correctly capture the production dynamics.     

Combination of Dependent Realizations Within the Gradual Deformation Method

Gradual deformation is a parameterization method that reduces considerably the unknown parameter space of stochastic models. This method can be used in an iterative optimization procedure for constraining stochastic simulations to data that are complex, nonanalytical functions of the simulated variables. This method is based on the fact that linear combinations of multi-Gaussian random functions remain multi-Gaussian random functions. During the past few years, we developed the gradual deformation method by combining independent realizations. This paper investigates another alternative: the combination of dependent realizations. One of our motivations for combining dependent realizations was to improve the numerical stability of the gradual deformation method. Because of limitations both in the size of simulation grids and in the precision of simulation algorithms, numerical realizations of a stochastic model are never perfectly independent. It was shown that the accumulation of very small dependence between realizations might result in significant structural drift from the initial stochastic model. From the combination of random functions whose covariance and cross-covariance are proportional to each other, we derived a new formulation of the gradual deformation method that can explicitly take into account the numerical dependence between realizations. This new formulation allows us to reduce the structural deterioration during the iterative optimization. The problem of combining dependent realizations also arises when deforming conditional realizations of a stochastic model. As opposed to the combination of independent realizations, combining conditional realizations avoids the additional conditioning step during the optimization process. However, this procedure is limited to global deformations with fixed structural parameters.     

Scenario Reduction Applied to Geostatistical Simulations

Having a large number of geostatistical simulations of a mineral or petroleum deposit provides a better idea of its upside potential and downside risk; however, large numbers of simulated realizations of a deposit may pose computational difficulties in subsequent decision-making phases. Hence, depending on the specific case, there can be a need to select a representative subset of conditionally simulated deposit realizations. This paper examines and extends an approach developed by the stochastic optimization community based on stochastic mathematical programming with recourse and is discussed here in the context of mineral deposits while it is possibly suitable for other earth science applications. The approach is based on measuring the “distance” between simulations and the introduced distance measure between simulated realizations of a mineral deposit is based on the metal above a given set of cutoff grades while a pre-existing mine design is available. The approach is tested on 100 simulations of the Walker Lake data with promising results.     

A Distance-based Prior Model Parameterization for Constraining Solutions of Spatial Inverse Problems

Spatial inverse problems in the Earth Sciences are often ill-posed, requiring the specification of a prior model to constrain the nature of the inverse solutions. Otherwise, inverted model realizations lack geological realism. In spatial modeling, such prior model determines the spatial variability of the inverse solution, for example as constrained by a variogram, a Boolean model, or a training image-based model. In many cases, particularly in subsurface modeling, one lacks the amount of data to fully determine the nature of the spatial variability. For example, many different training images could be proposed for a given study area. Such alternative training images or scenarios relate to the different possible geological concepts each exhibiting a distinctive geological architecture. Many inverse methods rely on priors that represent a single subjectively chosen geological concept (a single variogram within a multi-Gaussian model or a single training image). This paper proposes a novel and practical parameterization of the prior model allowing several discrete choices of geological architectures within the prior. This method does not attempt to parameterize the possibly complex architectures by a set of model parameters. Instead, a large set of prior model realizations is provided in advance, by means of Monte Carlo simulation, where the training image is randomized. The parameterization is achieved by defining a metric space which accommodates this large set of model realizations. This metric space is equipped with a “similarity distance” function or a distance function that measures the similarity of geometry between any two model realizations relevant to the problem at hand. Through examples, inverse solutions can be efficiently found in this metric space using a simple stochastic search method.     

An Initial Guess for the Levenberg–Marquardt Algorithm for Conditioning a Stochastic Channel to Pressure Data

A standard procedure for conditioning a stochastic channel to well-test pressure data requires the minimization of an objective function. The Levenberg–Marquardt algorithm is a natural choice for minimization, but may suffer from slow convergence or converge to a local minimum which gives an unacceptable match of observed pressure data if a poor initial guess is used. In this work, we present a procedure to generate a good initial guess when the Levenberg–Marquardt algorithm is used to condition a stochastic channel to pressure data and well observations of channel facies, channel thickness, and channel top depth. This technique yields improved computational efficiency when the Levenberg–Marquardt method is used as the optimization procedure for generating realizations of the model by the randomized maximum likelihood method.     

Conditional Simulation with Patterns

An entirely new approach to stochastic simulation is proposed through the direct simulation of patterns. Unlike pixel-based (single grid cells) or object-based stochastic simulation, pattern-based simulation simulates by pasting patterns directly onto the simulation grid. A pattern is a multi-pixel configuration identifying a meaningful entity (a puzzle piece) of the underlying spatial continuity. The methodology relies on the use of a training image from which the pattern set (database) is extracted. The use of training images is not new. The concept of a training image is extensively used in simulating Markov random fields or for sequentially simulating structures using multiple-point statistics. Both these approaches rely on extracting statistics from the training image, then reproducing these statistics in multiple stochastic realizations, at the same time conditioning to any available data. The proposed approach does not rely, explicitly, on either a statistical or probabilistic methodology. Instead, a sequential simulation method is proposed that borrows heavily from the pattern recognition literature and simulates by pasting at each visited location along a random path a pattern that is compatible with the available local data and any previously simulated patterns. This paper discusses the various implementation details to accomplish this idea. Several 2D illustrative as well as realistic and complex 3D examples are presented to showcase the versatility of the proposed algorithm.     

The FFT Moving Average (FFT-MA) Generator: An Efficient Numerical Method for Generating and Conditioning Gaussian Simulations

A fast Fourier transform (FFT) moving average (FFT-MA) method for generating Gaussian stochastic processes is derived. Using discrete Fourier transforms makes the calculations easy and fast so that large random fields can be produced. On the other hand, the basic moving average frame allows us to uncouple the random numbers from the structural parameters (mean, variance, correlation length, ... ), but also to draw the randomness components in spatial domain. Such features impart great flexibility to the FFT-MA generator. For instance, changing only the random numbers gives distinct realizations all having the same covariance function. Similarly, several realizations can be built from the same random number set, but from different structural parameters. Integrating the FFT-MA generator into an optimization procedure provides a tool theoretically capable to determine the random numbers identifying the Gaussian field as well as the structural parameters from dynamic data. Moreover, all or only some of the random numbers can be perturbed so that realizations produced using the FFT-MA generator can be locally updated through an optimization process.     

Evaluation of the Reduction in Uncertainty Obtained by Conditioning a 3D Stochastic Channel to Multiwell Pressure Data

A stochastic channel embedded in a background facies is conditioned to data observed at wells. The background facies is a fixed rectangular box. The model parameters consist of geometric parameters that describe the shape, size, and location of the channel, and permeability and porosity in the channel and nonchannel facies. We extend methodology previously developed to condition a stochastic channel to well-test pressure data, and well observations of the channel thickness and the depth of the top of the channel. The main objective of this work is to characterize the reduction in uncertainty in channel model parameters and predicted reservoir performance that can be achieved by conditioning to well-test pressure data at one or more wells. Multiple conditional realizations of the geometric parameters and rock properties are generated to evaluate the uncertainty in model parameters. The ensemble of predictions of reservoir performance generated from the suite of realizations provides a Monte Carlo estimate of the uncertainty in future performance predictions. In addition, we provide some insight on how prior variances, data measurement errors, and sensitivity coefficients interact to determine the reduction in model parameters obtained by conditioning to pressure data and examine the value of active and observation well data in resolving model parameters.     

Compensating for estimation smoothing in kriging

Smoothing is a characteristic inherent to all minimum mean-square-error spatial estimators such as kriging. Cross-validation can be used to detect and model such smoothing. Inversion of the model produces a new estimator—compensated kriging. A numerical comparison based on an exhaustive permeability sampling of a 4-ft² slab of Berea Sandstone shows that the estimation surface generated by compensated kriging has properties intermediate between those generated by ordinary kriging and stochastic realizations resulting from simulated annealing and sequential Gaussian simulation. The frequency distribution is well reproduced by the compensated kriging surface, which also approximates the experimental semivariogram well—better than ordinary kriging, but not as well as stochastic realizations. Compensated kriging produces surfaces that are more accurate than stochastic realizations, but not as accurate as ordinary kriging.     

Parallel simulation of groundwater non-point source pollution using algebraic multigrid preconditioners

The simulation of non-point source pollution in agricultural basins is a computationally demanding process due to the large number of individual sources and potential pollution receptors (e.g., drinking water wells). In this study, we present an efficient computational framework for parallel simulation of diffuse pollution in such groundwater basins. To derive a highly detailed velocity field, we employed algebraic multigrid (AMG) preconditioners to solve the groundwater flow equation. We compare two variants of AMG implementations, the multilevel preconditioning provided by Trilinos and the BoomerAMG provided by HYPRE. We also perform a sensitivity analysis on the configuration of AMG methods to evaluate the application of these libraries to groundwater flow problems. For the transport simulation of diffuse contamination, we use the streamline approach, which decomposes the 3D transport problem into a large number of 1D problems that can be executed in parallel. The proposed framework is applied to a 2,600-km² groundwater basin in California discretized into a grid with over 11 million degrees of freedom. Using a Monte Carlo approach with 200 nitrate loading realizations at the aquifer surface, we perform a stochastic analysis to quantify nitrate breakthrough prediction uncertainty at over 1,500 wells due to random, temporally distributed nitrate loading. The results show that there is a significant time lag between loading and aquifer response at production wells. Generally, typical production wells respond after 5–50 years depending on well depth and screen length, while the prediction uncertainty for nitrate in individual wells is very large—approximately twice the drinking water limit for nitrate.     

History Matching Geostatistical Model Realizations Using a Geometrical Domain Based Parameterization Technique

Reservoir characterization needs the integration of various data through history matching, especially dynamic information such as production or 4D seismic data. Although reservoir heterogeneities are commonly generated using geostatistical models, random realizations cannot generally match observed dynamic data. To constrain model realizations to reproduce measured dynamic data, an optimization procedure may be applied in an attempt to minimize an objective function, which quantifies the mismatch between real and simulated data. Such assisted history matching methods require a parameterization of the geostatistical model to allow the updating of an initial model realization. However, there are only a few parameterization methods available to update geostatistical models in a way consistent with the underlying geostatistical properties. This paper presents a local domain parameterization technique that updates geostatistical realizations using assisted history matching. This technique allows us to locally change model realizations through the variation of geometrical domains whose geometry and size can be easily controlled and parameterized. This approach provides a new way to parameterize geostatistical realizations in order to improve history matching efficiency.     

Dynamic data integration for structural modeling: model screening approach using a distance-based model parameterization

This paper proposes a novel history-matching method where reservoir structure is inverted from dynamic fluid flow response. The proposed workflow consists of searching for models that match production history from a large set of prior structural model realizations. This prior set represents the reservoir structural uncertainty because of interpretation uncertainty on seismic sections. To make such a search effective, we introduce a parameter space defined with a “similarity distance” for accommodating this large set of realizations. The inverse solutions are found using a stochastic search method. Realistic reservoir examples are presented to prove the applicability of the proposed method.     

考虑污染源强随机变化的感潮河流环境容量优化

考虑污染源强随机变化和感潮河流潮周期内动态水文条件对水质的影响,建立了优化污染负荷分配的流域水质管理模型。模型以总的允许排污量最大为目标函数,流域的水质控制点达标为约束条件。假设排污量是服从对数正态分布的随机变量,并且以潮周期内水质达标的概率作为衡量控制点达标的依据。采用遗传算法对该随机规划模型进行求解。研究结果表明,污染负荷优化分配结果能够满足随机条件下的水质达标率要求,并且与传统的确定性线性规划模型的分配结果相比有着明显差别。同时证实了遗传算法能够有效地解决复杂的随机规划模型。     

Uncertainty in fractal dimension estimated from power spectra and variograms

The reliability of using fractal dimension (D) as a quantitative parameter to describe geological variables is dependent mainly on the accuracy of estimated D values from observed data. Two widely used methods for the estimation of fractal dimensions are based on fitting a fractal model to experimental variograms or power-spectra on a log-log plot. The purpose of this paper is to study the uncertainty in the fractal dimension estimated by these two methods. The results indicate that both spectrum and variogram methods result in biased estimates of the D value. Fractal dimension calculated by these two methods for the same data will be different unless the bias is properly corrected. The spectral method results in overestimated D values. The variogram method has a critical fractal dimension, below which overestimation occurs and above which underestimation occurs. On the bases of 36,000 simulated realizations we propose empirical formulae to correct for biases in the spectral and variogram estimated fractal dimension. Pitfalls in estimating fractal dimension from data contaminated by white noise or data having several fractal components have been identified and illustrated by simulated examples.