bluuebunny/crossref_metadata_2025 huggingface dataset:
Enhancement of the CPR Preconditioner for Efficient Solution of the Adjoint Equation
Authors: Han Choongyong , Wallis John , Sarma Pallav , Li Gary , Schrader Mark L. , Chen Wen
Published: February 2011

It is well known that the adjoint approach is the most efficient approach for gradient calculation, which can be used with gradient-based optimization techniques to solve various optimization problems such as the production optimization problem, history-matching problem, etc. The adjoint equation to be solved in the approach is a linear equation formed with the "transpose" of the Jacobian matrix from a fully implicit reservoir simulator. For a large and/or complex reservoir model, generalized preconditioners often prove impractical for solving the adjoint equation. Preconditioned specialized for reservoir simulation such as Constrained Pressure Residual (CPR) exploit properties of the Jacobian matrix in order to accelerate convergence, so they cannot be directly applied to the adjoint equation. To overcome this challenge, we have developed a new two-stage preconditioner for efficient solution of the adjoint equation (named CPRA).

The CPRA preconditioner has been coupled with an Algebraic Multi Grid (AMG) linear solver and implemented in Chevron's in-house reservoir simulator. The AMG solver is well known for its outstanding capability to solve the pressure equation of complex reservoir models; and solving the linear system with the "transpose" of pressure matrix is one of the two-stages of construction of the CPRA preconditioner.

Through test cases, we have confirmed that the CPRA/AMG solver with GMRES acceleration solves the adjoint equation very efficiently with reasonable number of linear solver iterations. Adjoint simulations to calculate the gradients with the CPRA/AMG solver take about the same time at most as corresponding CPR/AMG forward simulations. Accuracy of the solutions has also been confirmed by verifying the gradients against solutions with a direct solver. A production optimization case study for a real field using the CPRA/AMG solver has further validated its accuracy, efficiency, and the capability to perform long term optimization for large, complex reservoir models at low computational cost.
Numerical Studies on Two-way Coupled Fluid Flow and Geomechanics in Hydrate Deposits
Authors: Kim J. , Yang D. , Moridis G. J. , Rutqvist J.
Published: February 2011

Coupled flow and geomechanics play an important role in gas hydrate reservoirs because the stiffness of the rock skeleton, porosity and permeability are directly influenced by changes of the fluid (water and gas) and solid (hydrate and water) phase saturations, and the deformation of the reservoir. The fluid and solid phases coexist, yielding a high nonlinearity for flow and mechanics, so the coupled problem for hydrates is exceptionally complicated.

In previous study, the stability of hydrate-bearing sediments was assessed by one-way coupled analysis, where the change of fluid properties affects mechanics, but there is no feedback from mechanics to flow. In this paper, we develop and test rigorous two-way coupling between fluid flow and mechanics, where the solutions from mechanics are used to solve the flow problem. We employ the fixed-stress split for a convergent sequential implicit scheme.

We have found noticeable differences between one- and two-way couplings for several cases. The nature of the elliptic boundary problem of quasi-static mechanics results in instantaneous compaction or dilation over the domain by loading from reservoir fluid production. This yields the pressure rise-up or drop at early time (low pressure diffusion), which changes the effective stress instantaneously and geological instability may occur. Also, the change of pressure or temperature affects the solid saturation, rock stiffness, porosity, and permeability fields, changing the fluid flow regime. These behaviors cannot be captured by one-way coupling.

In conclusion, we developed, verified, and demonstrated the applicability of the tightly coupled sequential approach. It pro- vides a rigorous two-way coupled simulator which is ready to be applied to large scale problems to hydrate-bearing sediments.
Robust Constrained Optimization of Short and Long-Term NPV for Closed-Loop Reservoir Management
Authors: Chen Chaohui , Li Gaoming , Reynolds A. C.
Published: February 2011

In this paper, we develop an efficient algorithm for production optimization under linear and nonlinear constraints and an uncertain reservoir description. The linear and nonlinear constraints are incorporated into the objective function using the the augmented Lagrangian method and the bound constraints are enforced using a gradient-projection trust-region method. Robust long-term optimization maximizes the expected life-cycle net-present value (NPV) over a set of geological models, which represent the uncertainty of reservoir description. As the life-cycle optimal controls may be in conflict with the operator's objective of maximizing short-time production, the method is adapted to maximize the expectation of short-term NPV over the next one or two years subject to the constraint that the life-cycle NPV will not be substantially decreased. The technique is applied to synthetic reservoir problems to demonstrate its efficiency and robustness. Experiments show that the field can not always achieve the optimal NPV using the optimal well controls obtained based on a single, but uncertain reservoir model, whereas the application of robust optimization significantly reduces this risk. Experimental results also show robust sequential optimization on each short-term period is not able to achieve an expected life-cycle NPV as high as the one obtained with robust long-term optimization.
Accelerating Reservoir Simulators using GPU Technology
Authors: Appleyard John R. , Appleyard Jeremy D. , Wakefield Mark A. , Desitter Arnaud L.
Published: February 2011

