Tracking Thermal Saturation Fronts by a High Level PC Programming Language
Authors: Erkal A. , Numbere D.T. 
Published: December 2007

Abstract
This paper presents a PC based alternative procedure for determining the water saturations within the hot water zone of a thermal project for use in analytical oil recovery calculations. Conventional analytical calculation of oil recovery under steam and hot water injection requires the tracking of the movements of the saturations advancing within the variable-temperature hot-water zone. This involves an adaptation of the Buckley-Leverett theory to this variable temperature zone after dividing it into a number of constant temperature or isothermal zones. If the number exceeds two, the calculation can become very tedious unless done with a computer program. FORTRAN programs have generally been used but they are not as intuitive or easy to use as modern PC based programming tools such as MathematicaTM. In this paper. MathematicaTM was used because it is relatively easy to program. easy to use, and is fast and robust. Unfortunately. depending on the number of isothermal zones, the time step size and the time at which the recovery calculation is desired, the tracking of saturations can tax the capabilities of even the most modern PCs. Therefore, an alternative method is introduced that eliminates the need for saturation tracking. This method calculates the instantaneous saturations within each isothermal zone at the time of interest. Oil recovery results by this method were found to be comparable to those by the saturation tracking method with considerable saving in computation time. Two examples are presented to demonstrate the utility of the method.

Introduction
For one-dimensional models, the analytical calculation of steamflood oil recovery with time requires calculating the following:the position of the steam front and its rate of advance or velocity.the temperature profile in the hot water zone and its rate of advance,the fluid saturation profile and its rate of advance.
The location of the steam front at any time can be calculated by the equations of Mandl and Volek who showed that before a critical time tc, the steam zone can be described by the equations of Marx and Langenheim. Beyond tc, Mandl and Volek present an approximate equation to calculate the steam front location with time. This equation was later improved by Prats and Vogiatzis and communicated to Myhill and Stegemeier who presented the solution in graphical form.
The temperature profile in the hot water zone can be calculated using the equation by Lauwerier for both hot waterflood and steamflood. However, for a steamflood, a simpler approximation assumes a linear temperature drop from the steam temperature to the cold reservoir temperature. The saturation profile in the steam and hot water zones can be calculated using the Buckley-Leverett theory provided that the non-isothermal hot water zone is first divided into an appropriate number of isothermal zones. Several procedures are available to calculate the saturation profile with time. This paper is concerned with the method whereby characteristic saturations are picked and their motions tracked. The tracking process is very time consuming and can easily exceed the capabilities of PCs' when small time steps are used together with many isothermal zones at long times. An alternative to tracking is used whereby instantaneous saturations in each temperature zone are calculated at any given time.
The saturations to be tracked can be picked in different ways. Farouq Ali determines the flood front saturation at the cold reservoir temperature (just as in a waterflood) and arbitrarily picks saturation values greater than this value ending at the maximum water saturation. This has the advantage of giving greater saturation definition where needed. However, it has the potential of using more saturations than is needed leading higher computation time on personal computers. Willman et. al. suggest a somewhat narrower range of saturations by finding as the starting point, that saturation in the cold zone that have the same velocity as the cold temperature at tc. Prats extends this further by finding a particular saturation for each isothermal zone referred as characteristic saturations. These are calculated as the saturation in each zone that has the same velocity as the cold temperature at tc. This way, there are only as many characteristic saturations to be tracked as there are isothermal zones as opposed to tracking an arbitrarily large number of saturations.
Even with these improvements in the choice of saturations, the tracking process of calculating their locations at each time step is time consuming even for personal computers.
P. 83