Solving multiphase equations from a dynamical systems point of view

We discuss the numerical solution of the partial differential equations (PDEs) that govern subsurface multiphase flow in porous media from a "dynamical systems" point of view. We study a two-phase problem in one spatial dimension (Green et al., 1970), discretized both temporally and spatially by finite differences. We use the formulation in Peaceman (1977), where at each time step we solve first a linear "pressure equation" directly, then a nonlinear "saturation equation" via an iterative process. We consider two iterative processes for the saturation equation: Newton's method and Picard linearization. Our present work focuses on the case where there are only two unknowns in the discretized saturation equation for which to solve at each iteration, a selection that simplifies analysis while permitting pertinent dynamics to be observed. Due to the strong nonlinearity in the saturation PDE, very interesting phenomena are observed during the iterative processes. For the Picard iterations, stable n-cycles for many different values of n are observed. Also, we find at each time step that there is often more than one solution for the discretized nonlinear saturation equation. To which of these solutions the Newton iterates converge depends, in quite a subtle manner, on the "initial guess" for the solution.