Cross-covariances and localization for EnKF in multiphase flow data assimilation

The ensemble Kalman filter has been successfully applied for data assimilation in very large models, including those in reservoir simulation and weather. Two problems become critical in a standard implementation of the ensemble Kalman filter, however, when the ensemble size is small. The first is that the ensemble approximation to cross-covariances of model and state variables to data can indicate the presence of correlations that are not real. These spurious correlations give rise to model or state variable updates in regions that should not be updated. The second problem is that the number of degrees of freedom in the ensemble is only as large as the size of the ensemble, so the assimilation of large amounts of precise, independent data is impossible. Localization of the Kalman gain is almost universal in the weather community, but applications of localization for the ensemble Kalman filter in porous media flow have been somewhat rare. It has been shown, however, that localization of updates to regions of non-zero sensitivity or regions of non-zero cross-covariance improves the performance of the EnKF when the ensemble size is small. Localization is necessary for assimilation of large amounts of independent data. The problem is to define appropriate localization functions for different types of data and different types of variables. We show that the knowledge of sensitivity alone is not sufficient for determination of the region of localization. The region depends also on the prior covariance for model variables and on the past history of data assimilation. Although the goal is to choose localization functions that are large enough to include the true region of non-zero cross-covariance, for EnKF applications, the choice of localization function needs to balance the harm done by spurious covariance resulting from small ensembles and the harm done by excluding real correlations. In this paper, we focus on the distance-based localization and provide insights for choosing suitable localization functions for data assimilation in multiphase flow problems. In practice, we conclude that it is reasonable to choose localization functions based on well patterns, that localization function should be larger than regions of non-zero sensitivity and should extend beyond a single well pattern.     

A probabilistic parametrization for geological uncertainty estimation using the ensemble Kalman filter (EnKF)

In the past years, many applications of history-matching methods in general and ensemble Kalman filter in particular have been proposed, especially in order to estimate fields that provide uncertainty in the stochastic process defined by the dynamical system of hydrocarbon recovery. Such fields can be permeability fields or porosity fields, but can also fields defined by the rock type (facies fields). The estimation of the boundaries of the geologic facies with ensemble Kalman filter (EnKF) was made, in different papers, with the aid of Gaussian random fields, which were truncated using various schemes and introduced in a history-matching process. In this paper, we estimate, in the frame of the EnKF process, the locations of three facies types that occur into a reservoir domain, with the property that each two could have a contact. The geological simulation model is a form of the general truncated plurigaussian method. The difference with other approaches consists in how the truncation scheme is introduced and in the observation operator of the facies types at the well locations. The projection from the continuous space of the Gaussian fields into the discrete space of the facies fields is realized through in an intermediary space (space with probabilities). This space connects the observation operator of the facies types at the well locations with the geological simulation model. We will test the model using a 2D reservoir which is connected with the EnKF method as a data assimilation technique. We will use different geostatistical properties for the Gaussian fields and different levels of the uncertainty introduced in the model parameters and also in the construction of the Gaussian fields.     

Assimilating spatially dense data for subsurface applications—balancing information and degrees of freedom

