Numerical Calculation of Equivalent Cell Permeability Tensors for General Quadrilateral Control Volumes

A new method for upscaling fine scale permeability fields to general quadrilateral-shaped coarse cells is presented. The procedure, referred to as the conforming scale up method, applies a triangle-based finite element technique, capable of accurately resolving both the coarse cell geometry and the subgrid heterogeneity, to the solution of the local fine scale problem. An appropriate averaging of this solution provides the equivalent permeability tensor for the coarse scale quadrilateral cell. The general level of accuracy of the technique is demonstrated through application to a number of flow problems. The real strength of the conforming scale up method is demonstrated when the method is applied in conjunction with a flow-based gridding technique. In this case, the approach is shown to provide results that are significantly more accurate than those obtained using standard techniques.     

Multiscale mixed/mimetic methods on corner-point grids

Multiscale simulation is a promising approach to facilitate direct simulation of large and complex grid models for highly heterogeneous petroleum reservoirs. Unlike traditional simulation, approaches based on upscaling/downscaling, multiscale methods seek to solve the full flow problem by incorporating subscale heterogeneities into local discrete approximation spaces. We consider a multiscale formulation based on a hierarchical grid approach, where basis functions with subgrid resolution are computed numerically to correctly and accurately account for subscale variations from an underlying (fine-scale) geomodel when solving the global flow equations on a coarse grid. By using multiscale basis functions to discretise the global flow equations on a (moderately sized) coarse grid, one can retain the efficiency of an upscaling method and, at the same time, produce detailed and conservative velocity fields on the underlying fine grid. For pressure equations, the multiscale mixed finite-element method (MsMFEM) has been shown to be a particularly versatile approach. In this paper, we extend the method to corner-point grids, which is the industry standard for modelling complex reservoir geology. To implement MsMFEM, one needs a discretisation method for solving local flow problems on the underlying fine grids. In principle, any stable and conservative method can be used. Here, we use a mimetic discretisation, which is a generalisation of mixed finite elements that gives a discrete inner product, allows for polyhedral elements, and can (easily) be extended to curved grid faces. The coarse grid can, in principle, be any partition of the subgrid, where each coarse block is a connected collection of subgrid cells. However, we argue that, when generating coarse grids, one should follow certain simple guidelines to achieve improved accuracy. We discuss partitioning in both index space and physical space and suggest simple processing techniques. The versatility and accuracy of the new multiscale mixed methodology is demonstrated on two corner-point models: a small Y-shaped sector model and a complex model of a layered sedimentary bed. A variety of coarse grids, both violating and obeying the above mentioned guidelines, are employed. The MsMFEM solutions are compared with a reference solution obtained by direct simulation on the subgrid.     

A locally conservative variational multiscale method for the simulation of porous media flow with multiscale source terms

We present a variational multiscale mixed finite element method for the solution of Darcy flow in porous media, in which both the permeability field and the source term display a multiscale character. The formulation is based on a multiscale split of the solution into coarse and subgrid scales. This decomposition is invoked in a variational setting that leads to a rigorous definition of a (global) coarse problem and a set of (local) subgrid problems. One of the key issues for the success of the method is the proper definition of the boundary conditions for the localization of the subgrid problems. We identify a weak compatibility condition that allows for subgrid communication across element interfaces, a feature that turns out to be essential for obtaining high-quality solutions. We also remove the singularities due to concentrated sources from the coarse-scale problem by introducing additional multiscale basis functions, based on a decomposition of fine-scale source terms into coarse and deviatoric components. The method is locally conservative and employs a low-order approximation of pressure and velocity at both scales. We illustrate the performance of the method on several synthetic cases and conclude that the method is able to capture the global and local flow patterns accurately.     

An upscaling method for one‐phase flow in heterogeneous reservoirs. A weighted output least squares (WOLS) approach

In this paper we study the problem of determining the effective permeability on a coarse scale level of problems with strongly varying and discontinuous coefficients defined on a fine scale. The upscaled permeability is defined as the solution of an optimization problem, where the difference between the fine scale and the coarse scale velocity field is minimized. We show that it is not necessary to solve the fine scale pressure equation in order to minimize the associated cost‐functional. Furthermore, we derive a simple technique for computing the derivatives of the cost‐functional needed in the fix‐point iteration used to compute the optimal permeability on the coarse mesh. Finally, the method is illustrated by several analytical examples and numerical experiments. This revised version was published online in July 2006 with corrections to the Cover Date.     

Nested gridding and streamline-based simulation for fast reservoir performance prediction

