Iterations

In vapor-liquid equilibria, it is relatively easy to start the iteration because assumption of ideal behavior (Raoult s law) provides a reasonable zeroth approximation. By contrast, there is no obvious corresponding method to start the iteration calculation for liquid-liquid equilibria. Further, when two liquid phases are present, we must calculate for each component activity coefficients in two phases since these are often strongly nonlinear functions of compositions, liquid-liquid equilibrium calculations are highly sensitive to small changes in composition. In vapor-liquid equilibria at modest pressures, this sensitivity is lower because vapor-phase fugacity coefficients are usually close to unity and only weak functions of composition. For liquid-liquid equilibria, it is therefore more difficult to construct a numerical iteration procedure that converges both rapidly and consistently. 

The algorithm employed in the estimation process linearizes the constraint equations at each iterative step at current estimates of the true values for the variables and parameters. 

The primary purpose for expressing experimental data through model equations is to obtain a representation that can be used confidently for systematic interpolations and extrapolations, especially to multicomponent systems. The confidence placed in the calculations depends on the confidence placed in the data and in the model. Therefore, the method of parameter estimation should also provide measures of reliability for the calculated results. This reliability depends on the uncertainties in the parameters, which, with the statistical method of data reduction used here, are estimated from the parameter variance-covariance matrix. This matrix is obtained as a last step in the iterative calculation of the parameters. 

Many different manipulations of these equations have been used to obtain solutions. As discussed by King (1971), many of the older approaches work in terms of V/L, which has the disadvantage of being unbounded and which, in the classical implementation, leads to poorly convergent iterative calculations. A preferable arrangement of this equation system for solution is based on the ratio V/F, which must lie between 0 and 1. If we substitute in Equation (7-1) for L from Equation (7-2) and for y from Equation (7-4), and then divide by F, we obtain... 

However, each of these forms possesses a spurious root and has other characteristics (maxima or minima) that often give rise to convergence problems with common iterative-solution techniques. 

These systems are solved by a step-limited Newton-Raphson iteration, which, because of its second-order convergence characteristic, avoids the problem of "creeping" often encountered with first-order methods (Law and Bailey, 1967) ... 

Here the superscript (r) represents the iteration number and g is the Jacobian derivative matrix whose elements are... 

The scalar t is in the range 0 to 1 and provides step-limiting, or damping, when required to obtain a convergent iteration. 

There is justification for allowing t to increase beyond 1, and in many particular applications this may be desirable. Here a more conservative approach is used to reduce the chance of unstable iterations. 

The procedure would then require calculation of (2m+2) partial derivatives per iteration, requiring 2m+2 evaluations of the thermodynamic functions per iteration. Since the computation effort is essentially proportional to the number of evaluations, this form of iteration is excessively expensive, even if it converges rapidly. Fortunately, simpler forms exist that are almost always much more efficient in application. 

It is important to stress that unnecessary thermodynamic function evaluations must be avoided in equilibrium separation calculations. Thus, for example, in an adiabatic vapor-liquid flash, no attempt should be made iteratively to correct compositions (and K s) at current estimates of T and a before proceeding with the Newton-Raphson iteration. Similarly, in liquid-liquid separations, iterations on phase compositions at the current estimate of phase ratio (a)r or at some estimate of the conjugate phase composition, are almost always counterproductive. Each thermodynamic function evaluation (set of K ) should be used to improve estimates of all variables in the system. 

In application of the Newton-Raphson iteration to these objective functions [Equations (7-23) through (7-26)], the near linear nature of the functions makes the use of step-limiting unnecessary. 

At low or moderate pressures,a Newton-Raphson iteration is not required, and the bubble and dew-point pressure iteration can be, respectively. 

Convergence is usually accomplished in 2 to 4 iterations. For example, an average of 2.6 iterations was required for 9 bubble-point-temperature calculations over the complete composition range for the azeotropic system ehtanol-ethyl acetate. Standard initial estimates were used. Figure 1 shows results for the incipient vapor-phase compositions together with the experimental data of Murti and van Winkle (1958). For this case, calculated bubble-point temperatures were never more than 0.4 K from observed values. 

Equations (7-8) and (7-9) are then used to calculate the compositions, which are normalized and used in the thermodynamic subroutines to find new equilibrium ratios,. These values are then used in the next Newton-Raphson iteration. The iterative process continues until the magnitude of the objective function 1g is less than a convergence criterion, e. If initial estimates of x, y, and a are not provided externally (for instance from previous calculations of the same separation under slightly different conditions), they are taken to be... 

In the case of the adiabatic flash, application of a two-dimensional Newton-Raphson iteration to the objective functions represented by Equations (7-13) and (7-14), with Q/F = 0, is used to provide new estimates of a and T simultaneously. The derivatives with respect to a in the Jacobian matrix are found analytically while those with respect to T are found by finite-difference approximation... 