The recent advances in graphics processing units (GPU) and associated development environments suitable for scientific modeling has generated significant interest in the high performance computing arena. In this paper we investigate strategies to incorporate this new technology into an existing commercial reservoir simulator. The use of the GPU for solving linear systems is demonstrated and the algorithmic considerations required in order to exploit the hardware are discussed. The paper describes a massively parallel incomplete factorization which is used as a preconditioner in conjunction with the GMRES algorithm. This factorization balances parallelism with accuracy resulting in a method that is significantly faster than the current serial solver when both are implemented on current hardware, despite in general requiring more iterations to converge. The performance of the implementation is shown to be dependent on problem size and indicates that when fully loaded the GPU is capable of producing a factor of 10 speed-up in the linear solver compared with the CPU based serial solver. The algorithm can also be applied on cluster systems, using domain decomposition, and although no numerical results are yet available, there are grounds for anticipating that performance will scale well for sufficiently large problems. The paper also discusses the benefits of migrating other simulator components of the simulator, such as the Jacobian matrix assembly, to the GPU. It is shown that this improves the overall performance of the simulator considerably.
A Comparison Between Wavelet and Renormalization Upscaling Methods and Iterative Upscaling-Downscaling Scheme
Authors: Babaei Masoud , King Peter R.
Published: February 2011

For all the advances in computational power of modern computers the problem of upscaling the properties of geological reservoir models is an ongoing challenge for reservoir simulation. Unfortunately, the upscaling of parameters like absolute permeability in single phase flow produces errors due to discarding the subgrid complexities. For tackling this problem, reconstruction techniques have been developed to downscale the local dynamic regions of the reservoir. Nonetheless, the upscaling part of these methods still has a significant effect on the final results since a poor upscaling can lead to errors even if coarse scale solution is downscaled or refined. In this work we focus on application of a range of conventional upscaling methods including pressure-solver, Renormalization, renormalization based modifications and wavelet transmissibility upscaling.

We simulate incompressible flooding scenarios in two and three dimensions. The results are subsequently compared to a modified local-global upscaling method coupled with modified nested-gridding downscaling. CPU runtime, error in breakthrough time and production values are the parameters considered in the comparative analysis.

We show improvements by iteration in providing better coarse scale solutions compared to conventional upscaling methods. Finally, we hope that this study can present a reliable ranking of renormalization and wavelet methods compared to iterative reconstruction method in terms of the trade-off between precision and speed.
Two-Stage Algebraic Multiscale Linear Solver for Highly Heterogeneous Reservoir Models
Authors: Zhou H. , Tchelepi H. A.
Published: February 2011

Previous research on multiscale methods for subsurface flow aims to obtain an efficient multiscale solution to the fine-scale problems. Such multiscale solution is usually a good approximation to the fine-scale problem. However, it has been reported that the multiscale solution deteriorates for high aspect ratios and channelized structures of permeability. Moreover, the multiscale solution does not converge to the fine-scale solution unless some special techniques are used.

In this paper, we propose an efficient two-stage algebraic multiscale (TAMS) method that converges to the fine-scale solution. The TAMS method consists of two stages, namely global and local stages. In the global stage, a multiscale solution is obtained purely algebraically from the fine-scale matrix. In the local stage, a local solution is constructed from a local preconditioner such as Block ILU(0) (BILU) or Additive Schwarz (AS) method. Spectral analysis shows that the multiscale solution step captures the low-frequency spectra in the original matrix very well and when combined with a local preconditioner that represents the high-frequency spectra, the full spectra can be well approximated. Thus the TAMS is guaranteed to converge to the fine-scale solution. Moreover, the spectra of the TAMS method tend to cluster together, which is a favorable property for Krylov subspace methods to converge fast. We also show that local mass conservation can be preserved if a multiscale solution step with the finite-volume restriction operator is applied before the iterative procedure exits. This allows for the TAMS method to be used to build efficient approximate solutions for multiphase flow problems that need to solve transport equations.

We test the numerical performance of the TAMS method using several challenging large-scale problems (fine-scale grid in the magnitude of one million) with complex heterogeneous structures and high aspect ratios. Different choices in the TAMS algorithm are employed, including the Galerkin or finite-volume type of restriction operator, BILU or AS preconditioner for the second stage, and the size of blocks for BILU and AS. The performance of TAMS is comparable or superior to the state-of-the-art algebraic multigrid (AMG) preconditioner when some optimal choice in the TAMS method is adopted. Moreover, the convergence of the TAMS method is insensitive to problem sizes, and the CPU time is almost linear to problem sizes. These indicate the TAMS method is efficient and robust for large-scale problems.
Reservoir Modelling with Feature Selection: A Kernel Learning Approach
Authors: Demyanov V. , Foresti L. , Christie M. , Kanevski M.
Published: February 2011