Detailed reservoir models routinely contain 10⁶–10⁸ grid blocks. These models often cannot be used directly in a reservoir simulation because of the time and memory required for solving the pressure grid on the fine grid. We propose a nested gridding technique that efficiently obtains an approximate solution for the pressure field. The domain is divided into a series of coarse blocks, each containing several fine cells. Effective mobilities are computed for each coarse grid block and the pressure is then found on the coarse scale. The pressure field within each coarse block is computed using flux boundary conditions obtained from the coarse pressure solution. Streamline-based simulation is used to move saturations forward in time. We test the method for a series of example waterflood problems and demonstrate that the method can give accurate estimates of oil production for large 3D models significantly faster than direct simulation using streamlines on the fine grid, making the method overall approximately up to 1,000 times faster than direct conventional simulation.     

Constraint energy minimizing generalized multiscale finite element method in the mixed formulation

This paper presents a novel mass-conservative mixed multiscale method for solving flow equations in heterogeneous porous media. The media properties (the permeability) contain multiple scales and high contrast. The proposed method solves the flow equation in a mixed formulation on a coarse grid by constructing multiscale basis functions. The resulting velocity field is mass-conservative on the fine grid. Our main goal is to obtain first-order convergence in terms of the mesh size which is independent of local contrast. This is achieved, first, by constructing some auxiliary spaces, which contain global information that cannot be localized, in general. This is built on our previous work on the generalized multiscale finite element method (GMsFEM). In the auxiliary space, multiscale basis functions corresponding to small (contrast-dependent) eigenvalues are selected. These basis functions represent the high-conductivity channels (which connect the boundaries of a coarse block). Next, we solve local problems to construct multiscale basis functions for the velocity field. These local problems are formulated in the oversampled domain, taking into account some constraints with respect to auxiliary spaces. The latter allows fast spatial decay of local solutions and, thus, allows taking smaller oversampled regions. The number of basis functions depends on small eigenvalues of the local spectral problems. Moreover, multiscale pressure basis functions are needed in constructing the velocity space. Our multiscale spaces have a minimal dimension, which is needed to avoid contrast dependence in the convergence. The method’s convergence requires an oversampling of several layers. We present an analysis of our approach. Our numerical results confirm that the convergence rate is first order with respect to the mesh size and independent of the contrast.     

Mixed multiscale finite element methods using approximate global information based on partial upscaling

The use of limited global information in multiscale simulations is needed when there is no scale separation. Previous approaches entail fine-scale simulations in the computation of the global information. The computation of the global information is expensive. In this paper, we propose the use of approximate global information based on partial upscaling. A requirement for partial homogenization is to capture long-range (non-local) effects present in the fine-scale solution, while homogenizing some of the smallest scales. The local information at these smallest scales is captured in the computation of basis functions. Thus, the proposed approach allows us to avoid the computations at the scales that can be homogenized. This results in coarser problems for the computation of global fields. We analyze the convergence of the proposed method. Mathematical formalism is introduced, which allows estimating the errors due to small scales that are homogenized. The proposed method is applied to simulate two-phase flows in heterogeneous porous media. Numerical results are presented for various permeability fields, including those generated using two-point correlation functions and channelized permeability fields from the SPE Comparative Project (Christie and Blunt, SPE Reserv Evalu Eng 4:308–317, 2001). We consider simple cases where one can identify the scales that can be homogenized. For more general cases, we suggest the use of upscaling on the coarse grid with the size smaller than the target coarse grid where multiscale basis functions are constructed. This intermediate coarse grid renders a partially upscaled solution that contains essential non-local information. Numerical examples demonstrate that the use of approximate global information provides better accuracy than purely local multiscale methods.     

Nearly perfectly matched layer boundary conditions for operator upscaling of the acoustic wave equation

Acoustic imaging and sensor modeling are processes that require repeated solution of the acoustic wave equation. Solution of the wave equation can be computationally expensive and memory intensive for large simulation domains. One scheme for speeding up solution of the wave equation is the operator-based upscaling method. The algorithm proceeds in two steps. First, the wave equation is solved for fine grid unknowns internal to coarse blocks assuming the coarse blocks do not need to communicate with neighboring blocks in parallel. Second, these fine grid solutions are used to form a new problem which is solved on the coarse grid. Accurate and efficient wave propagation schemes also must avoid artificial reflections off of the computational domain edges. One popular method for preventing artificial reflections is the nearly perfectly matched layer (NPML) method. In this paper, we discuss applying NPML to operator upscaling for the wave equation. We show that although we only apply NPML to the first step of this two step algorithm (directly affecting the fine grid unknowns only), we still see a significant reduction of reflections back into the domain. We describe three numerical experiments (one homogeneous medium experiment and two heterogeneous media examples) in which we validate that the solution of the wave equation exponentially decays in the NPML regions. Numerical experiments of acoustic wave propagation in two dimensions with a reasonable absorbing layer thickness resulted in a maximum pressure reflection of 3–8%. While the coarse grid acceleration is not explicitly damped in our algorithm, the tight coupling between the two steps of the algorithm results in only 0.1–1% of acceleration reflecting back into the computational domain.     

