Fluid Phase Equilibria, 94 (1994) 93-108 93

Elsevier Science B.V.

Non-iterative Flash Calculation Algorithm in

Compositional Reservoir Simulation.

Peng Wang and Erling H. Stenby’

Engineering Research Center, IVC-SEP, Department of Chemical Engineering

The Technxal University of Denmark, Building 229, DK-2800, Lyngby

Keywords: flash, equation of state, reservoir simulation, EOR

(Received March 22, 1993; accepted in final form December 29, 1993)

ABSTRACT

This paper describes an efficient two- and three-phase flash calculation algorithm
for equation of state (EOS) compositional reservoir simulation problems. In this
approach, two sets of linear equations with the same coefficient matrix are established
on the basis of the phase equilibrium information at the previous time step. The
solution of these equations is the distribution of mole number of each component
between phases. The phase saturation can thus be obtained without solving the
Rachford-Rice equation.
A number of reservoir simulations have been performed to show the validity of the
new algorithm. Comparisons are made of these runs with those conducted by a
rigorous fugacity iteration method. No apparent difference in the results from the two
methods is observed. However, the CPU time can be significantly saved by using the
new approach, generally around 50%.
In this work an attempt is also made to employ the non-iterative approach as a tool
for providing the initial guess for the rigorous EOS iteration. A considerable reduction
of the CPU time is again observed.

INTRODUCTION

A reservoir simulator constructed on the basis of Darcy’s law is actually used to

model phase flow behavior for processes taking place in a reservoir. During
the simulation, the reservoir is described by many small volume elements in one, two
or three-dimensional grids, dependent on the complexity of the problem dealt with.
Each of the small volume elements is characterized by the proper reservoir rock

* Author to whom correspondence should be addressed.

0378-3812/94/$07.00 01994 Elsevier Science B.V. All rights reserved

SSDI 0378-3812 (94) 02470-L
94 P. Wung, E.H. ~tenby~F~~i~ Phase E~ui~ib~i~ 94 (1994) 93-108

properties and the displacement process is demonstrated by the variations of the phase
properties and the interaction between rock and fluid with time in these small blocks.
The phases in each of the volume elements are conserved by means of components.
Therefore, it is necessary to determine how the components are distributed among the
phases.
It is commonly assumed that, during the simulation, a mass transfer between phases
takes place much faster than the phase flow in the reservoir rock. A series of
thermodynamic equilibria for each volume element is established at any time of the
simulation. This leads to the phase equilibrium calculation being required at least once
for each time step and each volume element in the grid.
For the simulation of enhanced oil recovery (EOR) processes, such as miscible and
immiscible gas injection processes, a compositional simulator is normally used, in
which an equation of state (EOS) approach is commonly employed for the
determination of phase equilibrium. Many studies have proved that for a typical
compositions simulation the largest portion (-70%) of the total execution time is
often spent on phase behavior computation. This has restricted the application of the
compositional model to a full-field reservoir simulation. A robust, efficient and fast
phase equilibrium calculation has been urgently demanded to be developed. Great
progress has been made in this field in the last decade. So far it may be said that the
determination of the phase number and the convergence of the solution method for
the hydrocarbon mixtures at the given temperature and pressure are no longer the
problem (Micheisen, 1982; Trangenstein, 1987). However, even if the composi~onal
simulator employs the most advanced phase equilibrium calculation routine, the
computation relevant to the flash calculation still takes a significant percentage of the
total CPU time. This is mainly because the solution is still obtained by iteration.
The purpose of this work is to find a non-iterative approach to get the dis~ibution
of the component concentration between phases in the EOS reservoir simulations. In
this study it is shown that the phase equilibrium calculations can be conducted by
solving two sets of linear equations formed on the basis of the PVT information at
the previous time step. The coefficient matrices of these two sets of equations, which
are partial derivatives of component fugacity with respect to the mole number of the
component, are the same and decomposition is thus required once only. Furthermore,
there is no need for solving the Rachford-Rice equation for the phase saturation in
this procedure, since the solution to the linear equations is the mole number of each
component in each phase.

DESCRIPTION OF THE NON-ITERATIVE APPROACH

