Solver Preconditioning Using the Combinatorial Multilevel Method

The purpose of this paper is to report the first preliminary study of the recently introduced Combinatorial Multilevel (CML) method for solver preconditioning in large-scale reservoir simulation with coupled geomechanics. The CML method is a variant of the popular Algebraic Multigrid (AMG) method yet with essential differences. The basic idea of this new approach is to construct a hierarchy of matrices by viewing the underlying matrix as a graph and by using the discrete geometry of the graph such as graph separators and expansion. In this way, the CML method combines the merits of both geometric and algebraic multigrid methods. The resulting hybrid approach not only provides a simpler and faster set-up phase compared to AMG, but the method can be proved to exhibit strong convergence guarantees for arbitrary symmetric diagonally-dominant matrices. In addition, the underlying theoretical soundness of the CML method contrasts to the heuristic AMG approach, which often can show slow convergence for difficult problems. This new approach is implemented in a reservoir simulator for both pressure and displacement preconditioners in the multi-stage preconditioning technique. We present results based on several known benchmark problems and provide a comparison of performance and complexity with the widespread preconditioning schemes used in large-scale reservoir simulation. An adaptation of CML for unsymmetric matrices is shown to exhibit excellent convergence properties for realistic cases. Introduction Reservoir simulation, which mimics or infers the behavior of fluid flow in a petroleum reservoir system through the use of mathematical models, is a practice that is widely used in petroleum upstream development and production. Reservoir simulation was born as an efficient tool for reservoir engineers to better understand and manage assets. However, like any numerical simulation tool, reservoir simulation is inherently computational intensive and easily becomes inefficient if more grids, coupled physics, and/or complex geometry are necessary to accurately describe the complex phenomena occurring in the subsurface. Mathematically speaking, reservoir simulation solves a system of discretized partial differential equations (PDEs) which describe the underlying physics. Due to stability constraints, an implicit formulation is required at least for the pressure system. Details about the numerical analysis for choosing an implicit formulation (or more specifically, the backward Euler method) can be found in the classic literature of Aziz and Settari (1979). However, as a recent exception, Piault and Ding (1993) attempted a fully explicit scheme in a reservoir simulation on a massively parallel computer and showed acceptable results. They adopted the Dufort and Frankel scheme which is unconditionally stable but numerically inconsistent (Dufort and Frankel 1953). This scheme is of order of 2 2 x t   accuracy, which clearly implies the truncation error can be significant if t  does not approach 0 faster than x  . In essence, implicit formulation is the only unconditionally stable and consistent scheme and is adopted by all commercial reservoir simulators. As a result, a linear solver is inevitable for reservoir simulation due to this implicit formulation. There are four main streams of formulations applied in reservoir simulation: IMPES, fully-implicit, AIM, sequential implicit. Of these, fully-implicit is the most robust formulation but the resulting coupled system matrix is numerically challenging and computationally expensive. In the fully-implicit formulation, pressure, saturation/mass, and/or temperature are to be solved simultaneously. The generated system matrix is highly non-symmetric and not positive definite, which brings great challenges for applying robust and efficient preconditioners and liner solvers. This situation is further exacerbated for large-scale models with highly heterogeneous coefficients and unstructured gridding. Since, generally speaking, in black-oil simulation the solution of linear system ( b Ax  ) usually consumes up to 90% of the total execution time, linear solver performance enhancement means significant reservoir simulator speedup.