Olav Møyner

MRST-AD - an Open-Source Framework for Rapid Prototyping and Evaluation of Reservoir Simulation Problems

We present MRST-AD, a free, open-source framework written as part of the Matlab Reservoir Simulation Toolbox and designed to provide researchers with the means for rapid prototyping and experimentation for problems in reservoir simulation. The article outlines the design principles and programming techniques used and explains in detail the implementation of a full-featured, industry-standard black-oil model on unstructured grids. The resulting simulator has been thoroughly validated against a leading commercial simulator on benchmarks from the SPE Comparative Solution Projects, as well as on a real-field model (Voador, Brazil). We also show in detail how practitioners can easily extend the black-oil model with new constitutive relationships, or additional features such as polymer flooding, thermal and reactive effects, and immediately benefit from existing functionality such as constrained-pressure-residual (CPR) type preconditioning, sensitivities and adjoint-based gradients. Technically, MRST-AD combines three key features: (i) a highly vectorized scripting language that enables the user to work with high-level mathematical objects and continue to develop a program while it runs; (ii) a flexible grid structure that enables simple construction of discrete differential operators; and (iii) automatic differentiation that ensures that no analytical derivatives have to be programmed explicitly as long as the discrete flow equations and constitutive relationships are implemented as a sequence of algebraic operations. We have implemented a modular, efficient framework for implementing and comparing different physical models, discretizations, and solution strategies by combining imperative and object-oriented paradigms with functional programming. The toolbox also offers additional features such as upscaling and grid coarsening, consistent discretizations, multiscale solvers, flow diagnostics and interactive visualization.

A multiscale restriction-smoothed basis method for high contrast porous media represented on unstructured grids

A wide variety of multiscale methods have been proposed in the literature to reduce runtime and provide better scaling for the solution of Poisson-type equations modeling flow in porous media. We present a new multiscale restricted-smoothed basis (MsRSB) method that is designed to be applicable to both rectilinear grids and unstructured grids. Like many other multiscale methods, MsRSB relies on a coarse partition of the underlying fine grid and a set of local prolongation operators (multiscale basis functions) that map unknowns associated with the fine grid cells to unknowns associated with blocks in the coarse partition. These mappings are constructed by restricted smoothing: Starting from a constant, a localized iterative scheme is applied directly to the fine-scale discretization to compute prolongation operators that are consistent with the local properties of the differential operators.The resulting method has three main advantages: First of all, both the coarse and the fine grid can have general polyhedral geometry and unstructured topology. This means that partitions and good prolongation operators can easily be constructed for complex models involving high media contrasts and unstructured cell connections introduced by faults, pinch-outs, erosion, local grid refinement, etc. In particular, the coarse partition can be adapted to geological or flow-field properties represented on cells or faces to improve accuracy. Secondly, the method is accurate and robust when compared to existing multiscale methods and does not need expensive recomputation of local basis functions to account for transient behavior: Dynamic mobility changes are incorporated by continuing to iterate a few extra steps on existing basis functions. This way, the cost of updating the prolongation operators becomes proportional to the amount of change in fluid mobility and one reduces the need for expensive, tolerance-based updates. Finally, since the MsRSB method is formulated on top of a cell-centered, conservative, finite-volume method, it is applicable to any flow model in which one can isolate a pressure equation. Herein, we only discuss single and two-phase incompressible models. Compressible flow, e.g., as modeled by the black-oil equations, is discussed in a separate paper.

A multiscale two-point flux-approximation method

A large number of multiscale finite-volume methods have been developed over the past decade to compute conservative approximations to multiphase flow problems in heterogeneous porous media. In particular, several iterative and algebraic multiscale frameworks that seek to reduce the fine-scale residual towards machine precision have been presented. Common for all such methods is that they rely on a compatible primal-dual coarse partition, which makes it challenging to extend them to stratigraphic and unstructured grids. Herein, we propose a general idea for how one can formulate multiscale finite-volume methods using only a primal coarse partition. To this end, we use two key ingredients that are computed numerically: (i) elementary functions that correspond to flow solutions used in transmissibility upscaling, and (ii) partition-of-unity functions used to combine elementary functions into basis functions. We exemplify the idea by deriving a multiscale two-point flux-approximation (MsTPFA) method, which is robust with regards to strong heterogeneities in the permeability field and can easily handle general grids with unstructured fine- and coarse-scale connections. The method can easily be adapted to arbitrary levels of coarsening, and can be used both as a standalone solver and as a preconditioner. Several numerical experiments are presented to demonstrate that the MsTPFA method can be used to solve elliptic pressure problems on a wide variety of geological models in a robust and efficient manner.

