An Algorithm For Calculating Vapour-Liquid Equilibrium
Authors: Kao J. 
Published: December 2007

Abstract
An algorithm for calculating vapour-liquid equilibrium in multi-component systems is presented. The algorithm uses the Newton-Raphson method to solve the non-linear algebraic equations defining equilibrium. Iterative sequences have been organised to trace the progress of a calculation to ensure an efficient and progress of a calculation to ensure an efficient and successful calculation. Initialisation of the calculation procedure is automatic and relies only on parameters procedure is automatic and relies only on parameters characterising the components in a system. Important features discussed in this paper include a generalised flash calculation defined as a calculation of equilibrium conditions specified by a feed and any valid combination of two thermodynamic variables, and the incorporation of any suitable equations of state for evaluating liquid and vapour properties. An example problem is studied to illustrate the efficiency of the algorithm and a comparison is made between the Newton-Raphson method and the successive substitution method which is commonly used in phase-equilibrium calculations.

Introduction
Phase-equilibrium calculations are often required in the solution of problems involving simultaneous flow of liquid and vapour. Retrograde condensation in gas-condensate reservoirs and gas liberation in volatile oil reservoirs are examples of two-phase flow problems currently attracting interest. In simulating such flows, knowledge of the mutual influences between the flow of the reservoir fluid and the thermodynamics of equilibrium is essential. Some insight into the complex nature of these influences can be gained using compositional reservoir simulators. The basic function of these simulators is to determine numerical solutions of the partial differential equations governing the conservation partial differential equations governing the conservation of mass, momentum etc of each component in a system coupled with the algebraic equations defining vapour-liquid equilibrium throughout the system. During the calculation procedure, a solution of the equations defining equilibrium is required at every time wad space coordinate to evaluate, for example, local flow coefficients etc. The solution of these equations can be either computed internally in a simulator or supplied as multi-dimensional tables, for example, equilibrium ratios as functions of pressure, temperature and composition. In Both cases, a routine is required to determine vapour-liquid equilibrium, hereinafter abbreviated to VLE.
The equations defining VLE are non-linear and in general can only be solved by iterative methods. Most commonly used is the method of successive substitution which is based on that proposed by Prausnitz and Chueh. However, this method exhibits slow or no convergence at conditions where the gradients of the thermodynamic variables are large, for example, in the vicinity of dew-, bubble- and critical point. To overcome convergence problems, Fussell and Yanosik developed an iterative sequence based on the Newton-Raphson method. A close examination of their paper revealed that their analysis can be extended. This prompted further investigations, the results of which are described in this paper.
The contents of this paper are arranged as follows. The equations defining VLE are first presented. On inspection, it became apparent that a flash calculation can be generalised to a VLE calculation specified by a feed and any valid combination of two thermodynamic variables, An algorithm based on the Newton-Raphson method for generalised flash calculations is then described, followed by a description of an automatic starring procedure which relies only on parameters characterizing the components in a system. A system studied by Yarborough was selected as a basis on which the efficiency of the algorithm was assessed, and a comparison between the Newton-Raphson method and the successive substitution method was made.
Phase-equilibrium calculations require equations of state to evaluate relevant properties pertaining to liquid and vapour. With several equations of state, for example ref. 4–7, available to an investigator, an algorithm which incorporates any suitable equations of state for liquid and vapour is desirable. The analysis has therefore been kept general throughout this investigation to allow for the insertion of any suitable equations of state. The experimental results of Yarborough were also used to assess the numerical results computed from the equation of Soave and that of Peng and Robinson.