The degrees of freedom (DOF) in standard ensemble-based data assimilation is limited by the ensemble size. Successful assimilation of a data set with large information content (IC) therefore requires that the DOF is sufficiently large. A too small number of DOF with respect to the IC may result in ensemble collapse, or at least in unwarranted uncertainty reduction in the estimation results. In this situation, one has two options to restore a proper balance between the DOF and the IC: to increase the DOF or to decrease the IC. Spatially dense data sets typically have a large IC. Within subsurface applications, inverted time-lapse seismic data used for reservoir history matching is an example of a spatially dense data set. Such data are considered to have great potential due to their large IC, but they also contain errors that are challenging to characterize properly. The computational cost of running the forward simulations for reservoir history matching with any kind of data is large for field cases, such that a moderately large ensemble size is standard. Realization of the potential in seismic data for ensemble-based reservoir history matching is therefore not straightforward, not only because of the unknown character of the associated data errors, but also due to the imbalance between a large IC and a too small number of DOF. Distance-based localization is often applied to increase the DOF but is example specific and involves cumbersome implementation work. We consider methods to obtain a proper balance between the IC and the DOF when assimilating inverted seismic data for reservoir history matching. To decrease the IC, we consider three ways to reduce the influence of the data space; subspace pseudo inversion, data coarsening, and a novel way of performing front extraction. To increase the DOF, we consider coarse-scale simulation, which allows for an increase in the DOF by increasing the ensemble size without increasing the total computational cost. We also consider a combination of decreasing the IC and increasing the DOF by proposing a novel method consisting of a combination of data coarsening and coarse-scale simulation. The methods were compared on one small and one moderately large example with seismic bulk-velocity fields at four assimilation times as data. The size of the examples allows for calculation of a reference solution obtained with standard ensemble-based data assimilation methodology and an unrealistically large ensemble size. With the reference solution as the yardstick with which the quality of other methods are measured, we find that the novel method combining data coarsening and coarse-scale simulations gave the best results. With very restricted computational resources available, this was the only method that gave satisfactory results.     

Multiscale Parameterization of Petrophysical Properties for Efficient History-Matching

The prediction of fluid flows within hydrocarbon reservoirs requires the characterization of petrophysical properties. Such characterization is performed on the basis of geostatistics and history-matching; in short, a reservoir model is first randomly drawn, and then sequentially adjusted until it reproduces the available dynamic data. Two main concerns typical of the problem under consideration are the heterogeneity of rocks occurring at all scales and the use of data of distinct resolution levels. Therefore, referring to sequential Gaussian simulation, this paper proposes a new stochastic simulation method able to handle several scales for both continuous or discrete random fields. This method adds flexibility to history-matching as it boils down to the multiscale parameterization of reservoir models. In other words, reservoir models can be updated at either coarse or fine scales, or both. Parameterization adapts to the available data; the coarser the scale targeted, the smaller the number of unknown parameters, and the more efficient the history-matching process. This paper focuses on the use of variational optimization techniques driven by the gradual deformation method to vary reservoir models. Other data assimilation methods and perturbation processes could have been envisioned as well. Last, a numerical application case is presented in order to highlight the advantages of the proposed method for conditioning permeability models to dynamic data. For simplicity, we focus on two-scale processes. The coarse scale describes the variations in the trend while the fine scale characterizes local variations around the trend. The relationships between data resolution and parameterization are investigated.     

Impact of 3D Var GSI-ENKF hybrid data assimilation system

The hybrid two-way coupled 3DEnsVar assimilation system was tested with the NCMRWF global data assimilation forecasting system. At present, this system consists of T574L64 deterministic model and the grid-point statistical interpolation analysis scheme. In this experiment, the analysis system is modified with a two-way coupling with an 80 member Ensemble Kalman Filter of T254L64 resolution and runs are carried out in parallel to the operational system for the Indian summer monsoon season (June–September) for the year 2015 to study its impact. Both the assimilation systems are based on NCEP GFS system. It is found that hybrid assimilation marginally improved the quality of the forecasts of all variables over the deterministic 3D Var system, in terms of statistical skill scores and also in terms of circulation features. The impact of the hybrid system in prediction of extreme rainfall and cyclone track is discussed.     

Towards a hierarchical parametrization to address prior uncertainty in ensemble-based data assimilation

Ensemble-based methods are becoming popular assisted history matching techniques with a growing number of field applications. These methods use an ensemble of model realizations, typically constructed by means of geostatistics, to represent the prior uncertainty. The performance of the history matching is very dependent on the quality of the initial ensemble. However, there is a significant level of uncertainty in the parameters used to define the geostatistical model. From a Bayesian viewpoint, the uncertainty in the geostatistical modeling can be represented by a hyper-prior in a hierarchical formulation. This paper presents the first steps towards a general parametrization to address the problem of uncertainty in the prior modeling. The proposed parametrization is inspired in Gaussian mixtures, where the uncertainty in the prior mean and prior covariance is accounted by defining weights for combining multiple Gaussian ensembles, which are estimated during the data assimilation. The parametrization was successfully tested in a simple reservoir problem where the orientation of the major anisotropic direction of the permeability field was unknown.     