The Multiscale Finite-Volume Method on Stratigraphic Grids

Finding a pressure solution for large and highly detailed reservoir models with fine-scale heterogeneities modeled on a meter scale is computationally demanding. One way of making such simulations less compute intensive is to employ multiscale methods that solve coarsened flow problems using a set of reusable basis functions to capture flow effects induced by local geological variations. One such method, the multiscale finite-volume (MsFV) method, is well studied for 2D Cartesian grids but has not been implemented for stratigraphic and unstructured grids with faults in 3D. We present an open-source implementation of the MsFV method in 3D along with a coarse partitioning algorithm that can handle stratigraphic grids with faults and wells. The resulting solver is an alternative to traditional upscaling methods, but can also be used for accelerating fine-scale simulations. To achieve better precision, the implementation can use the MsFV method as a preconditioner for Arnoldi iterations using GMRES, or as a preconditioner in combination with a standard inexpensive smoother. We conduct a series of numerical experiments in which approximate solutions computed by the new MsFV solver are compared with fine-scale solutions computed by a standard two-point scheme for grids with realistic permeabilities and geometries. On one hand, the results show that the MsFV method can produce accurate approximations for geological models with pinchouts, faults, and non-neighboring connections, but on the other hand, they also show that the method can fail quite spectacularly for highly anisotropic systems in a way that cannot efficiently be mitigated by iterative approaches The MsFV method is, in our opinion, not yet sufficiently robust to be applied as a black-box solver for models with industrystandard complexity. However, extending the method to realistic grids is an important step on the way towards fast and accurate multiscale solution of large-scale reservoir models. In particular, our open-source implementation provides an efficient framework suitable for further experimentation with partitioning algorithms and MsFV variants. Introduction Multiscale methods have been proposed as a way of bridging the gap in resolution between geological models (cell sizes: centimeters to decimeters in the vertical direction, meters to tens of meters in horizontal direction) and dynamic simulation models (cell sizes: meters to tens of meters in vertical direction, tens of meters to hundred of meters in horizontal direction). As an alternative to traditional upscaling techniques, multiscale methods (Efendiev and Hou 2009) can resolve fine-scale quantities with reduced computational complexity (Kippe et al. 2008) on highly detailed reservoir models by creating basis functions that relate localized flow problems on the fine scale (geological model) to a global flow problem on a much coarser scale (dynamic model). Herein, we consider the multiscale finite-volume (MsFV) method (Jenny et al. 2003) in which the basis functions in the MsFV method are computed using a dual-grid formulation with unitary pressure values at each vertex of the coarse blocks. The method has been extended to provide qualitatively correct approximations to wide variety of problems in subsurface flow, including density-driven flow (Lunati and Jenny 2008; Künze and Lunati 2012), compressible flow (Lunati and Jenny 2006; Hajibeygi and Jenny 2009; Zhou and Tchelepi 2008), three-phase and compositional flow (Lee et al. 2008; Hajibeygi and Jenny 2013), well modeling (Wolfsteiner et al. 2006; Jenny and Lunati 2009), and fractured systems (Hajibeygi et al. 2011b; Sandve et al. 2013). To overcome accuracy problems observed for strongly heterogeneous cases, and generally be able to systematically reduce fine-scale residuals, an iterative method was introduced by Hajibeygi et al. (2008) and further developed in (Hajibeygi and Jenny 2011). The key idea is to use the MsFV coarse-scale operator to reduce low-frequency errors and a combination of locally computed correction functions and inexpensive smoothers to reduce high-frequency errors. Likewise, the method has been formulated in an algebraic manner (Zhou and Tchelepi 2008; Lunati et al. 2011) and applied as a preconditioner (Lunati et al. 2011; Zhou and Tchelepi 2012; Wang et al. 2012). The close relationship between MsFV and domain-decomposition methods is discussed in more detail in (Nordbotten and Bjørstad 2008; Sandvin et al. 2011; Nordbotten et al. 2012; Sandvin et al. 2013). A particular advantage of the MsFV method, compared with multigrid and domain-decomposition methods, is that the multiscale

The Application of Flow Diagnostics for Reservoir Management

