A New Method for Solving the Implicit Reservoir Simulation Matrix Equation

Abstract

Summary The dominant computational cost associated with implicit reservoir simulation arises in solving the matrix equation at each Newton iteration. This paper describes a new solution method based on a reduced matrix equation constructed using singular value decomposition. The reduced matrix equation is used as the basis of a preconditioner. The resulting preconditioned iteration generally converges in two to four iterations. The new solver can reduce solve times substantially. Introduction Reservoir simulation proceeds in a series of timesteps, with each timestep requiring one or more Newton iterations. Generally, the most time-consuming operation in a Newton iteration is solving the matrix equation. In IMPES1,2, this solution generally accounts for between one-third and one-half of the total solution time. In sequential implicit calculations3, the same is true. However, in fully implicit calculations, the matrix equation solve is much more expensive. In fact, it often dominates the total computer time, costing several times as much as all other computations combined. This is the usual case in compositional models. The solution of the fully implicit matrix equation takes more time than solving the IMPES matrix equation because it is larger, containing more unknowns and matrix coefficients. The equation set has an unknown and mass balance equation for each component being modeled plus a pressure unknown paired with a volume constraint equation4. One of the unknowns can be eliminated easily, leaving a total of one unknown per component. In the most commonly-used black-oil formulation, there are three components. In compositional models, there typically are around ten components. The solve cost grows faster than linearly with the number of unknowns at each gridblock. If the method used to solve the IMPES matrix equation is simply generalized to handle multiple components, solving for three components per gridblock might cost roughly an order of magnitude more than solving the single-unknown IMPES matrix equation, and solving for ten components per gridblock might cost roughly an order of magnitude more than solving for three. This paper describes a way to reduce the implicit matrix equation such that it has far fewer implicit unknowns and coefficients than the original matrix equation but yields nearly the same solution. In a typical model, the reduced matrix equation has one unknown (pressure) at most gridblocks, two or three unknowns at most remaining gridblocks, and approximately the original number of unknowns at the few remaining gridblocks. Because of the reduction in the number of unknowns and the accompanying reduction in the number of matrix coefficients, the solution takes a fraction of the time that it would take to solve the full set of equations. The reduced matrix equation is used as a preconditioner in an iterative solution of the original matrix equation. The reduction in unknowns is based on singular value decomposition, or SVD5. Specifically, SVD is used to decompose the flow submatrices in each node?s column within the global matrix. These submatrices govern the effect of this node's solution variables on flow between it and its neighbors. They are the off-diagonal submatrices in the node?s column and the diagonal submatrix after subtracting the capacitance terms out of it. In general, flow between two nodes depends on pressures and masses in the two nodes. The result of the SVD indicates that the effects of masses on this flow can often be approximated accurately using significantly fewer variables than there are masses. The reduction in unknowns is obtained by dropping (i.e., setting to zero) the smallest singular values. The accuracy of the preconditioner depends on the singular values dropped at each node. The size of the largest dropped singular value is a reasonable indicator of this accuracy. This makes it possible to choose an appropriate threshold singular value size. An unknown and equation (in addition to the pressure unknown and corresponding equation) are retained for each singular value larger than the threshold. When an appropriate threshold is used, two to four solutions of the preconditioner are generally sufficient. When only a few singular values are needed, the precise number needed does not depend strongly on how many components there are. For example, if one singular value is needed, this is generally true whether the number of components is three or twenty. As a result, the more components there are, the more benefit this reduction should provide.