Statoil is implementing “Fast Model Update (FMU),” an integrated and automated workflow for reservoir modeling and characterization. FMU connects all steps and disciplines from seismic depth conversion to prediction and reservoir management taking into account relevant reservoir uncertainty. FMU delivers an ensemble of geologically-consistent and history-matched model realizations that together characterizes the reservoir uncertainty. FMU facilitates management of the reservoir (e.g. field development plan, well planning, and drainage strategies) in both early-phase as well as mature projects, properly taking geological uncertainty (including structure, facies, and rock/fluid properties) into account. The focus of the paper is to demonstrate how FMU is used in an algorithm for robust optimization of wells (e.g. well targets, infill wells, drilling priority, rate control). In the current paper we demonstrate some early results where the drilling sequence of wells is optimized under geological uncertainty, using an ensemble of models conditioned on all available data. The final product provides an optimized drilling schedule and drilling time. The paper includes examples from one synthetic and one real field application.
In this paper we construct coupled gas flow models in shale-matrix and hydraulic fractures within the framework of the reiterated homogenization procedure in conjunction with the treatment of fractures as (n-1) interfaces (
Flow equations are fundamental building blocks of any reservoir simulator and changing the formulation after the flow equations are interfaced with other building blocks such as solvers, facility networks, and well management may be virtually impossible. In a general purpose simulator it is also desirable to use the same formulation for black-oil, compositional, and thermal simulations. The challenge is therefore to find a formulation that works well all the time and to do so at a time when the simulator is at the beginning of its life-cycle and does not yet have all the features that are needed for thorough testing. This paper presents a general flow equation framework where the flow equations are expressed in terms of general material balance equations. The approach isolates the flow equations from the fluid model and makes the formulation independent of fluid type. Additionally, all derivatives are expressed in terms of generic gradients. Architecturally, this decouples the material balance equations from the choice of independent variables and isolates the part of the code that is formulation specific to a thin layer interfacing the flow equations with the property calculations. The end result is a modular flow equation framework which facilitates the use of any formulation with any fluid type.
The ensemble smoother with multiple data assimilation (ES-MDA) has been shown to outperform EnKF for both synthetic and field problems. Specifically, ES-MDA gives better data matches than EnKF, maintains the correct geology and appears to provide a better quantification of uncertainty than EnKF. However, if ES-MDA (or EnKF) is applied to update the rock-property fields where the underlying geological model is a facies model, then the boundaries between facies are not preserved. Here, we couple an ES-MDA update of the permeability field with an update of the distribution of facies for cases where both the distribution of geological facies and the distribution of permeability within each facies are unknown. To do this, the permeability field is represented as a Gaussian mixture model (GMM), where the permeability within each facies is represented by a different Gaussian probability distribution. ES-MDA is applied to update the permeability field in the normal way, but after each ES-MDA iteration, the facies value in each reservoir simulator gridblock is updated by calculating the probability of each possible facies with respect to the GMM. In addition, the updated permeability distributions for each facies are remapped to the original Gaussian distribution. To keep the facies distribution consistent with the underlying geological model, which in this work is based on multi-point statistics (MPS), every several ES-MDA iterations, the facies distribution is regenerated by using the facies probability map as soft data and by using certain permeability values as the hard data to avoid destroying the data match. For the example considered in this paper, the procedure is able to provide good data matches as well as posterior facies maps and permeability fields that reflect the main geological features of the true model.
Thermal-reactive-compositional flow simulation in porous media is essential to model unconventional thermal oil recovery processes for extra-heavy hydrocarbon resources, e.g., the In-situ Conversion Process (ICP) for oil-shale production. Computational costs can be very high for such a complex system, which may make simulation studies prohibitively time consuming for large field-scale applications on fine grids. On the other hand, significant errors are introduced with the use of coarse-scale models. In this paper, we developed an innovative multipoint multiscale modeling method to effectively capture the fine-scale reaction rates in coarse-scale ICP-simulation models. In our multiscale method, coupled thermal-reactive-compositional flow equations are solved only on the coarse-scale, with the kinetic parameters (frequency factors) calculated based on fine-scale reaction rates. We perform the temperature downscaling by solving the heat diffusion equation in local regions subject to temperature-gradient boundary conditions obtained from a multipoint evaluation on the coarse grid. A multiscale treatment is also developed for the heater well model. Coarse-scale heater well indices are calculated from fine-scale well models using downscaled temperatures. The newly developed multiscale method is applied to realistic cross-sectional ICP-pattern models with a vertical production well and multiple horizontal heater wells operated subject to a time-varying power. It is shown that the multiscale model delivers results that are in close agreement with the fine-scale reference results for all quantities of interest. Despite the fact that the multiscale method is implemented at the simulation-deck level, using the flexible scripting and monitor functionalities of our proprietary simulation package, significant computational improvements are achieved for all cases considered.
Geomechanical models of commercial reservoirs comprise hundreds of thousands to millions of grid cells. These lead to large systems of linear equations. The solution of these systems is time consuming and, in most cases, dominates the total runtime of geomechanical simulations (70 to 80% of the total runtime). In this paper, we will present a solution for accelerating the linear solver by combining a state of the art deflation algorithm with GPU technology to achieve a speed-up of 4 to 5 times for geomechanical simulations. In this paper, we will include details of both the deflation algorithm for geomechanics and our GPU implementation. We will conclude by presenting the numerical results obtained.
We present a method for porous media flow in the presence of complex fracture networks. The approach uses the Mimetic Finite Difference method (MFD) and takes advantage of MFD's ability to solve over a general set of polyhedral cells. This flexibility is used to mesh fracture intersections in two and three-dimensional settings without creating small cells at the intersection point. We also demonstrate how to use general polyhedra for embedding fracture boundaries in the reservoir domain. The target application is representing fracture networks inferred from microseismic analysis.
Algebraic multigrid (AMG) methods for directly solving coupled systems of partial differential equations (PDEs) have been extensively used in various types of numerical simulations in engineering. A necessary condition for its efficient applicability is the simulation process being driven by elliptic components. In reservoir simulation the pressure, described by Darcy's law, is known to drive the process and hence System-AMG should be applicable and outperform classical solvers.
In the context of adaptive and fully implicit reservoir simulations, the linearization of balance equations results in linear systems of equations that can be challenging. This makes it crucial to exploit all physical information from the full system to construct a robust AMG strategy and to extend it to more complex simulations than Black-Oil.
At the same time, the full set of information helps to get the best out of AMG. Just as with multigrid in general, System- AMG provides a framework for combining algorithmic modules rather than a fixed solution algorithm. The adaptation of the solution strategy to the concrete class of applications is the key to obtain the best performance.
Finally, System-AMG does not only allow for choosing more efficient algorithmic strategies, but also for exploiting parallelism in an optimized way, regardless of the simulation code being parallelized or not. Because the linear solver time is typically far dominating, even serial simulators can immediately and substantially benefit from MPI parallelism.
This paper considers the use of spatial-temporal (tensor) decompositions for the compact representations of saturation patterns in reservoir models. Reservoir flow patterns, in the sense of evolution of saturation patterns over time, can be considered to drive the economic performance of reservoirs as they drive the ultimate recovery. This makes the reservoir flow pattern a natural dissimilarity measure between models in the context of production optimization. We show that the application of multilinear algebra techniques allows the construction of low-complexity representations of the essential saturation patterns. The reservoir flow patterns are stored in large-scale multidimensional arrays, and tensor decompositions can be effectively used to describe the spatial-temporal behavior of the reservoir flow patterns. The dimensionality of the reservoir flow patterns can be substantially reduced, showing that a small number of spatial-temporal basis functions are required to characterize the dominant features. When applying the tensor decompositions to flow profiles of an ensemble of realizations they can be used for clustering models with similar dynamical properties, by allowing a fast calculation of a flow-relevant dissimilarity measure between realizations. For large ensembles of realizations they can lead to considerable computational advantages in (robust) optimization. This is illustrated by using a gradient-based technique to maximize Net Present Value (NPV) in water flooding for an ensemble of realizations. Besides in reducing ensemble sizes the resulting tools have potential use in constructing reduced order models.
Nonlinear convergence problems in numerical reservoir simulation can lead to unacceptably large computational time and are often the main impediment to performing simulation studies of large-scale problems. We analyze the nonlinearity of the discrete transport (mass conservation) equation for immiscible, incompressible, two-phase flow in porous media in the presence of viscous and buoyancy forces. Although simulation problems are multi-dimensional with large numbers of cells and variables, we find that the essence of the nonlinear behavior can be understood by studying the discretized (numerical) flux function for the interface between two cells. The numerical flux is expressed in terms of the saturations of the two cells. Discontinuities in the first derivative of the flux function (referred to as kinks) and inflection lines are identified as the cause of convergence difficulty. These critical features (kinks and inflections) change the curvature of the numerical flux function abruptly, and can lead to overshoots, oscillations, or divergence in Newton iterations.
Based on our understanding of the nonlinearity, a nonlinear solver is developed, referred to as the Numerical Trust Region (NTR) solver. The solver is able to guide the Newton iterations safely and efficiently through the different saturation ‘trust-regions’ delineated by the kinks and inflections. Specifically, overshoots and oscillations that often lead to convergence failure are avoided. Numerical examples demonstrate that our NTR solver has superior convergence performance compared with existing methods. In particular, convergence is achieved for a wide range of timestep sizes and Courant-Friedrichs-Lewy (CFL) numbers spanning several orders of magnitude. Our proposed numerical solution strategy that is based on the numerical flux extends the previous work by Jenny et al. (