A Multidomain Approach to Multiscale Compositional Reservoir Simulation

Seismic data are routinely used to infer in situ properties of earth materials on many scales, ranging from global studies to investigations of surficial geological formations. While inversion and imaging algorithms utilizing these data have improved steadily, there are remaining challenges that make detailed measurements of the properties of some geologic materials very difficult. For example, the determination of the concentration and orientation of fracture systems is prohibitively expensive to simulate on the fine grid and, thus, some type of coarse-grid simulations are needed. In this paper, we describe a new multiscale finite element algorithm for simulating seismic wave propagation in heterogeneous media. This method solves the wave equation on a coarse grid using multiscale basis functions and a global coupling mechanism to relate information between fine and coarse grids. Using a mixed formulation of the wave equation and staggered discontinuous basis functions, the proposed multiscale methods have the following properties. • The total wave energy is conserved. • Mass matrix is diagonal on a coarse grid and explicit energy-preserving time discretization does not require solving a linear system at each time step. • Multiscale basis functions can accurately capture the subgrid variations of the solution and the time stepping is performed on a coarse grid. We discuss various subgrid capturing mechanisms and present some preliminary numerical results.

Flow fields in a centrifugal fan for an indoor unit of an air-conditioner were calculated with finite element method-based large eddy simulation (LES) with the aim of predicting fan performance and aerodynamic noise in this study. The numerical simulation code employed throughout the LES was called FrontFlow/blue (FFB). We compared 10M grid [coarse grid] and 60M grid [fine grid] calculation results for investigation of influence of grid resolution. In the fine grid, the number of grid elements in blade-to-blade direction, and of region between the shroud and the bell mouth increased in particular. By calculating with the fine grid, calculated distributions of absolute velocities at blade exit reasonably agreed with experimental results. Because of this, maximum absolute velocity by fine grid near hub decreased as compared to those by coarse grid. Calculated sound pressure level by fine grid was therefore smaller than that by coarse grid, and the overestimation of sound pressure was suppressed by calculating with fine grid. This decrease of the absolute velocity was a first factor for the improvement of calculation accuracy. Moreover, number of captured streaks on the blade, hub, and shroud surfaces by fine grid increased as compared to those by coarse grid. As a result, size of streak by fine grid became smaller than that by coarse grid. Static pressure fluctuations by fine grid on the blade, hub, and shroud surfaces therefore reduced as compared to those by coarse grid. Aerodynamic noise was related to static pressure fluctuations according to Curle’s equation. This reduction of static pressure fluctuations was therefore a second factor for improvement of calculation accuracy.

Running a fine grid model with 107 - 109 of cells is possible using a supercomputer with 103 - 106 of CPUs but may not be always cost-effective. The most cost-effective way is to use a coarse grid model that is much smaller but with static/dynamic profiles very close to the fine grid model. This paper proposes a new layer optimization and upscaling method with the aim for creating a consistent coarse grid model. Unlike the industry's existing layer optimization and upscaling methods, the proposed method performs layer optimization and upscaling fully integrated with the Lorenz coefficient and curves (LCC). Coarse grid layers and their permeabilities are created by minimizing the difference between fine and coarse grid LCCs. The process consists of static and dynamic optimizations. The former is measured by LCC while the latter by pressure, GOR, and water-cut. A new LCC-based permeability upscaling method is developed to preserve the fine grid multiphase flow behaviors. A satisfactory coarse grid model is achieved when both static and dynamic criteria are met. The proposed method has been successfully applied to a giant carbonate oil field in the Caspian Sea that consists of a matrix dominated platform and a fracture/karst dominated rim. Due to the field's complex geology and high H2S content (15%), a dual porosity, dual permeability compositional model has been created to model compositional sour crude flow within and between the matrix and fracture/karst features. The reservoir drive mechanisms are fluid expansion, miscible gas injection and aquifer drive. The reservoir is undersaturated and has an abnormally high initial reservoir pressure. The fine-grid static model contains 104 million cells (370×225×625×2) and the optimized upscaled coarse-grid dynamic model has 8.3 million cells (370×225×50×2). The upscaled model can be run efficiently on the company's existing HPC infrastructure with a maximum of 64 CPUs. Excellent matches of the Lorenz coefficient maps for reservoir total/zones and Lorenz curves at all wells between the fine and coarse grid models have been achieved. Matches on the dynamic variables, e.g., pressure, gas breakthrough time, and GOR growth, in all producers are within the defined acceptable tolerances. The high quality of the static and dynamic matches between the coarse- and fine-grid models confirms that the reservoir properties of the coarse-grid model is very close to the fine-grid model and can be used a base model for history matching and uncertainty analysis.

Convergence analysis of two-grid fixed stress split iterative scheme for coupled flow and deformation in heterogeneous poroelastic media

