A New Front-Tracking Method for Reservoir Simulation

Abstract

Summary This paper presents a new numerical method for solving saturation equationswithout stability problems and without smearing saturation fronts. A reservoirsimulator based on this numerical method is under development. A set of testproblems is used to compare the simulation results of the new simulator withthose of an existing finite-difference simulator (FDS). Introduction The standard nonlinear, partial-differential equations that describe flow inporous media can be separated into a pressure equation and saturationequations. If the diffusion term is ignored, the saturation equations describea physical problem where sharp discontinuities in the physical data arepossible. Finite-difference methods used to solve these equations typicallyexhibit numerical dispersion. They also show numerical stability problems sothat very short timesteps may be required. Use of an implicit problems so thatvery short timesteps may be required. Use of an implicit formulation can reducethis limitation of the timestep length, but this will not solve the numericaldispersion problem. Different methods are suggested in the literature for representing thesaturation fronts more accurately. These methods usually center on the methodof characteristics and usually are limited to miscible flow problems. Glimm etal. developed a front-tracking simulator (FTS) based problems. Glimm et al. developed a front-tracking simulator (FTS) based on a method that follows thediscontinuities separately and solves the saturation equation in the remainingpart of the reservoir with a infinite- difference method. New methods in the field of hyperbolic conservation laws have led toalternative solution procedures for the saturation equations. A simulator basedon an implicit pressure, explicit saturation (IMPES) formulation and thesemethods is under development. The goal of this development work is a 3D, three-phase simulator. In this reservoir simulator, the pressure equation issolved by a infinite-element method (FEM). The grid for the pressure equationcan therefore be fitted to the reservoir geometry and the pressure equation cantherefore be fitted to the reservoir geometry and the geometry of the sharpdiscontinuities in the saturations with great flexibility. The linear equationsystem is solved with a preconditioned conjugate-gradient method. The method for solving the saturation equations is based on an approximationof the fractional-flow function by a piecewise linear function. This methodallows the saturation equations to be solved without stability problems andrestrictions on the timestep length caused by the Courant-Friedrichs-Levy (CFL)condition. A set of test problems has been used to compare the simulation results ofthe new simulator with those of an existing FDS. The test cases show the FTS tobe without the grid and numerical dispersion effects observed in a standardfive-point FDS- The FTS is also computationally more efficient thanfinite-difference methods. The simulator has been used in North Sea fieldsimulation problems. Front-Tracking Method Traditional numerical methods have a tendency to smear discontinuities owingto numerical dispersion. This results in a less accurate estimate of theposition of the fluids and their motion in the reservoir. By introducing thediscontinuities (fronts) as separate entities, their positions can be simulatedmore accurately. The crucial point is to positions can be simulated moreaccurately. The crucial point is to represent the physical discontinuities(saturation fronts) correctly. In front-tracking, the sharp change insaturation is replaced by a step function-i.e., a discontinuous function. Thissaturation jump is called a front, and the method handles the frontsseparately. In two space dimensions, each front is represented as a piece-wise linearcurve. The speed of these fronts can be computed, and the saturation equationcan be solved away from the fronts by a standard numerical method based on agrid, as proposed in Ref. 6. The new numerical method in this paper is based ona piecewise linear approximation of the flow function (flux function) definedin Eq. 4. This approximation leads to a saturation profile that is piecewiseconstant. The saturation therefore is represented by a step function with a setof discontinuities (fronts). A physical saturation profile normally will becontinuous, with some sharp transition zones that can be called realdiscontinuities (fronts). Our numerical method therefore distinguishes itselffrom other front-tracking algorithms in that the saturation consists solely ofa step function with a set of discrete fronts. As the number of frontsin-creases, however, the solution converges toward a continuous solution in thedomains where the physical solution is continuous. physical solution iscontinuous. A sound mathematical foundation for the solution of the numericalproblem has been developed that leads to the solution of the Riemann problemhas been developed that leads to the solution of the Riemann problemsimplemented in the simulator. The independence of a grid for the problemsimplemented in the simulator. The independence of a grid for the solution ofthe saturation equation also results in a method that is not limited by the CFLcondition for the length of the timestep. Because the saturation equation often can be seen as locally ID (flow alonga characteristic), a mathematical theory for the solution of this kind ofproblem is a first building block in the theory. A more thorough explanation isgiven later under 1D Riemann Problem. The mathematical results have been generated to higher space dimensions. Amore thorough explanation is given under Riemann Problem in 2D and 3D. Refs. 10through 12 are the basis of the new reservoir simulator described in thispaper. When several phases are present in the reservoir, the saturation equation isa system of hyperbolic conservation laws. The mathematical theory for systemsof conservation laws has not been sufficiently developed to study multiphaseflow in porous media. A particular model for three-phase flow has beenproposed, however, and is discussed under Triangular Model. Several new problems arise when the method is applied to real data. At theboundary between two geological layers, the flux/fractional-flow function maychange because of gravity segregation. This problem is solved theoretically andis currently being implemented in the reservoir simulator. The background forthis solution is explained under Triangular Model. Equations The equations that govern the flow in a porous media, Darcy'slaw and themass-conservation equations, are found in the Appendix. By combining Darcy'slaw and the mass-conservation equations, we get a set of coupled, nonlinearpartial-differential equations. These can be separated into a pressure equationand a saturation equation. The well terms q, and q. are pressure equation anda saturation equation. The well terms q, and q. are set to zero. Theproduction/injection terms are modeled as a boundary condition on an internalboundary. We also set the capillary pressure, P, to zero so that the phasepressures are equal to the average pressure. We assume a black-oilformulation.