Combining sensitivities and prior information for covariance localization in the ensemble Kalman filter for petroleum reservoir applications

Sampling errors can severely degrade the reliability of estimates of conditional means and uncertainty quantification obtained by the application of the ensemble Kalman filter (EnKF) for data assimilation. A standard recommendation for reducing the spurious correlations and loss of variance due to sampling errors is to use covariance localization. In distance-based localization, the prior (forecast) covariance matrix at each data assimilation step is replaced with the Schur product of a correlation matrix with compact support and the forecast covariance matrix. The most important decision to be made in this localization procedure is the choice of the critical length(s) used to generate this correlation matrix. Here, we give a simple argument that the appropriate choice of critical length(s) should be based both on the underlying principal correlation length(s) of the geological model and the range of the sensitivity matrices. Based on this result, we implement a procedure for covariance localization and demonstrate with a set of distinctive reservoir history-matching examples that this procedure yields improved results over the standard EnKF implementation and over covariance localization with other choices of critical length.     

土壤水分同化系统的敏感性试验研究

利用1998年7月6日至8月9日青藏高原GAME-Tibet试验区MS3608站点的4cm、20cm和100cm的土壤水分观测数据同化SiB2模型输出的表层、根区和深层土壤水分,探讨了一个基于集合卡尔曼滤波和简单生物圈模型的单点土壤水分同化方案。分析和评价了集合大小、同化周期、模型误差、背景场误差以及观测误差对同化系统性能的影响。结果表明:①增加集合数目可以减小土壤水分同化系统的误差,但同时又降低了运行效率;②对于集合卡尔曼滤波,初始场的估计是否准确对同化系统性能影响不大;③模型误差和观测误差的准确估计可以提高土壤水分的估计精度;④利用数据同化的方法对土壤水分的估计有显著提高。     

Non-parametric Bayesian networks for parameter estimation in reservoir simulation: a graphical take on the ensemble Kalman filter (part I)

Reservoir simulation models are used both in the development of new fields and in developed fields where production forecasts are needed for investment decisions. When simulating a reservoir, one must account for the physical and chemical processes taking place in the subsurface. Rock and fluid properties are crucial when describing the flow in porous media. In this paper, the authors are concerned with estimating the permeability field of a reservoir. The problem of estimating model parameters such as permeability is often referred to as a history-matching problem in reservoir engineering. Currently, one of the most widely used methodologies which address the history-matching problem is the ensemble Kalman filter (EnKF). EnKF is a Monte Carlo implementation of the Bayesian update problem. Nevertheless, the EnKF methodology has certain limitations that encourage the search for an alternative method.For this reason, a new approach based on graphical models is proposed and studied. In particular, the graphical model chosen for this purpose is a dynamic non-parametric Bayesian network (NPBN). This is the first attempt to approach a history-matching problem in reservoir simulation using a NPBN-based method. A two-phase, two-dimensional flow model was implemented for a synthetic reservoir simulation exercise, and initial results are shown. The methods’ performances are evaluated and compared. This paper features a completely novel approach to history matching and constitutes only the first part (part I) of a more detailed investigation. For these reasons (novelty and incompleteness), many questions are left open and a number of recommendations are formulated, to be investigated in part II of the same paper.     

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.     

基于集合卡尔曼平滑算法的土壤水分同化

为研究观测资料稀少情况下土壤质地及有机质对土壤水分同化的影响,发展了集合卡尔曼平滑(Ensemble Kalman Smooth, EnKS)的土壤水分同化方案。利用黑河上游阿柔冻融观测站2008年6月1日至10月29日的观测数据,使用EnKS算法将表层土壤水分观测数据同化到简单生物圈模型(Simple Biosphere Model 2, SiB2)中,分析不同方案对土壤水分估计的影响,并与集合卡尔曼滤波算法(EnKF)的结果进行比较。研究结果表明,土壤质地和有机质对表层土壤水分模拟结果影响最大而对深层的影响相对较小;利用EnKF和EnKS算法同化表层土壤水分观测数据,均能够显著提高表层和根区土壤水分估计的精度,EnKS算法的精度略高于EnKF且所受土壤质地和有机质的影响小于EnKF;当观测数据稀少时,EnKS算法仍然可以得到较高精度的土壤水分估计。     

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

