This paper describes a new iterative technique for the solution of linear systems of equations arising in petroleum reservoir simulation. The procedure exploits a powerful, vectorizable preconditioner to accelerate generalized conjugate residual methods. Speed and robustness of the new method compare favorably with previously published schemes. These attributes will be demonstrated in actual applications.

We describe a multistage parallel linear solver framework developed as part of the Intersect (IX) next-generation reservoir simulation project. The object-oriented framework allows for wide flexibility in the number of stages, methods and preconditioners. Here, we describe the specific components of a two-stage CPR[1] (Constraint Pressure Residual) scheme designed for large-scale parallel, structured and unstructured linear systems. We developed a highly efficient in-house Parallel Algebraic Multigrid (PAMG) solver as the first stage preconditioner. For the second stage, we use a parallel ILU-type scheme. This new and powerful combination of CPR and PAMG was the result of detailed analysis of the linear system of equations associated with reservoir simulation. Using several difficult reservoir simulation problems, we demonstrate the robustness and excellent parallel scalability of the IX linear solver. For the field case studies, the IX linear solver with CPR and PAMG is at least five times faster than an established and widely used industrial linear solver. The performance advantage of the IX linear solver over traditional reservoir simulation linear solvers increases with both problem size and the number of processors. Introduction Different types of grid may be used for reservoir flow simulation to model geometrically complex, highly detailed models and/or deviated or multi-lateral wells[2]. Grid types are often labeled based on their structure. Examples of simulation grids include:structured Cartesian,structured stratigraphic,multi-block stratigraphic,PEBI (Perpendicular Bisector), andgenerally unstructured. Hybrid grids that combine various types can also be used. It is now widely recognized that complete flexibility in representing complex and highly detailed simulation models can be achieved using generally unstructured grids[3].In recent years, significant efforts have focused on building multi-purpose reservoir flow simulators that can deal with geometrically complex and highly detailed structured and unstructured reservoir models[4,5,6]. These relatively large-scale efforts are being pursued because, for nearly three decades, the reservoir simulation community has focused on building robust and efficient reservoir simulators for structured grid problems. Today, the ability to routinely simulate a wide spectrum of practical black-oil problems on (effectively) structured models with O(105) gridblocks is widespread. However, the performance of traditional reservoir simulators typically deteriorates significantly with problem size and the number of processors. This is because the algorithms and software implementations were not designed for scalable, parallel computation. A scalable algorithm is one whose computational complexity (i.e., the number of operations to reach solution) is proportional to the number of unknowns; moreover, the algorithm should also have a convergence rate that is independent of problem size or the number of processors. In numerical solution algorithms, there is often a tradeoff between convergence rate and degree of parallelism. As a result, to obtain a useful measure of parallel efficiency, the best scalar (uniprocessor) algorithm should be used as reference. Scalable methods are needed because the size of problems of interest continues to grow quite significantly, and we want to avoid methods with computational complexities of O(Na) with an athat is (much) larger than unity.

This paper discusses the application of preconditioned generalized conjugate gradient acceleration to fully implicit thermal simulation. The preconditioning step utilizes incomplete Gaussian elimination (IGE) to form an approximate factorization of the Jacobian matrix.

It is well known that the adjoint approach is the most efficient approach for gradient calculation, and it can be used with gradient-based optimization techniques to solve various optimization problems, such as the production-optimization problem and the history-matching problem. 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. Preconditioners specialized for reservoir simulation, such as constrained pressure residual (CPR), exploit properties of the Jacobian matrix to accelerate convergence, so they cannot be applied directly to the adjoint equation. To overcome this challenge, we have developed a new two-stage preconditioner for efficient solution of the adjoint equation by adaptation of the CPR preconditioner (named CPRA: CPR preconditioner for adjoint equation).The CPRA preconditioner has been coupled with an algebraic multigrid (AMG) linear solver and implemented in Chevron's extended applications reservoir simulator (CHEARS V R ). The AMG solver is well known for its outstanding capability to solve the pressure equation of complex reservoir models; solving the linear system with the "transpose" of the 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 generalized minimal residual (GMRES) acceleration solves the adjoint equation very efficiently with a reasonable number of linear-solver iterations. Adjoint simulations to calculate the gradients with the CPRA/AMG solver take approximately the same amount of time (at most) as do the 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. Mathematical Derivation of the CPRA PreconditionerDefinition of the Optimization Problem. The optimization problem discussed in the preceding requires a sequence of control vectors u n be found for n ¼ 0, 1, …, N -1, where n is the control-step index and N is the total number of control steps, to maximize (or minimize) a performance measurement J(u 0 , …, u N-1 ) (Sarma et al. 2005). The

History matching of oil and gas reservoirs can be accelerated by directly calculating the gradients of observed quantities (e.g., well pressure) with respect to the adjustable reservoir parameters (e.g., permeability). This leads to a set of linear equations which add a significant overhead to the full simulation run without gradients. Direct Gauss elimination solvers can be used to address this problem by performing the factorization of the matrix only once and then reusing the factored matrix for the solution of the multiple right hand sides. This is a limited technique, however. Experience has shown that problems with greater than few thousand cells may not be practical for direct solvers because of computation time and memory limitations. This paper discusses the implementation of a multiple right hand side iterative linear equation solver (MRHS) for a system of adjoint equations to significantly enhance the performance of a gradient simulator. This was accomplished by saving the inverse of the preconditioning matrix M from the first gradient solution and utilizing it in subsequent solutions. Reduced System/ILU(O) was chosen for the gradient solution since this preconditioner is one of the more robust. In this way, the significant computation time required to calculate the preconditioner is only needed for the first gradient solution, while all subsequent gradient calculations take only a fraction of the CPU time of the first solution. Numerical experiments were conducted on models of up to 10000 cells and the new MRHS iterative solver was compared with standard red-black line successive overrelaxation (RBLSOR) and Direct(D4) solvers. For the problems which could be compared, overall speedups with the new MRHS solver varied from a factor of 3-6 over other solvers while gain in storage was about factor of 6 over the direct solver. Introduction Before simulating the future performance of a reservoir it is conventional to history match the model. History matching is the process by which a numerical reservoir model is validated by modelling the reservoir S past performance. The results of the model, for example well pressures, water cuts and gas-oil ratios, are compared with the measured values of these quantities. If the differences are acceptably small and the reservoir description is reasonable then the model is valid and is said to be history matched. If the differences are too large then the model is adjusted until a history match is obtained. Quantification of the error is achieved by defining an objective function which is the sum of the squares of the differences between calculated and measured variables, such as a shut-in bottom hole pressures. The objective function is minimized using a non-linear regression algorithm. This requires the derivative, or gradient, of the calculated variables with respect to the history matching parameters. P. 249