Flow-diagnostics, as referred to herein, are computational tools based on controlled numerical flow experiments that yield quantitative information regarding the flow behavior of a reservoir model in settings much simpler than would be encountered in the actual field. In contrast to output from traditional reservoir simulators, flow diagnostic measures can be obtained within seconds. The methodology can be used to evaluate, rank and/or compare realizations or strategies, and the computational speed makes it ideal for interactive visualization output. We also consider application of flow diagnostics as proxies in optimization of reservoir management workflows. In particular, based on finite volume discretizations for pressure, time-offlight (TOF) and stationary tracer, we efficiently compute general Lorenz coefficients (and variants) which are shown to correlate well with simulated recovery. For efficient optimization, we develop an adjoint code for gradient computations of the considered flow diagnostic measures. We present several numerical examples including optimization of rates, well-placements and drilling sequences for two and three phase synthetic and real field models. Overall, optimizing the diagnostic measures imply substantial improvement in simulation-based objectives.

The multiscale restriction smoothed basis method for fractured porous media (F-MsRSB)

A novel multiscale method for multiphase flow in heterogeneous fractured porous media is devised. The discrete fine-scale system is described using an embedded fracture modeling approach, in which the heterogeneous rock (matrix) and highly-conductive fractures are represented on independent grids. Given this fine-scale discrete system, the method first partitions the fine-scale volumetric grid representing the matrix and the lower-dimensional grids representing fractures into independent coarse grids. Then, basis functions for matrix and fractures are constructed by restricted smoothing, which gives a flexible and robust treatment of complex geometrical features and heterogeneous coefficients. From the basis functions one constructs a prolongation operator that maps between the coarse- and fine-scale systems. The resulting method allows for general coupling of matrix and fracture basis functions, giving efficient treatment of a large variety of fracture conductivities. In addition, basis functions can be adaptively updated using efficient global smoothing strategies to account for multiphase flow effects. The method is conservative and because it is described and implemented in algebraic form, it is straightforward to employ it to both rectilinear and unstructured grids. Through a series of challenging test cases for single and multiphase flow, in which synthetic and realistic fracture maps are combined with heterogeneous petrophysical matrix properties, we validate the method and conclude that it is an efficient and accurate approach for simulating flow in complex, large-scale, fractured media.

A simulation workflow for large-scale CO2 storage in the Norwegian North Sea

Large-scale CO2 injection problems have revived the interest in simple models, like percolation and vertically-averaged models, for simulating fluid flow in reservoirs and aquifers. A series of such models have been collected and implemented together with standard reservoir simulation capabilities in a high-level scripting language as part of the open-source MATLAB Reservoir Simulation Toolbox (MRST) to give a set of simulation methods of increasing computational complexity. Herein, we outline the methods and discuss how they can be combined to create a flexible tool-chain for investigating CO2 storage on a scale that would have significant impact on European CO2 emissions. In particular, we discuss geometrical methods for identifying structural traps, percolation-type methods for identifying potential spill paths, and vertical-equilibrium methods that can efficiently simulate structural, residual, and solubility trapping in a thousand-year perspective. The utility of the overall workflow is demonstrated using real-life depth and thickness maps of two geological formations from the recent CO2 Storage Atlas of the Norwegian North Sea.

A Multiscale Restriction-Smoothed Basis Method for Compressible Black-Oil Models

Simulation problems encountered in reservoir management are often computationally expensive because of the complex fluid physics for multiphase flow and the large number of grid cells required to honor geological heterogeneity. Multiscale methods have been proposed as a computationally inexpensive alternative to traditional fine-scale solvers for computing conservative approximations of the pressure and velocity fields on high-resolution geo-cellular models. Although a wide variety of such multiscale methods have been discussed in the literature, these methods have not yet seen widespread use in industry. One reason may be that no method has been presented so far that handles the combination of realistic flow physics and industrystandard grid formats in their full complexity. Herein, we present a multiscale method that handles both the most wide-spread type of flow physics (black-oil type models) and standard grid formats like corner-point, stair-stepped, PEBI, as well as general unstructured, polyhedral grids. Our approach is based on a finite-volume formulation in which the basis functions are constructed using restricted smoothing to effectively capture the local features of the permeability. The method can also easily be formulated for other types of flow models, provided one has a reliable (iterative) solution strategy that computes flow and transport in separate steps. The proposed method is implemented as open-source software and validated on a number of two and three-phase test cases with significant compressibility and gas dissolution. The test cases include both synthetic models and models of real fields with complex wells, faults, and inactive and degenerate cells. Through a prescribed tolerance, the solver can be set to either converge to a sequential or the fully implicit solution, in both cases with a significant speedup compared to a fine-scale multigrid solver. Altogether, this ensures that one can easily and systematically trade accuracy for efficiency, or vice versa.

Spill-point analysis and structural trapping capacity in saline aquifers using MRST-co2lab