This paper describes some model improvements developed to allow new possibilities in reservoirs representation with simulators using a finite difference scheme. It includes: possibility of flow between any grid block and from one to many blocks;a multiple time step calculation which allows to solve flow equations with several time step lengths for different group of grid blocks. These facilities provide a flexible and performant way of handling irregular and complex grid patterns with different nested levels of refinements. They are particularly useful to accurately simulate: –unconfined pilots or areas with smaller wells spacing;–faulted reservoirs;–vanishing or discordant geological layers;–reservoirs connected to a general aquifer. The new approach has been tested on a black oil simulator, implicit in pressure, implicit in saturation. The solutions used for some important aspects of model design are outlined. These are the Input/Output data organisation, the inner model grid pattern representation, the matrix solver, the multiple time steps selector and the associated equation solver. Two examples illustrate the use of the modified black oil simulator, and results are presented. The first example concerns the simulation of a repeated five spot injection. This test is made to check the local mesh refinements efficiency. The second example is a full field reservoir study for which using the modified model is suitable. The performance of a produced oil reservoir connected to an aquifer is simulated. In both cases, results obtained with a composite grid and the modified model are compared to those obtained with a coarse grid and then a fine grid used on a conventional model. The following points are outlined: –the use of local grid refinements give results very closed to those obtained with the fine gridded conventional model;–on the computers used, the running time of the modified model simulations is much lower than the time needed for the conventional model. Moreover, it is shown that the modified model results remain unchanged when the multiple time step selector is used. A different time step length is then automatically calculated for the oil reservoir and the aquifer, thus allowing a noticeable extra gain of processing running time.

: To illustrate Brandt's multigrid method this report considers the Helmholtz equation in cases where the frequency is close to an eigenfrequency. In a second part the eigenvalue problem is discussed. The basic idea for the treatment of the inhomogeneous problem is contained in Brandt's publications. Difficulties are encountered because the coarse grid eigenvalues differ from those of the finer grid. As a consequence, approximations to a solution in the coarse grid may differ widely from the corresponding approximation in a finer grid, if the frequency defined by the inhomogeneous problem is close to one of the eigenfrequencies. This may be detrimental to the aim of the coarse grid step; namely to remove from an approximation in a fine grid long waves which cause divergence of fine grid iterations. This is remedied by eliminating long wave contributions to the residual directly in the next finer grid. The success of this procedure still depends upon the choice of the subspace from which the long wave corrections in the finer grid are taken. Even this modified method fails if the frequency pertaining to the inhomogeneous term is too close to an eigenvalue, unless the subspace contains the pertinent eigenfunction.

Constraint energy minimizing generalized multiscale discontinuous Galerkin method

Based on two-grid discretizations, local and parallel finite element algorithms are proposed and analyzed for the time-dependent Oseen equations. Using conforming finite element pairs for the spatial discretization and backward Euler scheme for the temporal discretization, the basic idea of the fully discrete finite element algorithms is to approximate the generalized Oseen equations using a coarse grid on the entire domain, and then correct the resulted residual using a fine grid on overlapped subdomains by some local and parallel procedures at each time step. By the theoretical tool of local a priori estimate for the fully discrete finite element solution, error bounds of the approximate solutions from the algorithms are estimated. Numerical results are also given to demonstrate the efficiency of the algorithms.

This paper describes the use of vorticity to optimize the generation of coarse simulation grids from finely gridded geological models, thereby reducing (or removing) the need for two phase upscaling. Vorticity has been used in the simulation of flow in porous media since de Josselin de Jong1 first discussed it in 1960. However this is the first time, to our knowledge, that it has been used in the area of upscaling. Vorticity maps are generated from simulations of single-phase flow in the finely gridded geological model of interest. A coarse simulation grid is then derived which is more refined in areas of high vorticity and coarser in areas of low vorticity. The method is demonstrated for simple, layered and lens heterogeneities as well as Model 1 from the SPE Comparative Study in Upscaling2. Good matches between fine and coarse grid simulations of waterflooding are obtained, compared with using a conventionally generated coarse grid model.

We present a locally mass conservative scheme for the approximation of two-phase flow in a porous medium that allows us to obtain detailed fine scale solutions on relatively coarse meshes. The permeability is assumed to be resolvable on a fine numerical grid, but limits on computational power require that computations be performed on a coarse grid. We define a two-scale mixed finite element space and resulting method, and describe in detail the solution algorithm. It involves a coarse scale operator coupled to a subgrid scale operator localized in space to each coarse grid element. An influence function (numerical Greens function) technique allows us to solve these subgrid scale problems independently of the coarse grid approximation. The coarse grid problem is modified to take into account the subgrid scale solution and solved as a large linear system of equations posed over a coarse grid. Finally, the coarse scale solution is corrected on the subgrid scale, providing a fine grid representation of the solution. Numerical examples are presented, which show that near-well behavior and even extremely heterogeneous permeability barriers and streaks are upscaled well by the technique.