Risk Management for Petroleum Reservoir Production: A Simulation-Based Study of Prediction

We consider numerical solutions of the Darcy and Buckley–Leverett equations for flow in porous media. These solutions depend on a realization of a random field that describes the reservoir permeability. The main content of this paper is to formulate and analyze a probability model for the numerical coarse grid solution error. We explore the extent to which the coarse grid oil production rate is sufficient to predict future oil production rates. We find that very early oil production data is sufficient to reduce the prediction error in oil production by about 30%, relative to the prior probability prediction.     

An application of Green’s function technique for computing well inflow without radial flow assumption

Well modeling plays an important role in numerical reservoir simulation. The main difficulty in well modeling is the difference in scale between the wellbore radius and well gridblock dimension used in the simulation. The Peaceman equation is widely used in reservoir simulation to match gridblock pressure to the local solution of the diffusivity equation describing the flow near the well. However, this approach was developed under the assumption of radial flow. At the same time, the well inflow equation can be solved within the Green’s function (GF) formalism which allows the solution to be obtained without the assumption of radial flow. The GF solution can be presented as a series over the eigenvalues of the Laplace differential operator. However, this series converges conditionally and its direct summation is time-consuming. In Posvyanskii et al. (2008), a method for fast summation of such a series was proposed and successfully applied for analyzing the pressure build up curves. In this paper, we adopt the same technique for calculating the well indices for horizontal, slanted and partially penetrated wells. Additionally, the role of different boundary conditions is considered. The semi-analytical expressions for well indices are obtained and compared to the solution of the Peaceman equation. It is shown that in some cases, the difference between these solutions can be significant. The use of the obtained expression in numerical flow simulation allows well inflow to be modeled with high accuracy even on a coarse grid.     

Uncertainty quantification of two-phase flow problems via measure theory and the generalized multiscale finite element method

The problem of multiphase phase flow in heterogeneous subsurface porous media is one involving many uncertainties. In particular, the permeability of the medium is an important aspect of the model that is inherently uncertain. Properly quantifying these uncertainties is essential in order to make reliable probabilistic-based predictions and future decisions. In this work, a measure-theoretic framework is employed to quantify uncertainties in a two-phase subsurface flow model in high-contrast media. Given uncertain saturation data from observation wells, the stochastic inverse problem is solved numerically in order to obtain a probability measure on the space of unknown permeability parameters characterizing the two-phase flow. As solving the stochastic inverse problem requires a number of forward model solves, we also incorporate the use of a conservative version of the generalized multiscale finite element method for added efficiency. The parameter-space probability measure is used in order to make predictions of saturation values where measurements are not available, and to validate the effectiveness of the proposed approach in the context of fine and coarse model solves. A number of numerical examples are offered to illustrate the measure-theoretic methodology for solving the stochastic inverse problem using both fine and coarse solution schemes.     

A multiscale method for distributed parameter estimation with application to reservoir history matching

A method for multiscale parameter estimation with application to reservoir history matching is presented. Starting from a given fine-scale model, coarser models are generated using a global upscaling technique where the coarse models are tuned to match the solution of the fine model. Conditioning to dynamic data is done by history-matching the coarse model. Using consistently the same resolution both for the forward and inverse problems, this model is successively refined using a combination of downscaling and history matching until model-matching dynamic data are obtained at the finest scale. Large-scale corrections are obtained using fast models, which, combined with a downscaling procedure, provide a better initial model for the final adjustment on the fine scale. The result is thus a series of models with different resolution, all matching history as good as possible with this grid. Numerical examples show that this method may significantly reduce the computational effort and/or improve the quality of the solution when achieving a fine-scale match as compared to history-matching directly on the fine scale.     

Addressing Conditioning Data in Multiple-Point Statistics Simulation Algorithms Based on a Multiple Grid Approach

Multiple-point statistics (MPS) allows simulations reproducing structures of a conceptual model given by a training image (TI) to be generated within a stochastic framework. In classical implementations, fixed search templates are used to retrieve the patterns from the TI. A multiple grid approach allows the large-scale structures present in the TI to be captured, while keeping the search template small. The technique consists in decomposing the simulation grid into several grid levels: One grid level is composed of each second node of the grid level one rank finer. Then each grid level is successively simulated by using the corresponding rescaled search template from the coarse level to the fine level (the simulation grid itself). For a conditional simulation, a basic method (as in snesim) to honor the hard data consists in assigning the data to the closest nodes of the current grid level before simulating it. In this paper, another method (implemented in impala) that consists in assigning the hard data to the closest nodes of the simulation grid (fine level), and then in spreading them up to the coarse grid by using simulations based on the MPS inferred from the TI is presented in detail. We study the effect of conditioning and show that the first method leads to systematic biases depending on the location of the conditioning data relative to the grid levels, whereas the second method allows for properly dealing with conditional simulations and a multiple grid approach.     