Let us consider an isothermal, two-phase flash case in one block. Information such
as pressure and mole number of each component in both phases is known at the
previous time step. The task of flash calculation is just to determine the component
partitions between phases at the current time step. The differences of the conditions
over a time step are the pressure and the total mole number of each component,
P. Wang, E.H. StenbylFluid Phase Equilibria 94 (1994) 93-108 95

Thus the following equation may be written to account for these two effects on the
component partitions:

dns = ( 3$&J t

( i = 1, 2, .... N, )

The mole number of each component in a phase at the current time step can,
therefore, be readily obtained from:

n&v = ncold + dn” (2)

The phase molar saturation in absence of water is:

$ny_
s, = N i=lN (3)

g &L4
a=1 i=l

In the simulations, dP and dn, are already known before the flash. The derivatives
in Eq.(l) may be evaluated in terms of the following simple way.
For the fluid in equilibrium at the previous time step, the below two equations must
be obeyed:

I&-lnf,‘=O (4)

ntL+n,“=ni (5)

Taking the partial derivatives of the above two equations with respect to the
pressure, one may reach:

(6)

(7)
96 P. Wang, E.H. StenbylFluid Phase Equilibria 94 (1994) 93-108

Substituting Eq.(7) into Eq.(6), one may readily get:

This is a set of N, linear equations with N, unknowns, &iL/L3P.By solving these

equations, the impact of the pressure variation on the component partition can be
estimated. The contribution of variations of the total mole number of each component
to the component partition may be evaluated in a similar way.
Considering the partial derivatives of Eq.(l) with respect to the overall mole
number of component i, rearranging the resultant equation, we may yield:

(9)

( k = 1, 2, .... N, )

ani’ an: (10)

(~)P=+k&p
J J

where 6,j=1 fori=j;

6,, = 0 for i f j.

It appears that Eq.(9) is also a set of linear equations, but with N, x N, equations
and N, x N, unknowns. The right-hand sides of Eqs.(S) and (9) are easily estimated
from the EOS. The partial derivative of the component fugacity against component
mole number can also be analytically evaluated. Besides, it is obvious that the
coefficient matrices of Eqs.(S) and (9) are exactly the same. The decomposition is
therefore required once only in the solution.
When we know the mole number of each component in different phases, the mole
composition can readily be computed from:

(11)

The new Xi and y, are used to compute the component fugacity and the partial
derivatives for the flash in the next time step.
For an isothermal three-phase flash, one more equation is respectively added into
Eqs(8) and (9). The resultant equations are required to be solved simultaneously.
P. Wang, E.H. StenbylFluid Phase Equilibria 94 (1994) 93-108 97

It appears that there is no iteration involved in the new procedure except in the
initial simulation in which the EOS iteration is still required to perform the phase
stability analysis to determine the phase number and the distribution of component
mole number between phases in each small volume element. It is also noteworthy that
the computation burden in this procedure comes mainly from the estimation of the
derivatives. In the volume-balance EOS simulators, however, these derivatives are
already obtained through the solution of the pressure equation.

DETERMINATION OF PHASE NUMBER

When using this new algorithm, we have to determine the phase number of the
mixture with current total composition and pressure. For a two-phase mixture, this
problem has been well addressed in the literature when the K-values are known. Here
we describe how to determine the current phase number for a previous three-phase
mixture.
An advanced approach is the Gibbs free energy minimization technique. This
method allows the phase number to be directly addressed by analyzing the Gibbs
energies of the mixture at different phase numbers. Therefore, in order to determine
a three-phase mixture the Gibbs energy of that mixture at two- and single-phase state
is usually required. Then a third phase is dropped out of the system to see whether
it could reduce the energy of the system. If it does, a three-phase flash is conducted.
Nelson (1987) proposed a relatively simple approach of doing so, which is based on
a multi-phase formulation of the well-known bubble and dew-point criteria and
applied to the cases in which a good initial guess of the K-values are available. There
is an iteration procedure involved in this approach.
In this work, we use the mole number of each phase to determine the phase number
under the current conditions. The phase with a phase mole number less than zero will
be dropped. If nnewv < 0, for example, the system will be in liquid-liquid state at the
current time step. If mole numbers of two phases are below zero, the system will be
in single-phase region with the state of phase whose mole number is larger than zero.

VALIDATION OF THE NEW METHOD

We first test the abilities of the new approach to simulate differential liberation
(DL) and constant volume depletion (CVD) experiments. During the simulations, the
rigorous EOS iteration method is used to perform the phase equilibrium calculation
at the first pressure stage. In the rest of the stages, the new method is switched on.
Many examples including gas condensate and crude oil have been tested. The
typical results are shown in Figures 1 and 2. In Figure 1 is shown the K-values
computed corn the new procedure and the EOS iteration for a West Texas oil (Lim
et al., 1992). An estimation of phase molar saturation for this example is also made
using the new approach, Eq.(5), and the EOS iteration. The comparison is shown in
Figure 2.
98 P. Wang, E.H. StenbylFluid Phase Equilibria 94 (1994) 93-108

From these two plots it is apparent that the proposed approach can well match the
phase equilibrium calculations from the rigorous iteration of the EOS. No distinct
difference between the K-values computed from the two methods is observed.
Besides, the phase saturations at the different pressures can also be accurately
estimated from Eq.(5) for all the examples tested. However, the Newton-Raphson
iteration of this property in the Rachford-Rice equation is totally avoided in the new
approach.

APPLICATION TO THE RESERVOIR SIMULATIONS

The proposed non-iterative flash approach has been implemented into an EOS
compositional reservoir simulator, UTCOMP (Chang, 1990). The simulations of
water-alternating-gas (WAG) cycling of a oil reservoir and dry gas cycling of a gas
condensate reservoir have been conducted, in which only two-phase flash is required.
The new three-phase flash is validated in terms of simulating three CO,-injection
processes in which three-hydrocarbon-phase with water flow occurs. During the
simulation, the PR EOS is used to carry out the flash at the beginning of the
simulation. The new approach is then switched to perform the flash calculations for
two- and three-phase blocks in the rest of the simulation period. A rigorous phase
stability analysis, which is based on the PR EOS for determination of the phase
number, is carried out for the only single-hydrocarbon-phase blocks adjacent to the
multi-phase blocks at the previous time step. While the single-hydrocarbon-phase
blocks remote from the multi-phase blocks are assumed to keep their phase number
at the current time step.

l **** PR model
- This work

- This work
l **=* PR model

0.007 ,,,,#(I,, (, l,,,,,,,,,,,,,,,l,r

80 85 90 95 100 105 1 0
Pressure (atm)
Figure 1X-value computed from the
Figure 2. Molar saturationof gas
PR EOS and the new approach for a
phase computed from the PR EOS and
West Texas oil.
the new approach for a West Texas oil.
P. Wang, E.H. StenbylFluid Phase Equilibria 94 (1994) 93-108 99

TWO-PHASE CASE

Slim-Tube Simulation

A slim-tube displacement experiment is used to determine the minimum dynamic

miscibility pressure (MMP) for the system containing oil in place and injected
solvent. At the pressure around MMP, the fluid in the tube may be considered in the
critical region. Thus, the simulations of this displacement can probably be used to test
whether the proposed flash procedure is capable of handling the critical phase
behavior.
The fluid of the SPE fifth project (Killough and Kossack, 1987) is selected for this
test. During the simulation, the slim tube is represented by a 30 grid-block system
with an identical length of 2.0 ft. Figure 3 shows the oil recoveries computed from
the new method and the PR EOS iteration at the different pressures when 1.2 standard
pore volume (PV) solvent has been injected. The MMP of about 3000 psi is found
by the both runs. The oil recovery curve with 1.2 PV gas injection against the
pressure also agrees well from these two flash approaches.

IFifth

The purpose of this project is to test the abilities of the compositional model to
simulate the WAG injection process into a volatile oil reservoir. Three different WAG
cycling and production, injection scenarios are designed for this comparative project.
One three-dimensional 7 x 7 x 3 grid system is used for all runs. Detailed information
on this project is given in the original paper.
These three cases have been simulated by UTCOMP with the new flash approach.
Figures 4 to 6 show the results for case 1. It appears that the non-iterative flash
routine leads to results almost identical to those obtained from the PR EOS iteration.
105
.*.**PR model
_ - This work
Cv 95-
rl -

120 140 160 180 200 220 240 26

Pressure (atm)
Figure 3. Slim-tube simulation from Figure 4. Comparison of average reservoir
the PR model and the new approach. pressure for case one of the fifth SPE proJect.
100 P. Wang, E.H. StenbylFluad Phase Equilibria 94 (1994) 93-108

,
*****PI2 EOS

0 3 6 9 12 15 18 O”” Time (years)

B
Time (years)
Figure 5 Comparison of producmg GOR Figure 6 Comparison of 011recovery for case
for case one of the fifth SPE project. one of the fifth SPE project

For the other two cases, we observed the same degree of agreement between the
run with the new approach and that with the EOS iteration.
A displacement of dry gas cycling of a gas condensate is also selected to show the
validity of the new approach. The detailed specifications of this process is given in
elsewhere (Wang and Stenby, 1993). The simulation of the example clearly shows an
excellent agreement between the two runs with different flash routines.

Influence of Time Steu Size

In the new approach, we assume that the partial derivatives involved in Eqs.(l) and
(2) are maintained constant within each time step. An investigation of the influence
of the time step size on the simulation accuracy should be made in order to verify this
assumption. The fifth SPE project is selected for this test. It is found that the results
obtained from the runs with different sizes of time step are not different. Both runs
agree well with the one by use of the rigorous iteration of the PR model.

THREE-PHASE CASE

Three-hydrocarbon-phase equilibria have been experimentally observed at

temperatures below about 120 “F in some oil-CO, systems (Orr et al. 1981; Turek et
al. 1984; Larsen et al. 1989). In order to predict the reservoir performance in a CO1
flooding process accurately, a simulator which is able to deal with four-phase flow
(water is always present in the reservoir) must be available. Apparently, three-phase
flash calculation is the core of such simulator.
Several CO, displacements have been simulated by the simulator, UTCOMP, with
the new three-phase flash routine and determination of the phase number. Some
typical results will be shown in the following section.
P. Wang, E.H. StenbylFluid Phase Equilibria 94 11994) 93-108 101

