Ind. Eng. Chem. Res. 1993, 32, 1528-1530
1528
RESEARCH NOTES Polynomial Objective Functions for Flash Calculations: Binary, Ternary, and Quaternary Systems John H. Warrent and Michael A. Adewumi' The Pennsylvania State University, 202 Mineral Science Building, University Park, Pennsylvania 16802
Polynomial objective functions are developed for the vapor fraction in vapor-liquid equilibrium flash calculations for binary, ternary, and quaternary systems. The form of the polynomial suggests possible expansion to N-component systems. An extensive testing of these objective functions was performed using experimental data for binary systems. Comparative analysis with traditional Rachford-Rice objective function demonstrates better predictive capability and improved computational strategy for the new form. Comparison between predictions and experimental data shows a good match, even for cases where the Rachford-Rice objective function fails to yield realistic results. Introduction The use of vapor-liquid equilibrium (VLE) flash calculations is prevalent in many types of process simulations across a wide spectrum of industries, including chemical, petrochemical, petroleum, gas, etc. Primarily, flash calculations constitute the predictive method used in establishing the phase fractions and compositions of the two contiguous phases, vapor and liquid, in a state of equilibrium. The basic problem formulation is based on simple material balance and the equilibrium constraints. From this formulation evolves an objective function, usually a nonlinearly implicit function of some phase fraction. The latter has been defined in a number of ways in the literature but the most widely-used definition is the ratio of vapor to the total system (both in terms of moles). Invariably the implicit nonlinear equation is solved for the phase fraction using an iterative technique. One of the standard solution methods is the Newton-Raphson technique or a quasi-Newton method. Several other solution techniques have also been suggested in the literature. It is a widely-known fact that the major cost in many process simulations is incurred from repetitive performance of flash calculations. This is true for compositional reservoir simulation and compositional hydrodynamic modeling of gas/gas condensate pipelines, as well as many other process simulations (Coats, 1980; Mucharam and Adewumi, 1990). Any substantial improvement resulting from either the formulation of the objective function or the formulation of the solution algorithm can substantially reduce computational costs. The other motivation arises from the need to improve VLE predictive performance for high-pressure systems, particularly in the vicinity of the critical point and more specifically in the retrograde condensation region. It is precisely for these reasons that many attempts have been reported in the literature at improving solution algorithms. Some attention has been devoted to possible reformulation of the objective function, but to a lesser degree than to new computational algorithms. The focus of this work is on the objective function itself.
* Author to whom correspondence should be addressed. +
Present address: U.S.Army, Fort Lee, VA.
T h e Rachford-Rice Objective Function The Rachford-Rice objective function (Rachford and Rice, 1952) has emerged, almost since its publication, as the standard one used in isothermal flash calculations. N
f ( 4=
z i ( K i - 1)
21+
a(Ki - 1)
=O
(1)
where a is the vapor fraction, Z i and Ki are the overall mole fraction and equilibrium ratio of the ith component in the mixture, and N is the number of components in the mixture. Other forms of the objective function have alsoappeared [e.g., Holland and Davison (195711, but none has gained as much popularity as the Rachford-Rice objective function. A few other attempts have also been reported at developing appropriate bounds within which the search for the true vapor fraction should be made to prevent convergenceproblems or even to improve convergencerate. The latest attempt along this line is the work of Leibovici and Neoscill(1992). In addition to providing bounds, they also proposed the use of a modified Rachford-Rice objective function. This modification is made by imposing a multiplier function, which is dependent on the vapor fraction bounds. Although, unlike the Rachford-Rice objective function, this function cannot be derived solely based on the physics of the problem (Le., via material balances and equilibrium constraints), their example suggestssome improvement in convergencerate. Although their work pertains to a mathematical transformation of the Rachford-Rice function, it does not preserve the governing physics of the problem. The focus of this work is to reformulate the Rachford-Rice objective function in a polynomial form, whose behavior is well-known from the long-time study of polynomials by mathematicians. Perhaps because of the familiarity with the behavior of polynomials in general, our limited testing has found this new form, which preserves the formulation integrity of the Rachford-Rice function, to be better in terms of convergence rate. I t also helps to eliminate ambiguities in functional behavior, since the behavior of polynomial functions is predictable.
0888-5885/93/2632-1528$04.00/00 1993 American Chemical Society
Ind. Eng. Chem. Res., Vol. 32, No. 7, 1993 1629
The Polynomial Form of the Objective Function In developing the polynomial objective function, the Rachford-Rice equation (eq 1)is used as the starting point. Thus, all the principles and assumptions associated with that equation apply to our development with no additional assumptions. Basically the summation term of eq 1 is expanded termwise for binary, ternary, and quaternary systems. Algebraic manipulations are then undertaken to reduce the polynomial to ita simplest form. The algebra involved for the binary system is very simple, but beyond that, it is laborious. We present the detailed algebra for only the binary system and the final equations for ternary and quaternary systems. Binary Systems. Starting with the Rachford-Rice objective function and setting Ci = Ki - 1,we rewrite eq laS
N f(a) =
zici
-- arci- 0
&1+
For a binary system, setting N = 2, eq 2 expands to f(a) =
ZlCl ; z2c2 - o 1+aC1 1+aC,
Multiplying through by (1 rearranging yields
+
aC1)(1
while more than one of them may be real. Only one of these roots has physical significance, and that is the admissible one.
Heuristic Validation of the Polynomial Functions Analysis was performed to test whether or not the polynomial function for quaternary systems, the largest mixture (in terms of number of components) studied, whould reduce to the corresponding equations for ternary and binary systems. The analyses are described below. Reduction of Ternary System. Invoking the constraint that the sum of all mole fractions is unity, we write
czi N
=1
1=1
For a quaternary system, N = 4. We proceed to reduce this to a ternary system by setting 24 = x4 = y4 = 0. At first glance, this would seem to present a problem in assigning a value to K4. This, however, is not a problem, once L”6pital‘s rule is invoked. Recalling that K4 = y4/ x4, we write
(3)
+
aC2) and
(4) or (5)
One should observe that, for binary systems, one can express a as an explicit function of overall composition and K-values. Similar algebraic exercises lead to polynomial expressions in a for ternary and quaternary systems. Omitting the detailed algebra involved, which would take several pages, the final forms of these equations are reported below. Ternary Systems. 3
Thus, for the component not present, we can assign the value of its equilibrium ratio as unity. With that, K4 = 1 and C4 = 0. Multiplying eq 7 through by nf=,Ci yields a3C,C,C3C4 + a2[(1- Z1)(C,C3C,) + (1- Z,)(C,C,C,)I + a2N1- Z3)(C,C,C,) + (1- Z4)(C,C,C3)1 + a[z,C,(C, + c3 + c4) + z2c,(cl + c3 + c,)]+ (r[z3C3(cl + c2 + C4) + Z4C4(C1+ c, + C3)1 + zlc, + z2c, + z3c3+ z4c4 = 0 (10) Applying the constraints of z4 = 0 and C4 = 0 to eq 10 yields a2(C,C,C3)
+ a[z,C,(C, + C3) + Z,C,(C, + C3) + Z3C3(C1+ C,)]
+ Z1Cl + z2c, + z3c3= 0
(11)
Dividing this equation byn&,C, and writing it in the compact form yields 3
Quaternary Systems. 4
= CY3+
1-zi
+
a2zci
=O (7)
It is apparent from these equations that the degree of the polynomial objective function is always of one less than the number of components. All these equations have one advantage in common: they can all be solved analytically for a without the need for an iterative technique. Since they are polynomial functions in a,they also suggest that more than one root exists for the equation. Specifically, for an N-component mixture, (N - 1)roots exist. It should be noted that some of these roots may be complex,
Equation 12 is the same as eq 6 for ternary systems. Reduction to Binary Systems. Using similar analysis, but setting z3 = z4 = y3 = y4 = x3 = x4 = 0, we can show that K3 = K4 = 1and thus C3 = C4 = 0. Starting with eq 10, we impose these conditions to obtain aC1C, + zlc,
+ zzc,= 0
(13)
which, after dividing through by fl:.lCi yields eq 5 for binary systems. The advantage of these forms of the objective function vis-a-visthe original Rachford-Rice objective function is apparent. For binary, ternary, and quaternary systems, the resulting objective functions are first-order polynomials, (linear), second-order polynomial (quadratic), and third-order polynomial (cubic),respectively. Hence, they
1530 Ind. Eng. Chem. Res., Vol. 32, No. 7, 1993
can be solved for the vapor fraction without any iteration. The appropriate root is the one that satisfies the physics of the problem; namely, there should be only one value for the vapor fraction in the two-phase region and ita numerical value should be between 0 and 1.
Computational and Descriptive Improvements of the Polynomial Function A rather extensive testing of the polynomial objective function was undertaken. We chose to conduct the test using binary systems simply because of the ready availability of data on such systems. The well-knownand oftenused Institute of Gas Technology data (Bloomer et al., 1953)for methane-ethane mixtures was used. The original Peng-Robinson equation of state (Peng and Robinson, 1976) was used for calculating fugacities and hence equilibrium ratios. The binary interaction parameter was set to zero. The initial set of K-values was generated using Wilson’s (1969) correlation. A wide range of temperature, pressure, and composition was used in this testing. A total of 1600 flash calculations each were performed using the polynomial function and the Rachford-Rice function. Newton-Raphson technique was used in conjunction with the latter. Not surprisingly, even for the simple binary system, the computational times (on a SUN 4/280 computer) differ by a factor of 2 in favor of the polynomial objective function. More interesting is the fact that while the polynomial objective function yields appropriate solutions for all 1600 cases tested, the Rachford-Rice function experienced divergence problems in 65 cases. These cases were those where the pressure and temperature conditions are close to the phase boundaries (dew and bubble point curves). This is perhaps because of the large value of the derivative of the Rachford-Rice function near the phase boundaries. This observation was made by Rachford and Rice. When one considers the fact that a typical reservoir simulation problem (or any similar process) involves a huge number of flash calculations, typically of the order of 10million (Leiboviciand Neoschil, 1992),the savingsin computational overhead is significant. It is also interesting to note that, in gas condensate simulation in gas pipelines, the system usually operates in the region close to the dew point curve and hence is susceptible to possible divergence problems when the Rachford-Rice function is used as experienced in our testing. Conclusions The Rachford-Rice objective function commonly used for flash calculations is recast into polynomial forms for binary, ternary, and quaternary mixtures. The degree of the resulting polynomials is one less than the number of components. Hence, for binary, ternary, and quaternary systems, vapor fraction can be solved for explicitly without
iterations. This obviously results in significant computational time saving. Moreover, extensive testing of the new functions demonstrates their superiority over the Rachford-Rice form in the region close to the phase boundaries. A mathematical analysis demonstrates the validity of the four-component polynomial function by ensuring ita collapseto the functions for ternary and binary systems. The form of the functions suggest the possibility of generalization for N-component systems in which the resulting polynomial is of degree (N - 1). Further investigation continues in this respect.
Acknowledgment Partial support of this work by the US. Army and the Gas Research Institute is gratefully acknowledged.
Nomenclature Ki = equilibrium ratio of the ith component in a mixture N = number of components in a mixture zi = mole fraction of the ith component in a mixture a = vapor fraction of the system
Literature Cited Bloomer, 0. T.; Gami, D. C.; Parent, J. D. Physical-Chemical Properties of Methane-Ethane Mixture; Institute of Gas Technology Research Bulletin 22; Illinois Institute of Technology: Chicago, 1953. Boston, J. F.; Britt, H. I. A Radically Different Formulation and Solution of the Single-StageFlash Problem. Comput.Chem.Eng. 1978,2, 109. Coats, K. H. An Equation of State Compositional Model. SPE J. 1980,20, 363. Holland, C. D.; Davison,R. R. SimplifyFlash Distillation Calculations. Pet. Refin. 1957,36, 183. Leibovici, C. F.; Neoschil, J. A New Look at the Rachford-Rice Equation. Fluid Phase Equilib. 1992, 74, 303. Mucharam, L.; Adewumi, M. A. Compositional Multiphaee Hydrodynamic Modelingof Gas/Gas-Condensate Dispersed Flow in Gas Pipelines. SPE Prod. Eng. 1990,5,85. Obut, S. T.; Ertekin, T.; Geisbreicht, R. A. A Versatile Phase Equilibrium Package for CompositionalSimulation. Proceedings, SPE California Regional Meeting, Oakland, CA, April 2-4,1986; SPE: Brookfield, CT, 1986;No. SPE 15083. Peng, D.-Y.; Robinson, D. B. A New Two-Constant Equation of State. Ind. Eng. Chem. Fundam. 1976,15,59. Rachford, H. H. Jr.;Rice, J. D. Procedure for Use of Electronic Digital Computers in Calculating Flash Vaporization Hydrocarbon Equilibrium. Pet. Trans, AIME 1952,195, 327. Wilson, G. M. A Modified Redlich-Kwong Equation of State, Applications to General Physical Data Calculations. Proceedings, AIChE National Meeting, Cleveland, OH, May 4-7,1969; AIChE: New York, 1969;paper 15C.
Received for review December 2, 1992 Revised manuscript received April 22, 1993 Accepted April 30,1993