A method for performing nested grid calculations with a large-eddy simulation code is described. A common numerical method is used for all meshes, and the grid architecture consists of a single outer or coarse grid, and nested or fine grids, which overlap in some common region. Inter-grid communication matches the velocity, pressure and potential temperature fields in the overlap region. Resolved and sub-grid scale (SGS) turbulent fluxes and kinetic energy on the fine grid are averaged to the coarse grid using a conservation rule equivalent to Germano's identity used to develop dynamic SGS models. Simulations of a slightly convective, strong shear planetary boundary layer were carried out with varying surface-layer resolutions. Grid refinements in the (x, y, z) directions of up to (5, 5, 2) times were employed. Two-way interaction solutions on the coarse and fine meshes are successfully matched in the overlap region on an instantaneous basis, and the turbulent motions on the fine grid blend smoothly into the coarse grid across the grid interface. With surface-layer grid nesting, significant increases in resolved eddy fluxes and variances are found. The energy-scale content of the vertical velocity, and hence vertical turbulent fluxes, appear to be most influenced by increased grid resolution. Vertical velocity spectra show that the dominant scale shifts towards higher wavenumbers (smaller scales) and the magnitude of the peak energy is increased by more than a factor of 3 with finer resolution. Outside of the nested region the average heat and momentum fluxes and spectra are slightly influenced by the fine resolution in the surface layer. From these results we conclude that fine resolution is required to resolve the details of the turbulent motions in the surface layer. At the same time, however, increased resolution in the surface layer does not appreciably alter the ensemble statistics of the resolved and SGS motions outside of the nested region.

The miscible displacement of one incompressible fluid by another in a porous medium is governed by a system of two equations. One is elliptic form equation for the pressure and the other is parabolic form equation for the concentration of one of the fluids. Since only the velocity and not the pressure appears explicitly in the concentration equation, we use a mixed finite element method for the approximation of the pressure equation and mixed finite element method with characteristics for the concentration equation. To linearize the mixed-method equations, we use a two-grid algorithm based on the Newton iteration method for this full discrete scheme problems. First, we solve the original nonlinear equations on the coarse grid, then, we solve the linearized problem on the fine grid used Newton iteration once. It is shown that the coarse grid can be much coarser than the fine grid and achieve asymptotically optimal approximation as long as the mesh sizes satisfy h = H2 in this paper. Finally, numerical experiment indicates that two-grid algorithm is very effective.

New procedures for unstructured coarse-model generation are presented and applied. The underlying fine-grid model is considered to be unstructured, and the coarse-model cells are defined as groupings of fine-grid cells. The key flow quantity that must be computed for the coarse model is the upscaled transmissibility for each cell-to-cell connection. We introduce a global upscaling procedure for this computation. The method first requires several (minimum of three) global single-phase flow solutions. Appropriately defined linear combinations of these solutions are used to compute each upscaled transmissibility. This approach circumvents some of the limitations of existing (local and global) upscaling procedures. It also enables transmissibility to be quickly computed for a number of different coarse grids without performing any additional pressure solutions. Results are presented for an idealized two-phase flow problem. The fine grid contains nearly 200,000 cells, and coarse models of varying resolution are considered. Accurate results for total injector-producer flow rate are observed for all grid-resolution levels for the three different well configurations considered. Oil rate as a function of time is shown to improve in accuracy with increasing resolution, and is quite accurate for a model of about 10,000 cells.

On divergence- and gradient-preserving coarse-graining for finite volume primitive equation ocean models

A novel technique for upscaling of detailed geological reservoir descriptions is presented. The technique aims at reducing both numerical dispersion and homogenization error, generated due to incorporating a coarse computational grid and assigning effective permeability to coarse grid blocks respectively. In particular we consider implicit-pressure explicit-saturation (IMPES) scheme where homogenization error impacts the accuracy of the coarse grid solution of the pressure equation. To reduce the homogenization error, we employ the new vorticity-based gridding that generates a nonuniform coarse grid with high resolution at high vorticity zones. In addition, to control numerical dispersion, we use Dual Mesh Method (DMM), which uses different grids to solve pressure and saturation equations. The coarse grid generated from vorticity is used for computation of pressure and the reference fine grid is used for updating saturation explicitly. The most strong point of the method is that dual mesh method has been incorporated onto non-uniform grid structure. This combination removes the need to solve the full fine grid for two-phase flow modeling, resulting in a less computationally demanding and more accurate upscaling technique. To evaluate the method, we run two-phase simulation using different 2D test cases. Our results indicate that the non-uniform DMM is more accurate than the uniform DMM due to using vorticity-based grids. The speed-up achieved in the computation is significant depending on the complexity of model and degree of upscaling.