Geological realism is an increasingly important aspect of automated history matching and uncertainty quantification. While there are many algorithms that can automate the process of generating multiple history matched reservoir models, there is no guarantee that the changes made to the reservoir description to achieve a history match preserve geological characteristics.

This paper proposes a novel data-driven approach to modelling the spatial distribution of petrophysical properties in reservoirs. The proposed approach applies an advanced kernel-based approach – Multiple Kernel Learning (MKL) – to integrate multi-scale geological features. Measured reservoir properties and prior knowledge are fused in an intelligent way in order to improve the realism of reservoir models.

We show how we can use MKL to combine prior geological knowledge with seismic data. Geologically tailored features (e.g. large scale continuous fluvial sand bodies, small scale heterogeneity, spatial variability around facies boundaries, etc.) are identified from data using convolution filters to capture changes in gradients, orientations, sizes of geological bodies. Such spatial features are then used as separable input patterns for the MKL model to compute petrophysical reservoir properties through regression in high-dimensional kernel Hilbert space. MKL derives the internal model relations from data without enforcing predefined dependencies by optimising a weighted combination of kernels to fit to the conditioning data.

We used Multiple Kernel Learning to achieve realistic reproduction of geological properties of a complex fluvial reservoir system. It provides flexible control over the model properties by kernel weighting of spatial patterns, which reflect different geological features and data interpretations. MKL is also able to maintain the right balance between model complexity and goodness of fit, yielding robust history matched models.

Finally, Multiple Kernel Learning was integrated into a Bayesian uncertainty quantification framework for a synthetic reservoir to obtain uncertain production forecasts based on the ensemble of multiple history matched models.
A Dirichlet Neumann Representation Method for Simulating Flow in Reservoirs
Authors: Yang Yahan , Li Dachang , Parashkevov Rossen , Wu Xiao-Hui , Usadi Adam
Published: February 2011

Reservoir simulation models multiphase flow in reservoirs. Since geological models usually are heterogeneous and may contain millions of cells, upscaling is often a necessary step to reduce the model size to make them suitable for simulation. When grid coarsening is done to obtain fast simulation results, simulation grid cells can be quite large and may also contain localized geological features. Obviously, the upscaled model, though easier to simulate, is only an approximation of the original model. Errors introduced in the upscaling process generally may be small when flow field is relatively smooth but could be large otherwise. A critical technical challenge in upscaling is to be able to coarsen the reservoir model for fast computational performance while maintaining a high degree of accuracy for the results.

In this paper we propose a new Dirichlet Neumann representation method, DNR for short, for upscaling and simulating flow in reservoirs. With DNR, expressions are derived for flow rates as linear functions of multiple discrete pressure values along the boundary of each coarse block and pressure value at the center of the block. The number of pressure values at the boundary is flexible and may be chosen so that those pressure points provide an adequate representation of the possible pressure profiles and flow distributions during the simulation. Because of this, the error induced by artificial boundary conditions in traditional upscaling methods is greatly reduced. Test results show that DNR produces substantially more accurate results when compared to existing approaches using similar number of unknown variables, especially for reservoirs which contain high permeability streaks or are highly channelized. This method makes it possible to coarsen reservoir models while preserving a high degree of reliability for simulation results. DNR also allows straightforward reconstruction of fine-scale flow solutions and so can be used as part of a multiscale computational procedure.
Accurate Cell-Centered Discretizations for Modeling Multiphase Flow in Porous Media on General Hexahedral and Simplicial Grids
Authors: Wheeler Mary F. , Xue Guangri , Yotov Ivan
Published: February 2011

We introduce an accurate cell-centered method for modeling Darcy flow on general quadrilateral, hexahedral, and simplicial grids. We refer to these discretizations as the multipoint flux mixed finite element (MFMFE) method. The MFMFE method is locally conservative with continuous fluxes and can be viewed within a variational framework as a mixed finite element method with special approximating spaces and quadrature rules. We study two versions of the method: with a symmetric quadrature rule on smooth grids and a non-symmetric quadrature rule on rough grids. The framework allows for handling hexahedral grids with non-planar faces defined via trilinear mappings from the reference cube. Moreover, the MFMFE method allows for local elimination of the velocity, which leads to a cell-centered pressure system. Theoretical and numerical results demonstrate first order convergence on rough grids. Second order superconvergence is observed on smooth grids. We also discuss a new splitting scheme for modeling multiphase flows that can treat higher order transport discretizations for saturations. We apply the MFMFE method to obtain physically consistent approximations to the velocity and a reference pressure on quadrilateral or hexahedral grids, and a discontinuous Galerkin method for saturations. For higher order saturations, we propose an efficient postprocessing technique that gives accurate velocities in the interior of the gridblocks. Computational results are provided for flow in highly heterogenous reservoirs, including different capillary pressures arising from different rock types.
Simulating Two-phase Flow in Porous Media with Anisotropic Relative Permeabilities
Authors: Keilegavlen E. , Nordbotten J. M. , Stephansen Annette
Published: February 2011