The use of the ensemble smoother (ES) instead of the ensemble Kalman filter increases the nonlinearity of the update step during data assimilation and the need for iterative assimilation methods. A previous version of the iterative ensemble smoother based on Gauss–Newton formulation was able to match data relatively well but only after a large number of iterations. A multiple data assimilation method (MDA) was generally more efficient for large problems but lacked ability to continue “iterating” if the data mismatch was too large. In this paper, we develop an efficient, iterative ensemble smoother algorithm based on the Levenberg–Marquardt (LM) method of regularizing the update direction and choosing the step length. The incorporation of the LM damping parameter reduces the tendency to add model roughness at early iterations when the update step is highly nonlinear, as it often is when all data are assimilated simultaneously. In addition, the ensemble approximation of the Hessian is modified in a way that simplifies computation and increases stability. We also report on a simplified algorithm in which the model mismatch term in the updating equation is neglected. We thoroughly evaluated the new algorithm based on the modified LM method, LM-ensemble randomized maximum likelihood (LM-EnRML), and the simplified version of the algorithm, LM-EnRML (approx), on three test cases. The first is a highly nonlinear single-variable problem for which results can be compared against the true conditional pdf. The second test case is a one-dimensional two-phase flow problem in which the permeability of 31 grid cells is uncertain. In this case, Markov chain Monte Carlo results are available for comparison with ensemble-based results. The third test case is the Brugge benchmark case with both 10 and 20 years of history. The efficiency and quality of results of the new algorithms were compared with the standard ES (without iteration), the ensemble-based Gauss–Newton formulation, the standard ensemble-based LM formulation, and the MDA. Because of the high level of nonlinearity, the standard ES performed poorly on all test cases. The MDA often performed well, especially at early iterations where the reduction in data mismatch was quite rapid. The best results, however, were always achieved with the new iterative ensemble smoother algorithms, LM-EnRML and LM-EnRML (approx).     

An efficient matrix-free algorithm for the ensemble Kalman filter

In this work, we present an efficient matrix-free ensemble Kalman filter (EnKF) algorithm for the assimilation of large data sets. The EnKF has increasingly become an essential tool for data assimilation of numerical models. It is an attractive assimilation method because it can evolve the model covariance matrix for a non-linear model, through the use of an ensemble of model states, and it is easy to implement for any numerical model. Nevertheless, the computational cost of the EnKF can increase significantly for cases involving the assimilation of large data sets. As more data become available for assimilation, a potential bottleneck in most EnKF algorithms involves the operation of the Kalman gain matrix. To reduce the complexity and cost of assimilating large data sets, a matrix-free EnKF algorithm is proposed. The algorithm uses an efficient matrix-free linear solver, based on the Sherman–Morrison formulas, to solve the implicit linear system within the Kalman gain matrix and compute the analysis. Numerical experiments with a two-dimensional shallow water model on the sphere are presented, where results show the matrix-free implementation outperforming an singular value decomposition-based implementation in computational time.     

History matching on the Imperial College fault model using parallel tempering