A comparison of multiscale methods for elliptic problems in porous media flow

We review and perform comparison studies for three recent multiscale methods for solving elliptic problems in porous media flow; the multiscale mixed finite-element method, the numerical subgrid upscaling method, and the multiscale finite-volume method. These methods are based on a hierarchical strategy, where the global flow equations are solved on a coarsened mesh only. However, for each method, the discrete formulation of the partial differential equations on the coarse mesh is designed in a particular fashion to account for the impact of heterogeneous subgrid structures of the porous medium. The three multiscale methods produce solutions that are mass conservative on the underlying fine mesh. The methods may therefore be viewed as efficient, approximate fine-scale solvers, i.e., as an inexpensive alternative to solving the elliptic problem on the fine mesh. In addition, the methods may be utilized as an alternative to upscaling, as they generate mass-conservative solutions on the coarse mesh. We therefore choose to also compare the multiscale methods with a state-of-the-art upscaling method – the adaptive local–global upscaling method, which may be viewed as a multiscale method when coupled with a mass-conservative downscaling procedure. We investigate the properties of all four methods through a series of numerical experiments designed to reveal differences with regard to accuracy and robustness. The numerical experiments reveal particular problems with some of the methods, and these will be discussed in detail along with possible solutions. Next, we comment on implementational aspects and perform a simple analysis and comparison of the computational costs associated with each of the methods. Finally, we apply the three multiscale methods to a dynamic two-phase flow case and demonstrate that high efficiency and accurate results can be obtained when the subgrid computations are made part of a preprocessing step and not updated, or updated infrequently, throughout the simulation. The research is funded by the Research Council of Norway under grant nos. 152732 and 158908.     

Use of Border Regions for Improved Permeability Upscaling

A procedure for the improved calculation of upscaled grid block permeability tensors on Cartesian grids is described and applied. The method entails the use of a border region of fine-scale cells surrounding the coarse block for which the upscaled permeability is to be computed. The implementation allows for the use of full-tensor permeability fields on the fine and coarse scales. Either periodic or pressure–no flow boundary conditions are imposed over the extended local domain (target block plus border regions) though averaged quantities, used to compute the upscaled permeability tensor, are computed only over the target block region. Flow and transport results using this procedure are compared to those from standard methods for different types of geological and simulation models. Improvement using the new approach is consistently observed for the cases considered, though the degree of improvement varies for different models and flow quantities.     

Multiscale finite-volume method for density-driven flow in porous media

The multiscale finite-volume (MSFV) method has been developed to solve multiphase flow problems on large and highly heterogeneous domains efficiently. It employs an auxiliary coarse grid, together with its dual, to define and solve a coarse-scale pressure problem. A set of basis functions, which are local solutions on dual cells, is used to interpolate the coarse-grid pressure and obtain an approximate fine-scale pressure distribution. However, if flow takes place in presence of gravity (or capillarity), the basis functions are not good interpolators. To treat this case correctly, a correction function is added to the basis function interpolated pressure. This function, which is similar to a supplementary basis function independent of the coarse-scale pressure, allows for a very accurate fine-scale approximation. In the coarse-scale pressure equation, it appears as an additional source term and can be regarded as a local correction to the coarse-scale operator: It modifies the fluxes across the coarse-cell interfaces defined by the basis functions. Given the closure assumption that localizes the pressure problem in a dual cell, the derivation of the local problem that defines the correction function is exact, and no additional hypothesis is needed. Therefore, as in the original MSFV method, the only closure approximation is the localization assumption. The numerical experiments performed for density-driven flow problems (counter-current flow and lock exchange) demonstrate excellent agreement between the MSFV solutions and the corresponding fine-scale reference solutions.     

油藏精细地质模型网格粗化算法及其效果

在前人研究基础上, 根据DP(Dykstra-Parsons)系数能定量评价储层非均质性, 微网格块的渗透率值粗化后, 其等效渗透率的上、下限(C_min、C_max)能反映渗透率的各向异性的特点, 提出了一种运算速度快和相对有效的网格粗化算法。该算法能考虑到储层非均质性对不同方向渗透率值的影响, 且求解过程相对简单。应用该方法对鄂尔多斯盆地中部某油藏的陆相储层的精细地质模型进行了网格粗化计算, 然后在粗化后的模型上进行油藏数值模拟研究, 同时针对研究区地质背景和产出流体微可压缩的物性特征, 首次利用流线模拟器对精细地质模型进行了油藏数值模拟研究, 并以此结果为标准, 对该网格粗化算法时效性进行了系统评价。分析表明, 该算法具有较快的计算速度和较高的可靠性, 是解决储层非均质强、物性差的陆相成因油藏精细油藏数值模拟的一种行之有效的手段。