In multi-phase flow simulation, the prevailing approach to discretizing flux terms treats the elliptic and the hyperbolic terms in the equations separately. With this concept, the flux is calculated analogously as for single phase flow, and later multiplied by the upstream mobility. This approach is valid when the mobility is a scalar quantity, which is the case for most traditional models. However, tensorial relative permeability (and thus tensorial mobilities) in general arise on all scales, as is seen in both laboratory and field experiments. Furthermore, upscaling methods almost invariably lead to anisotropic relative permeabilities. The saturation dependency of the fluid permeability tensor means that the upstream direction is no longer uniquely defined, which challenges common numerical schemes. In this work, we give examples of how the relative permeability has preferential flow directions on different scales. We then study how to incorporate tensorial relative permeability fields into control volume methods. In particular, we address two immediate challenges: Firstly, the non-determinacy of the upstream direction invalidates the use of upstream weighting for the saturation equation. We discuss the applicability of Godunov methods to handle flow cases where interfaces have either no or two upstream directions. Secondly, there is a marked increase in computational complexity associated with the pressure equation. We discuss the possibility of mitigating this computational complexity through mass lumping methods. The validity and computational efficiency of our approaches is discussed on theoretical terms and with the support of numerical implementations.
Automatic Differentiation Framework for Compositional Simulation on Unstructured Grids with Multi-Point Discretization Schemes
Authors: Zhou Y. , Tchelepi H. A. , Mallison B. T.
Published: February 2011

We present a flexible general-purpose reservoir simulation framework based on Automatic Differentiation (AD). The new AD-based simulator supports unstructured grids, employs a generalized Multi-Point Flux Approximation (MPFA) for spatial discretization, and uses a multi-level Adaptive Implicit Method (AIM) for time discretization. Given the discrete form of the governing nonlinear residual equations and a declaration of the independent variables, the AD library employs advanced expression templates with block data-structures to automatically generate compact computer code for the Jacobian matrix. Test results indicate that the construction of the Jacobian matrix with MPFA is efficient, and the overhead associated with treating a two-point flux approximation (TPFA) as a special case of MPFA is negligible. Our AIM implementation is designed to facilitate a systematic application of the method to new fluid models and variable formulations. The AD simulator allows for any combination of TPFA (Two-Point Flux Approximation), MPFA, FIM (Fully Implicit Method), and AIM. The generic and modular design is amenable to extension, both in terms of modeling additional flow processes and implementing new numerical methods. The AD-based modeling capability is demonstrated for highly nonlinear compositional problems using challenging large-scale reservoir models that include full-tensor permeability fields and non-orthogonal grids. The behaviors of TPFA and several MPFA schemes are analyzed for both FIM and AIM simulations. The implications of using MPFA and AIM on both the nonlinear and linear solvers are discussed and analyzed.
Domain Decomposition Methods Applied to Coupled Flow-Geomechanics Reservoir Simulation
Authors: Florez Horacio , Wheeler Mary , Rodriguez Adolfo A. , Monteagudo Jorge E.
Published: February 2011

In this work, well established Domain Decomposition techniques have been studied in order to carry out efficient simulations of coupled flow/geomechanics problems by taking full advantage of current parallel computer architectures. Different solution schemes can be defined depending upon transmission conditions among sub-domain interfaces. Three different schemes, i.e. Dirichlet-Neumann, Neumann-Neumann and Mortar-FEM, are tested and the advantages and disadvantages of each of them identified. This work will focus in the coupling of different meshes and/or physics on different domains by means of the Mortar-FEM plus the above DD-Schemes. Several examples of coupling of elasticity and poroelasticity in the context of reservoir compaction and subsidence are presented. In order to facilitate the implementation of complex workflows, we have implemented an advanced Python wrapper interface that allows programming capabilities. We have applied this platform to a variety of problems ranging from near-wellbore applications to field level subsidence calculations.
A Practical Approach for Assisted History Matching Using Grid Coarsening and Streamline-Based Inversion: Experiences in a Giant Carbonate Reservoir
Authors: Taware Satyajit , Friedel Torsten , Datta-Gupta Akhil
Published: February 2011