The history-matching inverse problem from petroleum engineering is analysed using the Imperial College fault model. This fault model produces a challenging inverse problem and is designed to show some of the problems which can occur whilst performing history-matching calculations on complicated geologies. It is shown that there can be multiple distinct geologies which match the history data. Furthermore, it is shown that the maximum-a-posteriori estimate does not correspond to the true geology in some cases. Both of these statements are corroborated via numerical examples where the parameter spaces are ?, ?³, ?⁷ and ?¹³. In addition, it is shown that the number of matches which agree with the data increases with dimension for these examples. It is also shown that the different matches can result in different reservoir management decision which, if incorrectly taken, would incur substantial financial penalties. All of these analyses are performed in a systematic manner, where it is shown that the standard algorithms can give a misleading answer. The history-matching problem is written as a minimisation problem, and it is shown that knowledge of all of the local minima is required. This presents significant computational issues as the resulting objective function is highly nonlinear, expensive to evaluate and multimodal. Previously used algorithms have been proved to be inadequate. Parallel tempering is a method which, if run for long enough, can find all the local minima. However, as the objective is expensive, a number of algorithm modifications had to be used to ensure convergence within a reasonable time. This new information is outlined in the paper. The algorithm as implemented produced results and new insights into this problem which were not suspected before. The results produced by this algorithm for the multimodal history-matching problem are superior to all other results of which we are aware. However, a considered amount of computation time was used within this paper, so this result does not infer that the algorithm cannot be improved upon. This algorithm not only produces good results but can be applied to all other history-matching problems. We have shown that this method provides a robust route of finding multiple local optima/solutions to the inverse problem, which is of considerable benefit to the petroleum industry. Furthermore, it is an entirely parallel algorithm which is becoming computationally feasible for other history-matching problems.     

Coarse-scale data assimilation as a generic alternative to localization

Coarse-scale data assimilation (DA) with large ensemble size is proposed as a robust alternative to standard DA with localization for reservoir history matching problems. With coarse-scale DA, the unknown property function associated with each ensemble member is upscaled to a grid significantly coarser than the original reservoir simulator grid. The grid coarsening is automatic, ensemble-specific and non-uniform. The selection of regions where the grid can be coarsened without introducing too large modelling errors is performed using a second-generation wavelet transform allowing for seamless handling of non-dyadic grids and inactive grid cells. An inexpensive local-local upscaling is performed on each ensemble member. A DA algorithm that restarts from initial time is utilized, which avoids the need for downscaling. Since the DA computational cost roughly equals the number of ensemble members times the cost of a single forward simulation, coarse-scale DA allows for a significant increase in the number of ensemble members at the same computational cost as standard DA with localization. Fixing the computational cost for both approaches, the quality of coarse-scale DA is compared to that of standard DA with localization (using state-of-the-art localization techniques) on examples spanning a large degree of variability. It is found that coarse-scale DA is more robust with respect to variation in example type than each of the localization techniques considered with standard DA. Although the paper is concerned with two spatial dimensions, coarse-scale DA is easily extendible to three spatial dimensions, where it is expected that its advantage with respect to standard DA with localization will increase.     

Gaussian Processes for history-matching: application to an unconventional gas reservoir

The process of reservoir history-matching is a costly task. Many available history-matching algorithms either fail to perform such a task or they require a large number of simulation runs. To overcome such struggles, we apply the Gaussian Process (GP) modeling technique to approximate the costly objective functions and to expedite finding the global optima. A GP model is a proxy, which is employed to model the input-output relationships by assuming a multi-Gaussian distribution on the output values. An infill criterion is used in conjunction with a GP model to help sequentially add the samples with potentially lower outputs. The IC fault model is used to compare the efficiency of GP-based optimization method with other typical optimization methods for minimizing the objective function. In this paper, we present the applicability of using a GP modeling approach for reservoir history-matching problems, which is exemplified by numerical analysis of production data from a horizontal multi-stage fractured tight gas condensate well. The results for the case that is studied here show a quick convergence to the lowest objective values in less than 100 simulations for this 20-dimensional problem. This amounts to an almost 10 times faster performance compared to the Differential Evolution (DE) algorithm that is also known to be a powerful optimization technique. The sensitivities are conducted to explain the performance of the GP-based optimization technique with various correlation functions.     

Ensemble filter methods with perturbed observations applied to nonlinear problems

