Many enhanced oil recovery (EOR) processes can be described using partial differential equations with parameters that are strongly non-linear functions of one or more of the state variables. Typically these nonlinearities result in solution components changing several orders of magnitude over small spatial or temporal distances. The numerical simulation of such processes with the aid of finite volume or finite element techniques poses challenges. In particular, temporally oscillating state variable values are observed for realistic grid sizes when conventional discretization schemes are used. These oscillations, which do not represent a physical process but are discretization artifacts, hamper the use of the forward simulation model for optimization purposes. To analyze these problems, we study the dynamics of a simple foam model describing the interaction of water, gas and surfactants in a porous medium. It contains sharp gradients due to the formation of foam. The simplicity of the model allows us to gain a better understanding of the underlying processes and difficulties of the problem. The foam equations are discretized by a first-order finite volume method. Instead of using a finite volume method with a standard interpolation procedure, we opt for an integral average, which smooths out the discontinuity caused by foam generation. We introduce this method by applying it to the heat equation with discontinuous thermal conductivity. A similar technique is then applied to the foam model, reducing the oscillations drastically, but not removing them. @en