Flow simulation of multimillion grid cell models with hundreds of wells and decades long production history can be extremely time-consuming. This often limits the applicability of assisted history matching techniques. A pragmatic solution to this problem is grid coarsening which is now embedded in many commercial reservoir simulators.

Instead of up-scaling geological models in external packages, grid cells are automatically amalgamated within the simulator while preserving flux distribution and reducing the total number of active cells. The resulting speedup can be significant, often only at small loss of accuracy. Both characteristics are essential elements of any inversion technique in a multimillion grid cell environment. For water-flood history matching, we have utilized a commercial finite-volume simulator and the streamline-based generalized travel time inversion whereby water-cut behavior is matched by adjusting inter-well permeabilities.

To apply the assisted history matching technique on high resolution models, a flux reconstruction method is devised which makes full use of the benefits of automatic grid coarsening. It approximates fluxes in the original geological grid by recalculating transmissibility at the fine scale and redistributing coarse scale fluxes accordingly. The inversion is thus conducted on the fine scale grid, while the forward simulation model uses the coarse grid. While solving for an inexpensive coarse model, the calibration takes into account fine grid resolution. Upscaling is not required and the degree of coarsening can be adjusted during the matching process.

The proposed method was successfully tested on a supergiant carbonate oilfield with about hundred wells and large-scale water injection. The history match improved dramatically at relatively low numerical cost, which also allowed for investigating multiple sensitivities. The results were verified against non-coarsened model. The significant increase in efficiency makes this a potential method of choice for cases, where previously assisted history matching techniques could not be deployed due to excessive run times.
An Expanded Well Model for Accurate Simulation of Well-Reservoir Interactions
Authors: Karimi-Fard M. , Durlofsky L. J.
Published: February 2011

We present a new framework for modeling wells in reservoir simulators. The approach is based on an expanded well model in which the well region is expanded geometrically to include portions of the reservoir. The method proceeds by first defining an underlying fine-scale model, in which the well and any key near-well features such as hydraulic fractures are fully resolved using an unstructured grid. The geometry, grid structure and transmissibilities in the coarse-model well region, and the linkage between the well and reservoir regions, are determined from the solution of a single-phase steady-state pressure equation on the underlying fine grid. The method is compatible with any reservoir simulator that allows general (unstructured) cell-to-cell connections. The expanded well modeling procedure is applied to challenging three-dimensional problems involving production in a low-permeability reservoir, tight-gas production via a hydraulically-fractured well, and production in a gas-condensate reservoir. In the first two cases, where it is possible to simulate the fine-grid unstructured model, results using the expanded well model closely match the reference solutions, while standard approaches lead to significant error. In the gas-condensate example, which involves a nine-component compositional model, the reference solution is not computed, but the solution using the expanded well model is shown to be physically reasonable while standard coarse-grid solutions show large variation under grid refinement.
Joint History Matching of Well Data and Surface Subsidence Observations Using the Ensemble Kalman Filter: A Field Study
Authors: Wilschut F. , Peters E. , Visser K. , Fokker P. A. , van Hooff P. M.
Published: February 2011

The number of reported applications of the Ensemble Kalman Filter (EnKF) for history matching reservoir models is increasing steadily for various reasons. Here, we report on exploiting the capability of EnKF to handle observations from different sources simultaneously. While traditionally only well data are matched, we use surface subsidence observations together with well data. Surface subsidence results from compaction of the reservoir rock through the mechanical response of the subsurface. Compaction is caused by decreasing pore pressures during reservoir depletion. Therefore, the subsidence data contains information about dynamic pressure distributions in the reservoir.

The joint history matching of well data and surface subsidence observations was applied to the Roswinkel gas field. This field has been operated for 25 years, during which nine leveling campaigns generated a valuable data set of subsidence data. An important feature of the reservoir was the uncertainty about its compartmentalization, due to a large number of possibly sealing faults in the anticlinal structure. Therefore, instead of uncertain rock properties, fault transmissibilities were estimated in this study.

In a previous study on Roswinkel, a compaction field was estimated by inverting subsidence measurements. The results indicated several sealing faults in the reservoir, dividing the field into different compartments with independent pressure histories. The post-inversion history match of production data, however, was unsatisfactory. We have now been able to show that estimating the driving parameters, in casu the fault transmissibilities in the reservoir, can be achieved in a consistent way with both production and surface subsidence data for a synthetic case. Furthermore, for the actual Roswinkel case, the joint history match clearly reveals the discreapancy between the data from both sources and the current reservoir simulation model.

The joint history match of land surface movement data together with well production data is a success for EnKF as a flexible method for history matching. But more importantly, it demonstrates the potential of using complementary sources of information for improved reservoir characterization, and the possibility of estimating the driving parameters.
Application of Dynamic Gridding Techniques to IOR/EOR-Processes
Authors: van Batenburg D. W. , Bosch M. , Boerrigter P. M. , de Zwart A. H. , Vink J. C.
Published: February 2011