Again, Equations (7-8) and (7-9) are then used to calculate new compositions. These compositions, normalized, and the new value for T are utilized in thermodynamic subroutine calls to find equilibrium ratios and enthalpies for use in the next iteration. 

Convergence of the iteration requires the norm of the objective vector 1g to be less than the convergence criterion, e. The initial estimates used, if not provided externally, are, in addition to Equation (7-28)... 

The convergence rate depends somewhat on the problem and on the initial estimates used. For mixtures that are not extremely wide-boiling, convergence is usually accomplished in three or four iterations,t even in the presence of relatively strong liquid-phase nonidealities. For example, cases 1 through 4 in Table 1 are typical of relatively close-boiling mixtures the latter three exhibit significant liquid-phase nonidealities. 

Case Flash Number Type Components Mole Fraction Pressure (bar) Temperature (K) Pressure (bar) Temperature V (K) F Mole Fractions Liquid Vapor No. of Iterations... 

Flash calculations for these mixtures usually require four to eight iterations. Cases 5 and 6 in Table 1 have feeds of this type, including noncondensable components in Case 6. Within the limits of the thermodynamic framework used here, no case has been encountered where FLASH has required more than 12 iterations for satisfactory convergence. 

Liquid-liquid equilibrium separation calculations are superficially similar to isothermal vapor-liquid flash calculations. They also use the objective function. Equation (7-13), in a step-limited Newton-Raphson iteration for a, which is here E/F. However, because of the very strong dependence of equilibrium ratios on phase compositions, a computation as described for isothermal flash processes can converge very slowly, especially near the plait point. (Sometimes 50 or more iterations are required. )... 

For liquid-liquid separations, the basic Newton-Raphson iteration for a is converged for equilibrium ratios (K ) determined at the previous composition estimate. (It helps, and costs very little, to converge this iteration quite tightly.) Then, using new compositions from this converged inner iteration loop, new values for equilibrium ratios are obtained. This procedure is applied directly for the first three iterations of composition. If convergence has not occurred after three iterations, the mole fractions of all components in both phases are accelerated linearly with the deviation function... 

Each iteration requires only one call of the thermodynamic liquid-liquid subroutine LILIK. The inner iteration loop requires no thermodynamic subroutine calls thus is uses extremely little computation effort. 

As the feed composition approaches a plait point, the rate of convergence of the calculation procedure is markedly reduced. Typically, 10 to 20 iterations are required, as shown in Cases 2 and 6 for ternary type-I systems. Very near a plait point, convergence can be extremely slow, requiring 50 iterations or more. ELIPS checks for these situations, terminates without a solution, and returns an error flag (ERR=7) to avoid unwarranted computational effort. This is not a significant disadvantage since liquid-liquid separations are not intentionally conducted near plait points. 

Outside the two-phase region, ELIPS yields a value of 0 for E/F on the R-phase side and 1 for E/F on the E-phase side. Con-, vergence to these values again requires about eight or fewer iterations, except near the plait-point region where convergence is somewhat slower. 

In the highly nonlinear equilibrium situations characteristic of liquid separations, the use of priori initial estimates of phase compositions that are not very close to the true compositions of these phases can lead to divergence of iterative computations or to spurious convergence upon feed composition. 

The special estimates used in BLIPS, which are essentially pure phases of the "solvent components" (98%, with 2% of the other solvent) are chosen to avoid these problems in this iterative procedure. 

From Equation (35), an iteration function can be developed in the form... 

Convergence of this iteration is influenced by initial estimates for the true mole fractions, zThe following rules have been found to lead to rapid convergence in all cases. 

These initial estimates are used in the iteration function. Equation (37), to obtain values of the 2 s that do not change significantly from one iteration to the next. These true mole fractions, with Equation (3-13), yield the desired fugacity... 

The subsequent representations are probably reliable within the range of data used (always less broad than 200° to 600°K), but they are only approximations outside that range. The functions are, however, always monotonic in temperature, to provide appropriate corrections when iterative programs choose temperature excursions outside the range of data. 

Subroutine MULLER. MULLER iteratively solves the equilibrium relations and computes the equilibrium vapor composition when organic acids are present. These compositions are used by subroutine PHIS2 to calculate fugacity coefficients by the chemical theory. 

Second card FORMAT(8F10.2), control variables for the regression. This program uses a Newton-Raphson type iteration which is susceptible to convergence problems with poor initial parameter estimates. Therefore, several features are implemented which help control oscillations, prevent divergence, and determine when convergence has been achieved. These features are controlled by the parameters on this card. The default values are the result of considerable experience and are adequate for the majority of situations. However, convergence may be enhanced in some cases with user supplied values. 

---