Abstract Geological carbon storage represents a substantial challenge for the subsurface geosciences. Knowledge of the subsurface can be captured in a quantitative form using computational methods developed within petroleum production. However, to provide good estimates of the likely outcomes over thousands of years, traditional 3D simulation methods should be combined with other techniques developed specifically to study large-scale, long-term migration problems, e.g., in basin modeling. A number of such methods have been developed as a separate module in the open-source Matlab Reservoir Simulation Toolbox (MRST). In this paper, we present a set of tools provided by this module, consisting of geometrical and percolation type methods for computing structural traps and spill paths below a sealing caprock. Using concepts from water management, these tools can be applied on large-scale aquifer models to quickly estimate potential for structural trapping, determine spill paths from potential injection points, suggest optimal injection locations, etc. We demonstrate this by a series of examples applied on publicly available datasets. The corresponding source code is provided along with the examples.

A Multiscale Method Based on Restriction-Smoothed Basis Functions Suitable for General Grids in High Contrast Media

A wide variety of multiscale methods have been proposed in the literature to improve simulation runtimes and provide better scaling to large models. With a few notable exceptions, the methods proposed so far are mostly limited to structured grids. We present a new multiscale restricted-smoothed basis (MsRSB) method that is designed to be applicable to stratigraphic and fully unstructured grids. Like many other multiscale methods, it is based on a coarse partition of an underlying fine grid with a set of prolongation operators (also called multiscale basis functions) that map from unknowns associated with the fine grid cells to unknowns associated with the coarse grid blocks. These mappings are constructed by restricted smoothing: starting from a constant, a localized iterative scheme is applied directly to the fine-scale discretization to give prolongation operators that are consistent with the local properties of the differential operators. The resulting method has three main advantages: First, there are almost no requirements on the geometry and topology of the fine and the coarse grids. Coarse partitions and good prolongation operators can therefore easily be constructed on complex models involving high media contrasts and unstructured cell connections introduced by faults, pinch-outs, erosion, local grid refinement, etc. Moreover, the coarse grid can easily be adapted to features in the geo-cellular model or any precomputed flow field to improve accuracy. Secondly, the method is accurate and robust when compared to existing multiscale methods. In particular, the method does not need to recompute local basis functions to account for transient behavior: dynamic mobility changes are incorporated by continuing previous iterations with a few extra steps. This way, the cost of updating the prolongation operators will be proportional to the amount of change in fluid mobility and one avoids tolerance-based updates. Finally, since the MsRSB method is formulated on top of a cell-centered, conservative, finite-volume method, it is applicable to any flow model in which one can isolate a pressure equation; in the paper, we discuss incompressible two-phase flow and compressible, three-phase, black-oil type models. For our fine-scale discretization, we use the standard two-point flux-approximation scheme, but the method could equally well have been formulated using a multipoint discretization. Several numerical examples are presented to highlight features of the method. First, we compare the MsRSB method with the multiscale finite-volume (MsFV) method for single-phase flow problems with petrophysical parameters from the SPE 10 benchmark. Then we perform several validation studies using two-phase flow geometry and petrophysical properties from simulation models of two real fields (Gullfaks and Norne), as well as a compressible gas-injection case described by the black-oil equations. Introduction The general movement of fluids in a hydrocarbon reservoir is induced by global forces like gravity and pressure differentials. The micro-scale displacement, however, is determined by small-scale flow paths throughout highly heterogeneous porous rocks. Flow modeling therefore needs to take into account processes taking place on a wide range of spatial and temporal scales. The traditional approach to solving such problems is to use upscaling or homogenization techniques to develop effective parameters that represent a sub-scale behavior in an averaged sense so that flow can be simulated on a coarser scale. Such methods have proved to be very effective for problems with scale separation, e.g., as seen in material science. However, porous rocks seldom exhibit clear scale separation and upscaling techniques are therefore not as robust and accurate as one would wish. Effective properties are generally process dependent, and because one needs to assume a specific set of localization conditions to compute effective properties, upscaling techniques tend to only produce reliable results for a limited range of flow scenarios. In an attempt to overcome some of the limitations of upscaling methods, so-called multiscale discretization methods have been proposed over the past two decades to solve second-order elliptic equations with strongly heterogeneous coefficients (Efendiev and Hou 2009). This includes methods such as the generalized finite-element methods (Babuka et al. 1994), finite-element methods (Hou and Wu 1997), numerical-subgrid upscaling (Arbogast 2002), multiscale mixed finite-element methods (Chen and Hou