IOR/EOR processes have moving fronts and require an accurate description at a scale much smaller than typical grid block sizes used in black-oil simulations. Furthermore, they may require extra fluid components and phases, or thermal properties, or sometimes chemical reactions between components. However, simulation of these processes on a high resolution grid throughout the reservoir is still not practical with today's computing power. One solution is to start on a coarse grid that is dynamically adjusted to provide sufficient resolution where needed.

The challenge in such an approach is in identifying where and when to adjust the grid. This paper presents a semi-implicit scheme, that allows for the evaluation of grid adaption criteria within a time step, and guarantees that grid refinements occur where the fronts have propagated at the end of each time step. No assumptions on front movement need to be made in advance. The scheme is robust and accurate, while still providing strict control on the number of blocks by minimizing the extent of refined zone(s). This is of particular importance for nested refinements as the number of blocks rapidly grows with each extra refinement level. The combined use of changes in properties in space and in time allows for the definition of criteria that always reliably position the refinements along the fronts.

Application of this dynamic gridding approach in simulation models with water-flooding, chemical-flooding, gas injection, and in-situ combustion demonstrates that efficient nested dynamic gridding can be implemented in a general purpose simulator to provide sufficient physical and spatial detail in meaningful field or pattern models for IOR/EOR. The key to fast and accurate dynamic gridding is the implicit evaluation of proper criteria for refining and coarsening.
Multiple-Objective Optimization Applied to Well Path Design under Geological Uncertainty
Authors: Schulze-Riegert R. , Bagheri M. , Krosche M. , Kück N. , Ma D.
Published: February 2011

Well placement and design under geological uncertainties defines a major risk for field development processes. Including alternative geological realizations in a manually optimized development plan has been cumbersome and time consuming. In practical cases, a single reference model was often used as a basis for the flow simulation model which included a major risk on underestimating geological uncertainties. Optimization approaches and related optimization methods applied in this field previously included an objective function definition which aggregates multiple objectives, e.g. cumulative well rates or economical indicators. These optimization methods do not allow distinguishing easily between correlated and uncorrelated objectives which defines a major limitation for efficiency.

This paper presents a novel approach individually linking objectives and design parameters to multiple optimization processes which operate in parallel on a single simulation case. This approach increases the efficiency of stochastic optimization processes for complex application scenarios with multiple objectives. The extension to an uncertainty workflow including multiple geological models in the optimization processes is discussed.

Well placement optimization of horizontal wells in independent target formations of the gas condensate field is used to verify the optimization methodology. Alternative geological models are included in the optimization workflow representing geological uncertainties. The specific objective of the case study is to find optimum well trajectories by maximizing cumulative gas production over the life cycle of the field. Performance indicators for all well candidates (design parameters) are related to individual objectives (regional production). Robustness criteria for well candidates are implicitly fulfilled on the basis of risk weighted geological realizations.
Automatic Decomposition for Parallel Reservoir Simulation
Authors: Tarman Mehmet , Wang Kefei , Killough John E. , Sepehrnoori Kamy
Published: February 2011

A large amount of research in domain decomposition methods for parallel computing has been carried out in the last twenty years. This paper presents an automatic domain decomposition framework for parallel reservoir simulation. It emphasizes practical issues and implementation rather than theoretical complexity.

At least three main factors affect the parallel simulator performance when decomposing the domain for the given total number of sub-grids: 1) the active cell distribution in the reservoir, 2) reservoir geometry, and 3) fault locations. This discussion primarily focuses on the active cell distribution to improve the performance of the parallel reservoir simulator by load balancing.

Four domain decomposition algorithms are implemented before determining the best decomposition for the given total number of sub-domains. These four algorithms applied to all factors of the total number of sub-domains. The ideal decomposition for active cells and the offset size between the ideal number and real number from each scenario are then calculated. The smallest offset size determines the best decomposition.

Several examples, such as a synthetic model, a GOM (Gulf of Mexico) model with over 550,000+ grid blocks and a synthetic African reservoir compositional model, illustrate the improved parallel simulator performance resulting from the automatic decomposition framework. For both cases, the reservoir models were decomposed by different algorithms, and the CPU time for each algorithm was summarized. The automatic decomposition framework provides an effective way to improve the parallel simulator performance.
Compositional Nonlinear Solver Based on Trust Regions of the Flux Function Along Key Tie-lines
Authors: Voskov D. V. , Tchelepi H. A.
Published: February 2011