The ensemble Kalman filter (EnKF) has been shown repeatedly to be an effective method for data assimilation in large-scale problems, including those in petroleum engineering. Data assimilation for multiphase flow in porous media is particularly difficult, however, because the relationships between model variables (e.g., permeability and porosity) and observations (e.g., water cut and gas–oil ratio) are highly nonlinear. Because of the linear approximation in the update step and the use of a limited number of realizations in an ensemble, the EnKF has a tendency to systematically underestimate the variance of the model variables. Various approaches have been suggested to reduce the magnitude of this problem, including the application of ensemble filter methods that do not require perturbations to the observed data. On the other hand, iterative least-squares data assimilation methods with perturbations of the observations have been shown to be fairly robust to nonlinearity in the data relationship. In this paper, we present EnKF with perturbed observations as a square root filter in an enlarged state space. By imposing second-order-exact sampling of the observation errors and independence constraints to eliminate the cross-covariance with predicted observation perturbations, we show that it is possible in linear problems to obtain results from EnKF with observation perturbations that are equivalent to ensemble square-root filter results. Results from a standard EnKF, EnKF with second-order-exact sampling of measurement errors that satisfy independence constraints (EnKF (SIC)), and an ensemble square-root filter (ETKF) are compared on various test problems with varying degrees of nonlinearity and dimensions. The first test problem is a simple one-variable quadratic model in which the nonlinearity of the observation operator is varied over a wide range by adjusting the magnitude of the coefficient of the quadratic term. The second problem has increased observation and model dimensions to test the EnKF (SIC) algorithm. The third test problem is a two-dimensional, two-phase reservoir flow problem in which permeability and porosity of every grid cell (5,000 model parameters) are unknown. The EnKF (SIC) and the mean-preserving ETKF (SRF) give similar results when applied to linear problems, and both are better than the standard EnKF. Although the ensemble methods are expected to handle the forecast step well in nonlinear problems, the estimates of the mean and the variance from the analysis step for all variants of ensemble filters are also surprisingly good, with little difference between ensemble methods when applied to nonlinear problems.     

Localization and regularization for iterative ensemble smoothers

Ensemble-based data assimilation methods have recently become popular for solving reservoir history matching problems, but because of the practical limitation on ensemble size, using localization is necessary to reduce the effect of sampling error and to increase the degrees of freedom for incorporating large amounts of data. Local analysis in the ensemble Kalman filter has been used extensively for very large models in numerical weather prediction. It scales well with the model size and the number of data and is easily parallelized. In the petroleum literature, however, iterative ensemble smoothers with localization of the Kalman gain matrix have become the state-of-the-art approach for ensemble-based history matching. By forming the Kalman gain matrix row-by-row, the analysis step can also be parallelized. Localization regularizes updates to model parameters and state variables using information on the distance between the these variables and the observations. The truncation of small singular values in truncated singular value decomposition (TSVD) at the analysis step provides another type of regularization by projecting updates to dominant directions spanned by the simulated data ensemble. Typically, the combined use of localization and TSVD is necessary for problems with large amounts of data. In this paper, we compare the performance of Kalman gain localization to two forms of local analysis for parameter estimation problems with nonlocal data. The effect of TSVD with different localization methods and with the use of iteration is also analyzed. With several examples, we show that good results can be obtained for all localization methods if the localization range is chosen appropriately, but the optimal localization range differs for the various methods. In general, for local analysis with observation taper, the optimal range is somewhat shorter than the optimal range for other localization methods. Although all methods gave equivalent results when used in an iterative ensemble smoother, the local analysis methods generally converged more quickly than Kalman gain localization when the amount of data is large compared to ensemble size.     

Shrinked (1???α) ensemble Kalman filter and α Gaussian mixture filter

State estimation in high dimensional systems remains a challenging part of real time analysis. The ensemble Kalman filter addresses this challenge by using Gaussian approximations constructed from a number of samples. This method has been a large success in many applications. Unfortunately, for some cases, Gaussian approximations are no longer valid, and the filter does not work so well. In this paper, we use the idea of the ensemble Kalman filter together with the more theoretically valid particle filter. We outline a Gaussian mixture approach based on shrinking the predicted samples to overcome sample degeneracy, while maintaining non-Gaussian nature. A tuning parameter determines the degree of shrinkage. The computational cost is similar to the ensemble Kalman filter. We compare several filtering methods on three different cases: a target tracking model, the Lorenz 40 model, and a reservoir simulation example conditional on seismic and electromagnetic data.