Table 1 lists an oil composition and relative permeability function used in this test.
Figure 7 shows the P-z diagram for this oil/CO, mixture at 100°F computed by the
PR model. The appearance of a third phase is found at pressures between 1100 and
1170 psi and CO, mole fraction higher than 75%.
The sample reservoir is 100 ft long, 100 ft wide, 10 ft thick and is represented by
a 11 x 11 grid with equal-size block in the X and Y directions. The initial reservoir
pressure is specified at 1127.5 psi. Reservoir temperature is 100°F. Pure CO, is
injected at one comer of the reservoir with a constant bottom hole pressure of 1127.5
psi. A constant production rate of 12 Bl-moles/day is applied to the production well
at the opposite comer. Water with a saturation of 47% is initially in place but
immobile in the later simulation.

TABLE 1. Characterization of an Amarume Stock Tank Oil

%Cal 0.1140 250.98 1365 00 20.624 0.6185 227.77 600

C 211 0.1199 147.43 1580.00 30.776 0.9171 420.24 700

Relative permeability model:

S, - 0.16 si
k,,,,=O ; k, = 0.574( I3 ; krj=0.278(o)
6.89
o’=V, LJ
0.37 .
Parameters of Leverett’s capillary pressure function (Leverett, 1941):

c =68psiaJj;;;i ; E =2 .
’ (Jwo= 42 dyne/cm
P - P
dyne[cm

The maximum simulation period of 512 days is specified. A limiting GOR of 50

Mscf/STB is used for the shut-in criteria. Figure 8 shows the average phase saturation
history simulated by the new approach and the EOS iteration. It can be seen that he
results from these two methods agree excellently.
The three-hydrocarbon-phase equilibrium is developed shortly after the
displacement. Figure 9 also show a very close match of the producing GOR estimated
from both approaches for this test case.
102 P. Wang, E.H. StenbylFlurd Phase Equilibria 94 (1994) 93-108

Distribution of three-hydrocarbon-phase blocks at three different times are

demonstrated in Figure 10. The phase number detection technique suggested in this
work yields the same result as estimated by the rigorous procedure of the Gibbs free
energy minimization. No three-phase blocks are missing during the simulation with
the non-iterative flash routine.

0 100 200 300 400 500

CO2 mole fraction Time (days)

Figure 7. P-z diagram for a CO,/oil system Figure 8. Average phase saturation from the
displaying VLLE at 100°F. EOS iteration and the new approach.

1 ’ I 1 I ’ ’ I-
0.0 0.2 0.4 0.6 0.0 1
PV injection
Figure 10. Distribution of three-phase
blocks from the EOS iteration and the new
Figure 9. Producing GOR from the EOS
approach.
iteration and the new approach.
P.Wang,E.H.StenbylFluidPhaseEquilibria
94(1994)93-108 103

COMPARISONBETWEEN THREE AND FOUR-PHASE SOLUTIONS

Since the computational expense with a four-phase flow simulation is at least

doubled as compared to the one of a three-phase solution, one may prefer to make the
runs with a three-phase flow option if there is only acceptable deviation between these
two solutions. Figures 11 and 12 show a comparison of the ultimate oil recovery
between three and four-phase simulations for 1-D and 2-D cases studied in this work.
In the three-phase solution, we ignore CO,-rich liquid phase. For 1-D case, the four-
phase flow simulation results in lower oil recovery during the early period of the
production, but the ultimate oil recovery is around 10% higher than in the three-phase
solution. However, the oil recovery predicted by these two simulations are matched
fairly well in the 2-D case for most of the project life. A 4% difference in the
ultimate oil recovery by both runs is observed.
With regard to the computational expense, simulations with the four-phase flow
option for the 1-D and 2-D runs of this work require an expense of 2.63 and 2.32
times as much, respectively, as in the case of the three-phase solution.

1.0 0.7 _p
- - Four-phase flow
: - - - Three-phase flow

- Four-phaae flow
. - - - Three-phase flow
0.0 I I I , I I 1, I , I , I I I , I I I , I I ,
0.0 0.2 0.4 0.6 0.0 1.0 1.2
PV injection 0.0 0.2 0.4 0.6 0.8 1.0 1
PV injection
Figure 11. Comparison of 011 recovery
Figure 12. Comparison of oil recovery
between three and four-phase solution for between three and four-phase solution for
1-D simulation. Z-D slmulatlon.

COMPUTATIONAL EFFICIENCY

Our purpose of introducing this non-iterative flash approach is to reduce the

computational burden of the compositional reservoir simulations. The results reported
earlier have indicated that all the simulations performed with the new approach are
confined in satisfactory agreement with those made with the EOS iteration method.
Here we present the experiences of execution time of this approach. The HP 9000/720
computer is used for all the runs.
104 P. Wang, E.H. StenbylFluid Phase Equilibna 94 (1994) 93-108

Two-Phase Simulation

It is found that, for the simulation of the fifth SPE project, the new procedure takes
around half, on average, the CPU time required by the run with the EOS iteration
method. For the gas injection process, a reduction of 30% of the CPU cost is
observed, even in this case the hydrocarbon fluid in around half the blocks in the later
simulation period is in the single-phase region.
If assuming that 65% of the total execution time is required by the flash for the
runs made in this work, we also evaluate the reduction of the CPU time of the flash
calculation by using the new approach. It is found that the computational burden in
this type of calculation can be reduced by around 70% in average for the SPE fifth
project, and 45% for the gas cycling process.

In the simulations with the EOS iteration method, it is usually accepted that the
phase equilibrium information at the previous time step can provide good initial guess
for the current flash calculation. In this work, we have also tried to employ the non-
iterative flash approach as a tool of providing the initial guess for the EOS iteration
method. The results of this test indicate that the computational time by means of this
technique can still be reduced significantly. Besides, it is found that the rigorous
iterative method of the EOS can generally reach the solution within two iterations
because the initial guess from the new approach has already been very close to the
true solution. It can further be inferred that more significant reduction of the CPU
time by use of this technique for the EOS iteration method can be expected for the
simulations performed on vector processors, because the maximum calculation rate
on the vector machines need all the flash calculations to reach the solution in about
the same iterations.
The relation between the number of
components and the execution time in
a compositional simulation has also
QQQDDPR EOS
been investigated in this work. This is u
l This work
done by dividing each component into
two components with equal content and
critical properties. Therefore the
simulation results are independent of
component number, while the
computation time increases as the
number of components increases. Case
4
one of the SPE fifth comparative
solution project is selected to carry out 5
OI 6 7 8 9 10 11 12 13
this study. The curves of the total CPU Component number
time against the component number by
the PR EOS and the new approach are F.lgure 13. Effect of number of components on
illustrated in Figure 13. CPU time for case one of the fifth SPE project
P. Wang, E.H. StenbylFluid Phase Equilibria 94 (1994) 93-108 105

It should be noted that the slope of the curve estimated by the simulation with the
PR EOS iteration is much steeper than the one from the new approach. A slight
influence of the number of components on the CPU time of the run with the new
approach is observed. Therefore, the cost of the compositional simulation with the
proposed approach will be far less than the one with the EOS model for the reservoir
fluid consisting of more components.

Three-Phase Simulation

In this section, we will use the 2-D CO,-injection case to show the saving of the
new three-phase flash approach.
In the simulation, a total number of time steps of 2715 is specified. The rigorous
two and three-phase flashes are performed by 211,344 and 43,083 times, respectively,
for the run with the conventional PR EOS iteration method. When the non-iterative
flash option is chosen, 48% of the rigorous two-phase flashes and 99% of the rigorous
three-phase flashes are saved, respectively. This gives us a benefit of CPU time of
38% compared to the one spent by the PR EOS iteration. It is easy to show that the
longer simulation time, the more saving can be expected when using the new flash
approach.
Finally, it may be worth noting that this new flash algorithm is easy to combine
with any EOS, even non-cubic equations or cubic equations with more complex
mixing rules. In literature some EOS of this type have been proposed to deal with the
complicated phase behavior around the mixture critical point. The CPU time in flash
calculation with these models or mixing rules by using the traditional method could
be increased by a few times compared to the one of the two-parameter cubic models
with the conventional mixing rule (Mathias and Benson, 1986). Much CPU time could
be expected to be saved if the phase equilibrium calculation from these models is
performed non-iteratively.

CONCLUSION

This paper presents a non-iterative phase equilibrium calculation approach for the
compositional reservoir simulation problems. The results show that this new approach
is efficient and robust. Unlike the K-value method in which the phase saturation is
still obtained by iteration, the proposed new routine does not involve any iteration for
the flash calculation. Besides, the new approach is easily implemented in the existing
simulator.
The new flash algorithm has been tested against many different simulation cases.
Satisfactory results are obtained. In general, the non-iterative method gives as good
results as the EOS model iteration. However, the execution time for the runs with the
new procedure can be significantly reduced, usually by a factor of 2. For the
simulation of immiscible displacement processes, more savings of the CPU time can
be expected.
106 P. Wang, E.H. StenbylFluid Phase Equilibria 94 (1994) 93-108

The computation time of the simulations with the rigorous iteration of the EOS for
flash calculations can also be considerably reduced by using the results of the non-
iterative approach as initial guess. In this case, the solution can generally be obtained
within two iterations.
Increase in the CPU time with component number for the simulation with the EOS
iterative method is much faster than that with the new approach. An increase in
component number used is accompanied by more benefit to be gained when using the
new approach.
The influence of time step size on the simulations with the new approach is also
studied in this work. No apparent difference in results is observed, even the size of
time step is enlarged by a factor of 5.
The effect of the four-phase flow may be significant in the simulations of low-
temperature fields in which the displacement pressures lie in the range of the
development of three-hydrocarbon phases.

ACKNOWLEDGEMENTS

The authors wish to thank the following:

Dr. Gary A. Pope and Dr. Kamy Sepehrnoori at the University of Texas, Austin, for
permission to use the compositional simulator UTCOMP. Dr. Claude F. Leibovici, Elf
Aquitaine Production, Pau, France, for his suggestions for testing the non-iterative
flash algorithm. The Danish Technical Research Council for financial support.

NOMENCLATURE

F total mole number of the mixture

fl fugacity of component i
K, equilibrium ratio for component i
MW molecular weight
N, total component number
n, total mole number of component i
P pressure
P, critical pressure
S phase molar saturation
V, critical volume
x, mole fraction of component i in liquid phase
Y, mole fraction of component i in vapor phase
z, total mole fraction of component i

w parachor, dyne1’4cm”‘4/mol
w acentric factor
P. Wang, E.H. StenbylFluid Phase Equilibrin 94 (1994) 93-108 107

Supersmpts

L liquid phase
V vapor phase

Subscripts

new current time step

old previous time step

REFERENCES

Chang, Y.-B., 1990. Development of an equation of state compositional simulator. Ph D.

dissertation, the Umversity of Texas at Austin,

Gautam, R and Seider, W D., 1979. Computation of phase and chemical equthbria: part I. Local
and constrained mmima in Gibbs free energy; part II. Phase-splitting. AIChE J, 25. 991- 999

Inoue, S., Mitsutshi, H , Morii, S. and Nanba, T., 1983. An evaluation of posstble CO, floodmg m
an oil reservoir. Paper presented at the Fourth Workshop of the IEA Collaborative Project on EOR,
Vienna, Austria, August 25-26.

Khan, S.A., Pope, G A. and Sepehrnoort, K , 1992. Fluid characterization of three-phase CO,/oil
mixtures Paper SPE/DOE 24130 presented at the SPE/DOE Eighth Symposium on EOR, Tulsa,
Oklahoma, April 21-24

Ktllough, J.E. and Kossack, C A., 1987. Fifth SPE Comparative Solution ProJect: Evaluatron of
Miscible Flood Simulators Paper SPE 16000 presented at the Ninth SPE Symposmm on Reservoir
Stmulation, San Antonio, Texas, Feb 1, 55-73

Larsen et al., 1989. Temperature Dependence of Ll - L2 - V Phase Behavior m COZ-Hydrocarbon

Systems, Sot Petro Eng J., 105-114.

Leverett, M.C., 1941. Capillary Behavior m Porous Solids, Trans., AIME, 142. 159-172.

Lim, M.T., Khan, S A., Pope, G.A and Sepehrnoori, K , 1992. Stmulation of carbon dioxide
flooding using horizontal wells. Paper SPE 24929 presented at the SPE 67th Annual Techmcal
Conference and Exhtbttton, Washmgton, DC. October 4-7

Mathias, P M. and Benson, M.S , 1986. Computational Aspects of Equattons of State Fact and
Fiction. AIChE Journal, Vol. 32, No 12. 2087-2090.

Mtchelsen, M., 1982 The isothermal flash problem. Part I. Stability Fluid Phase Equilibria, 9 1-19

Nanba, T., Kmchi, T and Tanakadate, T., 1986. Multiphase equilibrium calculations for CO,
flooding. Paper presented at the Seventh Workshop of the IEA Collaborative Project on EOR,
Hannover, Germany, Sept 16-17
108 P. Wang, E.H. StenbylFluid Phase Equilibria 94 (1994) 93-108

Nelson, P.A , 1987. Rapid phase determmation m multtple-phase flash calculations, Compu Chem
Eng., 11 581-591

Ogmo, K., 1988. Compositional simulation of carbon dioxide oil recovery experiments. M.S. Thesis,
The University of Texas at Austin, May

Orr, FM., Jr., Yu, A D and Lien, C L., 1981 Phase behavior of CO, and crude 011 m low
temperature reservoirs. Sot Pet Eng. J., 21, No.4 480-491.

Peng, D.-Y and Robmson, D B., 1976. A new two-constant equation of state. I&EC Fundam, 15
59-64.

Trangenstem, J A., 1987 Customized mmimrzation techmques for phase equilibrium computations
in reservoir simulation. Chem Eng. Sci., 42: 2847-2863

Turek, E A, Metcalfe, R S , Yarborough, L. and Robmson, R L., 1984 Phase equihbria in CO,-
multicomponent hydrocarbon systems. Experimental data and an improved prediction technique Sot
Pet Eng J., 24, No 3 308-324.

Young, L.C , 1991 Full-field compositional modellmg on vector processors, SPE Res Eng , Feb
107-114.