A new nonlinear solver that combines tie-line and fractional-flow information is developed for general-purpose compositional simulation. The solver, which takes advantage of the hyperbolic nature of the transport equations, is demonstrated for a compositional formulation that employs molar (or mass) variables. Since compositional recovery processes evolve along a few ‘key’ tie-lines, the flux functions (fractional flow curves) associated with these key tie-lines play a dominant role in the evolution of the solution. A two-stage nonlinear solver is proposed, and its applicability and efficiency for general-purpose compositional simulation is demonstrated. In the first stage, given the pressure and overall composition of the current nonlinear iteration, the key tie-lines are computed using a negative-flash procedure, and the nonlinear unknowns are updated accordingly. In the second stage, the flux functions associated with the computed tie-lines are segmented into trust regions (i.e., appearance and disappearance of a phase, and the inflection point of the S-shaped two-phase region), which are then used to guide the consistent evolution of the compositions. These trust regions delineate convex parts of the flux function, where convergence of the Newton iterations is guaranteed. Several challenging, multi-contact, miscible compositional problems, including the SPE 10 permeability field with the 6-component SPE 5 fluid, are used to test the robustness and efficiency of this new nonlinear solver. The convergence rate of the new nonlinear solver is always better than our standard safeguarded Newton method, which employs heuristics on maximum changes in the variables. We demonstrate that for aggressive time stepping (i.e., maximum throughput - so called CFL numbers - is hundreds to thousands of times larger than the IMPES limit), the new nonlinear solver converges within a few Newton iterations for the wide range of problems we studied.
Pressure Preconditioning Using Proper Orthogonal Decomposition
Authors: Astrid P. , Papaioannou G. V. , Vink J. C. , Jansen J. D.
Published: February 2011

We developed and implemented a new first-stage preconditioning method for large-scale reservoir simulation as an alternative to the popular Algebraic Multi Grid (AMG) method. We used Proper Orthogonal Decomposition (POD) to derive a reduced-order model for the linearized pressure equation in a proprietary reservoir simulator. A small set of pre-computed pressure solutions are used to transform the equation into a lower-order system that can be solved economically yet provides a relatively accurate estimation of the solution of the original high-order system. We present results for a two-phase (oil-water) reservoir model under water flooding conditions. Use of the POD-based preconditioner lead to a significantly faster convergence of the overall simulation compared to AMG, and reduced the linear solver times between approximately 70 % and 80%, i.e. we obtained accelerations with factors between three and five. The highest speed up was achieved for a case in which the flow rates in the eight injection wells were changed frequently during the entire production period, using ten POD basis functions obtained from short training simulations in which the injectors were operated only once each. These first results are very encouraging, especially because there is room to optimize the trade-off between the number of basis functions used in the POD method and the required number of linear iterations. However, further research is required to assess the applicability of the POD-based preconditioner to more complex cases including e.g. strongly compressible flow or compositional effects. The overhead required to pre-compute the POD solutions implies that the new method will be particularly attractive when many solutions of near-similar simulation models are required such as in computer-assisted history matching or flooding optimization.
An Ensemble Smoother for assisted History Matching
Authors: Skjervheim J.-A. , Evensen G. , Hove J. , Vabø J. G.
Published: February 2011

This paper compares two ensemble-based data-assimilation methods when solving the history-matching problem in reservoir-simulation models. The methods are the Ensemble Kalman Filter (EnKF) and the Ensemble Smoother (ES). Several publications have discussed the use of EnKF in petroleum applications while ES is now used for the first time for history matching. ES differs from EnKF by computing a global update in the space-time domain, rather than using recursive updates in time as in EnKF. Thus, the sequential updating of the realizations with associated restarts is avoided.

EnKF and ES provide identical solutions for state estimation with linear dynamical models. However, for nonlinear dynamical models, and in particular models with chaotic dynamics, EnKF is superior to ES, due to the fact that the recursive updates keep the model on track and close to the true solution. Thus, ES is not much used and EnKF has been the method of choice in most data assimilation studies where ensemble methods are used.

On the other hand, reservoir simulation models are rather diffusive systems when compared to the chaotic dynamical models that were previously used to test ES. If we can assume that the model solution is stable with respect to small perturbations in the initial conditions and the history-matching parameters, then ES should give similar results to EnKF, and ES will be a more efficient and much simpler method to implement and apply.

The technical advantages of using ES compared to EnKF are severe, especially when the methods are applied with complex real reservoir models. ES provides a significant reduction in simulation time. Furthermore, a more flexible parameterization is possible, which makes it easier to handle structural and geological model parameters in the history-matching process.

In this paper we compare EnKF and ES and show that ES indeed provide for an efficient ensemble-based method for history matching.
On the Coupling of Finite Volume and Discontinuous Galerkin for Reservoir Simulation Problems
Authors: Riviere Beatrice , Chidyagwai Prince , Mishev Ilya
Published: February 2011

Finite volume (FV) methods are still the most popular methods in practical reservoir simulation. They have good accuracy when used properly and are computationally very cheap. FV can be used on unstructured grids, but grid generation becomes a very difficult task. More advanced methods, discontinuous Galerkin (DG) for example, work on more general meshes, but are computationally more expensive and are still not accepted by the practitioners. We propose coupling of FV and DG methods that can improve the accuracy of the FV with reasonable increase of the computational cost.

We developed the algorithms for coupling DG and FV for model diffusion and convection diffusion problems on Voronoi /PEBI grids. We demonstrate with examples how DG can be used to alleviate the problems with grids in 2-D, gridding around pinch-outs and patch local refinement for areas with tensor coefficients.
A Loosely Coupled Hierarchical Fracture Model for the Iterative Multiscale Finite Volume Method
Authors: Hajibeygi H. , Karvounis D. , Jenny P.
Published: February 2011

The hierarchical fracture modeling approach is extended to a loosely coupled strategy in order to enhance its efficiency for highly conductive fractures and to easily integrate it into the iterative Multiscale Finite Volume (i-MSFV) framework. The fracture and matrix equations are decoupled and a sequential implicit approach is used to solve them. The i-MSFV method is then introduced for the efficient solution of the matrix equation which also involves fracture terms. Very important it that only one degree of freedom per connected fracture network appears in the coarse multiscale operator. This results in very efficient simulations for realistic problems, and also ensures total mass balance for each fracture network during the iterations. Local fracture functions are introduced in order to capture the fractures at the coarse scale accurately; similar to the complex well functions developed by Jenny and Lunati (2009). Numerical results show that the i-MSFV is an efficient conservative multiscale solver for simulations of flow in fractured porous media.
A New Method to Model Relative Permeability in Compositional Simulators to Avoid Discontinuous Changes Caused by Phase Identification Problems
Authors: Yuan Chengwu , Pope Gary A.
Published: February 2011

Simple methods such as using density during compositional simulations often fail to identify the phases correctly and this can cause discontinuities in the computed relative permeability values. The results are then physically incorrect. Furthermore, numerical simulators often slow down or even stop due to discontinuities. There are many important applications where the phase behavior can be single phase, gas-liquid, liquid-liquid, gas-liquid-liquid or gas-liquid-solid at different times in different grid blocks. Assigning physically correct phase identities during a compositional simulation turns out to be a very difficult problem that has resisted a general solution for decades. We know the intensive thermodynamic properties such as molar Gibbs free energy must be continuous assuming local equilibrium, but this condition is difficult to impose in numerical simulators because of the discrete nature of the calculations. An alternative approach is to develop a relative permeability model that is continuous and independent of the phase numbers assigned by the flash calculation. Relative permeability is a function of saturation, but also composition since composition affects the phase distribution in the pores i.e. the wettability. The equilibrium distribution of fluids in pores corresponds to the minimum in the Gibbs free energy for the entire fluid-rock system including interfaces. However, in general this relationship is very difficult to model from first principles. What we can easily do is calculate the molar Gibbs free energy (G) of each phase at reference compositions where the relative permeabilities are known or assumed to be known and then interpolate between these values using the G calculated during each time step of the simulation. Relative permeability values calculated this way are unconditionally continuous for all possible phase behavior changes including even critical points. We tested the new relative permeability model on a variety of extremely difficult simulation problems with up to four phases and it has not failed yet. We illustrate several of these applications.
CO2 Injection in Sub-salt Water Layers at 7000m Depth
Authors: Eftekhari Ali Akbar , Bruining J. , Wahanik Helmut , Marchesin Dan
Published: February 2011

Huge oil deposits (> 8 billion barrels) were recently discovered off the coast of Brazil below 2000 m of seawater and 1000–2000 m of reservoir, below 2000 m of salt layers. The high pressures (~ 700 bars) make these reservoirs excellent candidates for high pressure miscible gas injection, e.g., CO2 injection. For the geothermal gradients encountered in this area the temperatures are around 473 K. At high pressures it is possible that the mass density of carbon dioxide is higher than of oil, which will allow gravity stabilized injection. We consider the CO2-brine system, firstly for its storage as a greenhouse gas. Secondly, it is an initial step to describe displacement in the CO2-oil-brine system.

In this paper we study the phase behavior of carbon dioxide in a brine bearing layer at a temperature of 473 K and pressures between 600 and 800 bars. The quantification of the displacement process requires the thermodynamic equilibrium relations in the oil-water-salt-CO2 system at the prevailing pressures. Here we limit our interest to the thermodynamic behavior of the carbon dioxide-brine system. We show that it is possible to use the PRSV equation of state to estimate the composition of the vapor and aqueous phases. The PRSV is a modification of the Peng-Robinson equation of state. We use an activity-coefficient based mixing rule for the thermodynamic calculations. A volume shift procedure is applied to improve liquid density.