Ind. Eng. Chem. Res. 1997, 36, 4291-4298
4291
SEPARATIONS An Efficient Initialization Procedure for Simulation and Optimization of Large Distillation Problems Pierre Rabeau and Rafiqul Gani* Department of Chemical Engineering, Technical University of Denmark, DK-2800 Lyngby, Denmark
Claude Leibovici Elf Aquitaine Production, Centre Scientifique et Technique Jean Feger, 64018 Pau Cedex, France
An efficient initialization procedure for steady-state distillation simulation problems involving a large number of components is presented. The proposed initialization procedure employs “lumping” and “delumping” of the components present in the mixture. The problem of initialization is then reduced to a system of component mass balance equations for each plate. The procedure generates good initial estimates and is not time consuming. The method initially developed for petroleum fluids has proven to be efficient also for other hydrocarbon mixtures. This is shown with two series of examples of distillation problems, by comparison with traditional methods of initialization. Two optimization problems solved with a “lumped” system are also presented and compared to the optimum solution found with the original system. Introduction Although many computational schemes have been reported in the past, active research is still being carried out in developing efficient methods of solution for simulating steady-state behavior of distillation columns. The method of solution generally requires iterative calculations, so the computational time depends strongly on the initial estimates of the iterative variables. The computational time also depends on the size of the problem. Very large problems are generated due to a large number of components in the mixture to be separated and/or a large number of stages. In the petroleum industry, systems of more than 300 components are not uncommon since the crude oil feed to the distillation column is usually a blend of several crudes, each represented by 50 or more real and hypothetical components. The computational time may become prohibitive for such problems; hence, new time-reducing solution schemes are required. Wayburn and Seader (1983) proposed a mathematical model of great generality for distillation-based (interlinked or noninterlinked) separation processes and a robust method of solution. While their procedure addressed the important issues of flexibility and robustness, it did not address the question of computational time for large problems. In this paper, we address flexibility, robustness, and computational time, with special emphasis on large distillation simulation problems. Even in design computations, which are usually off-line, reduction in computational time can be useful since alternatives can be explored. Also, reductions in computational times should be machine independent. The goal should not be to use an inefficient solution procedure on a fast computer to obtain rapid solutions. Faster computers and efficient machines should provide * Author to whom all correspondence should be addressed. S0888-5885(97)00270-4 CCC: $14.00
additional computations which would not be feasible otherwise. In petroleum engineering, the usual way to reduce the size of the problem when simulating operations with distillation is to lump the components present in the oil mixture. Montel and Gouel (1984) proposed a lumping scheme based on the similarities of a few properties of all the components. An iterative clustering algorithm around mobile centers yields a classification into “pseudo”components optimum with respect to the considered equation of state. Leibovici (1993) has established a procedure for the estimation of properties of the pseudocomponents which is accurate and completely consistent with the thermodynamic model to be used. Even though simulations can be performed with this reduced system of pseudocomponents in a reasonable time, it should be noted that the results from simulation only deal with a simplified representation of the fluid (the system of pseudocomponents) and no information on the original system is provided here. Therefore, a large knowledge gap exists between the increasing amount of analytical data provided by modern laboratory equipments and the simplified simulation results. It is obvious that employing a detailed description of the mixture in process optimization or control will lead to more accuracy and quality of the resulting solution. However, these calculations are not done presently since it is impossible for reasons of either memory storage or computational time. When a simulation of a flash operation is performed, the detailed results can be obtained by using a “delumping” procedure, i.e., generating the detailed mixture results from the “lumped” solution. Leibovici et al. (1996) found a new way to “delump” a system which is completely consistent with the thermodynamic model to be used. The procedure is based on properties of the equilibrium constants involved in the phase equilibrium calculation. The new initialization procedure for steady-state © 1997 American Chemical Society
4292 Ind. Eng. Chem. Res., Vol. 36, No. 10, 1997
Figure 1. Distillation column configuration and an arbitrary stage in the column.
one equation per stage for energy balance: distillation models proposed in this present paper employs the “lumping” and “delumping” techniques which lead to a linear system of algebraic equations. Solution of this set of equations provides the initial estimates for the actual nonlinear problem (steady-state distillation). Steady-state distillation problems for mixtures with 20 < N < 320 components have been solved by employing this new procedure. It has been observed that the ratio of computational time needed for initialization to the computational time needed for N components rigorous distillation simulation decreases as N increases. Thus, while the order of reduction of computational times for N ) 20 is low (around 2-3), the corresponding order of reduction for N ) 320 is 8. Various distillation problems, with and without side products, pumparounds, and side-strippers, have been considered. In all cases, significant reductions in computational times to reach the rigorous solution of the distillation problem have been achieved. It has been observed that the initial estimates obtained with the new initialization procedure were almost as precise as the rigorous solution of the distillation problem. The details of the new method are presented together with the solution of industrially relevant problems. The results were obtained with the commercial simulator PRO/II (1992) as well as the in-house simulator SEPSIM (Andersen, 1985).
Model Description Figure 1 represents a typical equilibrium stage of a distillation column. For each stage, equilibrium conditions between the liquid and the vapor streams leaving the stage are assumed and the model equations are given by
(Li + SiL)hi + (Vi + SiV)Hi - Li+1hi+1 Vi-1Hi-1 - Qi ) 0 (1) N equations per stage for component material balance:
(
) (
)
SiL SiV lij 1 + + vij 1 + - li+1,j - vi-1,j - fij ) 0 Li Vi
(2)
N equilibrium relationships per stage: Vi vij - Kij lij ) 0 Li
(3)
The dependent variables Li and Vi are calculated from the summation equations: N
Li )
lij ∑ j)1
(4)
N
Vi )
vij ∑ j)1
Naphtali and Sandholm (1971) chose the liquid and vapor component flowrates (lij, vij) and the temperatures (Ti) for each plate as iterative variables, thereby giving a total of 2N + 1 variables. They found that, if the equations and variables are grouped in terms of a stage, the Jacobian matrix in the Newton-Raphson method for these equations has a block tridiagonal form which greatly simplifies the solution procedure. The method of Naphtali and Sandholm also has the advantage of great flexibility because it can handle absorption and distillation problems without requiring any additional changes. Many other algorithms have been proposed. They may be classified into two categories: equation
decoupling methods and simultaneous solution methods. The θ-method (Holland, 1975), the inside-out class of methods originated by Boston and Sullivan (1974), and the method of La´ng et al. (1991) developed especially for crude oil distillation fall into the former category, whereas the method of Naphtali and Sandholm, Christiansen et al. (1979), Hofeling and Seader (1978), and Letourneau et al. (1995) when pumparounds and/or recycles are present and Buzzi Ferraris (1983) fall into the second category. One point to note is that the most general initialization that can be done for distillation problems solved with any of these methods is to provide initial values for lij, vij, and Ti and the corresponding Kij, Li, and Vi.
) [
(
) ]
Vi-1 Vi SiL SiV - Ki-1,j li-1,j + 1 + + 1+ Kij lij Li-1 Li Vi Li li+1,j - fij ) 0 (6) Since the feed (fij) and side stream flowrates (SiV, SiL) are given in the problem definition, the problem consists of solving for each component a linear system of p algebraic equations in the liquid component flowrates lij. Then, verification of the solution obtained is made. The physical constraint concerning the positiveness of the component flowrates is checked
∀i,j
lij g 0
(7)
and the values of lij are normalized with respect to the total internal flowrates (Li) to verify
Method of Solution Let us consider a mixture of N components. In the lumping procedure described by Montel and Gouel (1984), components with similar properties (Tci, Pci, mi, and Mwi) are grouped in clusters. In this way, one can decide to generate M clusters (M < N), also called pseudocomponents. The critical and other physical properties of these pseudocomponents needed for the thermodynamic description of the phase equilibrium are calculated with the procedure developed by Leibovici (1993). With the reduced list of M pseudocomponents, a simulation is performed. Then, the first assumption is made: the total flowrates and temperature profiles (Li, Vi, Ti) found with this simulation are fully representative of the original problem, i.e., Li, Vi, and Ti are supposed to be equal in the simulation of the M pseudocomponents and of the N original components systems. The second assumption is that the equilibrium constants (Kij) of the N components system can be calculated from the ones of the M pseudocomponents system. Michelsen (1986) presented an efficient flash procedure for systems described by cubic equations of state with all binary interaction coefficients equal to zero. On the basis of this simplification, Leibovici et al. (1996) showed that the equilibrium constants lie on a geometrical locus when the mixture parameters of an equation of state can be expressed as a linear combination of pure-component parameters and molar compositions. For example, for a system described by a twoparameter cubic equation of state, the locus has the following expression:
ln Ki ) ∆C0 + ∆C1xai + ∆C2bi
(
Ind. Eng. Chem. Res., Vol. 36, No. 10, 1997 4293
(5)
where ai and bi are the equation of state purecomponent parameters. Therefore, knowing the ∆C coefficients, the K-factors for the original N components mixture can be determined by extrapolation and/or interpolation and use of the temperature profile in the evaluation of ai. This procedure, which is detailed in the appendix, was used for the second assumption. At this stage of the procedure, one knows the Kfactors and the internal and external total flowrates for the original N components mixture without any timeconsuming calculations. If the K-factors are introduced in the mass balance equations (2) through their definition in eq 3, the following equation is obtained:
N
Li )
lij ∑ j)1
(8)
The vapor component flowrates vij are also available directly from the equilibrium relationships (eq 3), and the same verification as for lij is made. Now, the whole set of variables is initialized for the rigorous solution of the distillation problem. Computational Procedure A step-by-step solution technique has been defined (Chart 1); it is given below and illustrated in Figure 2. Chart 1. New Method of Initialization Step 1. For the oil blend represented by N components, use a lumping procedure to generate M clusters (n + 1 < M , N), each representing a subset of components and possessing a full set of properties. Step 2. Perform a steady-state distillation simulation with the M clusters and save the results for Li, Vi, Ti, and Kij. Step 3. From the Kij for the M clusters and Ti obtained in step 2, calculate the Kij for the N components of the original mixture (Leibovici et al., 1996). Step 4. From the Kij calculated in step 3 and the Li and Vi obtained in step 2, calculate lij by solving component by component the system of linear equations (6). Step 5. From lij, Vi, Li, and Kij obtained in steps 2-4, calculate vij from eq 3.
Example Problems In order to try the initialization procedure, several test examples were used. Results from the ones which best represent all the results and highlight the main features of the initialization method are presented. Although the method was initially developed for applications involving petroleum mixtures, petrochemical mixtures have also been tested. The differences between these two types of mixtures are the number of components present in the mixture and the continuity of the mixture properties. A petrochemical mixture is considered as a “discrete mixture” due to the high discontinuity of its set of pure-component properties, whereas a petroleum mixture is considered to be a “continuous mixture”. Example 1 deals with a petrochemical mixture, and examples 2 and 3 deal with petroleum mixtures. Then, two optimization problems such as cost minimizations are presented, and the optimum solutions obtained with lumped and original systems are compared and discussed.
4294 Ind. Eng. Chem. Res., Vol. 36, No. 10, 1997 Chart 2. Clustering Algorithm for Distillation Step 1. From a shortcut distillation simulation (Smith, 1963) or other existing data, determine the recovery ratio for each component (Ri ) di/fi). Step 2. From the values of Ri found in step 1, generate four classes as follows:
1.00 g Ri g 0.98 w i ∈ C1 0.98 > Ri g 0.50 w i ∈ C2 0.50 > Ri > 0.02 w i ∈ C3 0.02 g Ri g 0.00 w i ∈ C4
Figure 2. Scheme of the new approach of initialization. Table 1. Output Results for the Dominant Components in Example 1 fluid rates, kmol/h feed liquid
distillate vapor
bottom liquid
n-pentane 2-methyl-2-butene cyclopentene cyclopentane benzene toluene ethylbenzene dicyclopentadiene others
29.58 20.45 28.88 23.78 203.69 64.78 39.39 29.61 103.32
29.58 20.45 28.87 23.74 0.72 0.00 0.00 0.00 62.11
0.00 0.00 0.01 0.04 202.97 64.78 39.39 29.61 41.20
total rates, kmol/h temperature, K pressure, atm
543.48 314.45 1.0
165.48 345.86 2.5
378.00 405.53 2.8
Petrochemical Mixture. Data from a petrochemical plant operating in the region of Normandy (France) are used in example 1. A 34-stage depentanizer column operating at a pressure of 2.5 atm and having one feed entering the column on stage 17 and two products is considered (stages numbered from bottom to top of the column). From chromatographic analysis, a 79-component mixture is considered as the feed of the column. Eight components of the feed are in higher amounts; they represent 80.44% (molar basis) of the mixture. Table 1 gives the steady-state condition in terms of pressure, temperature, and the component flowrates in the feed, in the distillate product, and in the bottom product for the eight dominant components. The output results for a 79-component mixture are available as supporting information. The location of the clusters was imposed instead of using the “automatic” dynamic clustering method of Montel and Gouel (1984) because numerical experiences suggested special caution when applying this method to distillation column simulation for “discrete” mixtures. Since the mixture composition changes along the length of the column, lumping two key components together, for example, will result in incorrect predictions. Thus, from various experiments and engineering insights, a step by step procedure for choosing the location of the clusters (Chart 2) was defined, used, and is recommended. Following this procedure, 10 clusters were made. The properties for each cluster were predicted by the method of Leibovici (1993). Simulations were performed with the rigorous distillation algorithm of SEPSIM, which is a simultaneous solution method. With the conventional initial estimates generated in this program, the computational time (IBM RS/6000) needed to reach the solution was 47.7 s (4 iterations), whereas with the new initialization
Step 3. Identify dominant components in the feed of the column and assign them as centers for new classes (included in C1, C2, C3, C4). Step 4. Assign each component to the “nearest” center according to a defined distance (the distance formulation given by Montel and Gouel (1984) could be used for example).
method, 23.3 s (2 iterations) was needed (note that only 2 s is used for the generation of the initial estimates). The time reduction was then a factor of 2, which is low because the problem was not difficult to converge and the conventional initial estimates were already close to the solution. Particular attention was given here to the quality of the initial estimates generated with the new method or the deviation between these values and the actual solution. The specifications for this simulation problem were a mass percentage of benzene equal to 0.5% in the distillate and a reflux ratio equal to 3. The bottom product flowrate found with the new initialization method was 376.5 kmol/h, and it was 378 kmol/h in the rigorous solution. The values of Li and Vi corresponding to simulations with the lumped system and the original components are available as supporting information. The average deviation for Li and Vi between the two simulations is 1.29 kmol/h, and the maximum deviation is 2.48 kmol/h. In Figure 3, the temperature profiles corresponding to these two simulations are plotted, the average deviation is 0.38 K, and the maximum deviation is 0.72 K. These comparisons suggest that the calculations could be stopped after the initialization proposed in this paper; i.e., consider the initial estimates as the actual solution. A significant time reduction will then be achieved; for this example the reduction is the ratio of the computational times for a 79-component simulation and for a 10-component simulation (factor of 24 in this case). This suggestion requires that the guidelines of Chart 2 are automatized in a way which guarantees a correct representation of the mixture in the range of temperature and compositions considered in the distillation problem. Petroleum Mixture. Moderately Large Problems. A set of three cases involving columns with different configurations were treated in example 2. In each case, a 20-tray column was considered; 132 components were involved in examples 2a and 2c and 172 components in example 2b; no side products were withdrawn in examples 2a and 2b, whereas two side products were withdrawn in example 2c; in each example 39 clusters were considered in the lumping procedure (Montel and Gouel, 1984). Table 2 gives the number of iterations and computational times for each example. These results are presented for four simulation modes: SIM A, the “normal” distillation simulation (conventional initialization) with N components; SIM B, simulation with M clusters; SIM C, simulation with N components and initial estimates (coming from SIM B) for Li, Vi, and Ti; SIM D, simulation with N components and initial
Ind. Eng. Chem. Res., Vol. 36, No. 10, 1997 4295
Figure 3. Temperature profiles for example 1. Table 2. Computational Times and Number of Iterations for Examples 2a-c simulation mode SIM A simulator
iter.
timea
SIM B
SIM C
SIM D
iter. time iter. time iter. time
SEPSIM CMO SEPSIM RIG
9 6
Example 2a 121 10 4 372 7 14
SEPSIM CMO SEPSIM RIG
11 6
Example 2b 638 10 4 2042 6 13
NA 3
NA 727
3 2
180 471
SEPSIM CMO SEPSIM RIG
12 5
Example 2c 167 12 5 382 6 13
NA 2
NA 90
4 2
59 90
a
NA 3
NA 127
3 2
44 85 Figure 4. Crude oil unit of example 3.
Simulations were run under UNIX system on an IBM RS/6000.
estimates for lij, vij, Ti, Li, and Vi generated with the new initialization method. With the rigorous distillation simulation option of SEPSIM (referred to here as SEPSIM RIG), the computational time reduction by applying the new method of initialization varied from 3.7 to 4.2. With the “constant molar overflow” version of the distillation algorithm of SEPSIM (SEPSIM CMO), the computational time reduction varied from 2.5 to 3.5. The same problems were also solved with PRO/II (1992) (using SURE algorithm), and good improvements (between factors of 2 and 4) were also obtained. Very Large Problem Example. A special version of the “inside-out” column (Leibovici, 1996) able to use the complete set of initial guesses for all iterative variables has been used for the simulation of example 3 dealing with a mixture of five crude oils represented by a total of 327 components (7 pure components + 320 petroleum cuts). As shown in Figure 4, a typical crude unit operating at atmospheric pressure, having the crude oil feed and a bottom feed of steam and involving three side strippers (each with a bottom feed of steam) and three recycling pumps, was considered. In the lumping procedure, the seven pure components were unchanged and the 320 petroleum cuts were lumped into 30 pseudocomponents using the dynamic clustering method of Montel and Gouel (1984). Decant water was treated as a third pure liquid phase in the calculations. By applying the new initialization procedure, the total computational time to reach the rigorous distillation solution was reduced by a factor of 8 for this problem. The computational times for lumping and delumping as well as for the distillation simulation of the lumped system of 37 components were negligible compared to the other computational times.
For all the examples of this Petroleum Mixture section, comparisons of the variables (lij, vij, Ti) between the actual solution of the problem and the initial estimates calculated with the new method showed an average relative deviation of between 0.5 and 1.5% and a maximum relative deviation between 3 and 5%. Therefore, the option to take the initial estimates as the actual solution of the problem appears again to be a promising alternative. Their use in optimization/control studies is highly recommended by this paper, since these results (estimates of lij, vij, Ti, Li, Vi and Kij) are obtained without much computational effort (in example 2c for instance, 13 s is necessary to calculate the initial estimates which differ by 0.5% from the rigorous solution which was found in 382 s; the corresponding time reduction is a factor of 29.38). This example shows that the initialization procedure could also be used as a first estimate in design-related computations. It should be noted, however, that for detailed design a final rigorous simulation would be needed. Since “good” initial estimates would be available, this step would be performed at significantly reduced computational time. Finally, we have seen that the initialization method is valid and efficient for both equation decoupling methods and simultaneous solution methods. Optimization. The solutions of two cost minimization problems are reported in this section. The emphasis is given on the computational aspects of the problem and not the actual operating cost. In both problems, a petrochemical mixture represented by 20 hydrocarbons (from C1 to C10) is considered. A 25-stage debutanizer column operating at a pressure of 15 atm and having one feed entering the column on stage 20 is considered (stages numbered from bottom to top of the column). Table 3 gives the steady-state condition in terms of pressure, temperature, enthalpy, molecular weight, and the component flowrates in the feed, in the distillate product and in the bottom product at the starting point (base case). The choice of the location of the clusters was done “manually” by taking into account engineering insight
4296 Ind. Eng. Chem. Res., Vol. 36, No. 10, 1997 Table 3. Output Results for the 20-Component Mixture Distillation (Base Case)
Table 6. Results for Optimization Problems 1 and 2 (a) Problem 1
fluid rates, kg‚mol/h feed liquid
distillate vapor
bottom liquid
5.000 10.000 40.000 50.000 260.000 250.000 140.000 100.000 100.000 80.000 30.000 20.000 5.000 20.000 30.000 30.000 20.000 30.000 20.000 20.000
5.00 9.999 39.999 49.863 238.257 241.170 8.018 27.624 0.000 0.000 0.000 0.000 0.000 0.004 0.000 0.065 0.000 0.000 0.000 0.000
0.000 0.001 0.001 0.137 21.743 8.830 131.982 72.376 100.000 80.000 30.000 20.000 5.000 19.996 30.000 29.935 20.000 30.000 20.000 20.000
1260.000 375.00 15.00 21.49 70.28
620.000 361.64 15.00 16.67 55.07
640.000 447.42 15.00 23.79 85.01
1. methane 2. ethane 3. ethylene 4. propane 5. butane 6. isobutane 7. pentane 8. isopentane 9. hexane 10. heptane 11. octane 12. nonane 13. decane 14. 2,2-dimethylbutane 15. 3-methylpentane 16. cyclopentane 17. cyclohexane 18. methylcyclopentane 19. cycloheptane 20. methylcyclohexane total rate, kg‚mol/h temperature, K pressure, atm enthalpy, M‚kJ/h molecular weight
Table 4. Components Present in Each Cluster for Systems 1 and 2 clusters
components present in the cluster
system 1
system 2
system 1
system 2
1 2 3 4
1 2 3 4 5 6
1-4 5 and 6 7 and 8 9-20
1-4 5 6 7 8 9-20
Table 5. Details of the Cost Minimization Problems problem
1
2
f ) |Qreb| + |Qcond| 1 < FTRA < 25 350 < T(feed) < 450 (1) recovery of iC4 + nC4 in the bottom ) 2% (2) recovery of iC5 + nC5 in the top ) 2%
objective function f ) Qreb design variables 1 < FTRA < 25 specifications of the column
and the guidelines given in Chart 2. Two lumped systems (systems 1 and 2) were considered because the specifications of the debutanizer given in the calculations were directly related to the key components. The compositions of the clusters of systems 1 and 2 are shown in Table 4. For this problem, the initialization method is applied in the following way: the optimization problem is solved with the lumped system and then the optimal solution is delumped, producing the required initialization for the final rigorous distillation. In the first optimization problem, the optimum feed tray location is investigated in order to minimize the operating costs; the objective function solely depends on the reboiler heat duty. In the second optimization problem, the optimum feed tray location and temperature are investigated in order to minimize the operating costs; the objective function depends on the reboiler and condenser heat duty. The details of the problems are shown in Table 5 and the results in Table 6. These problems were treated with the simulator PRO/II (CHEMDIST algorithm for the distillation columns simulations).
mixture
computational time, s
optimum feed tray locationa
system 1 system 2 original
14.66 20.76 135.00
12 12 12
objective function, M‚kJ/h starting final gain 61.91 66.43 67.13
36.86 39.46 39.29
25.05 26.97 27.84
(b) Problem 2 objective function, compu- optimum optimum M‚kJ/h tational tray feed mixture time, s location temperature starting final gain system 1 22.08 system 2 25.32 original 161.00 a
12 13 12
399.92 402.44 399.40
104.02 51.84 52.18 113.14 57.34 55.80 115.27 57.99 57.28
Rounded numbers.
In both problems considered here, small differences were noticed between the results obtained with the original and the lumped systems. The results are nearly fully reproduced when working with the lumped systems. This shows the benefit of also using the initialization procedure in optimization. Although the two examples chosen are simple, much more complex configurations involving separation processes with interlinked columns can be considered. In that case, several lumped systems at different points of the flowsheet may be necessary. This should not introduce any new computationally expensive steps to the method of solution. Conclusions A rapid and safe procedure for initialization of steadystate distillation calculations has been developed. The procedure is particularly attractive for cases where the mixture properties can be considered to be “continuous”. The computational times needed to converge the rigorous distillation model for the full system can be reduced by up to a factor of 8 through the initialization procedure. Note that the corresponding reduction of the absolute computational time (which is very high for these problems) is of greatest interest. Furthermore, since the estimated values obtained with the initialization procedure are very close to the actual solution, they can be used, in many cases, directly as the actual solution, thereby leading to much higher time reductions. It has been shown that the method is also applicable to smaller discrete mixtures, such as the ones found in the petrochemical industry, but with special care (engineering insights related to distillation) when choosing the clusters. The computational time reduction can be very profitable for optimization studies. The results show the future potential of the initialization method in terms of design, simulation, and optimization of complex flowsheets involving several distillation columns and recycles. It should be noted that, even though computation examples for simulation involves only a simple column, the advantage of using the proposed initialization procedure should be obvious when such columns are part of a flowsheet with recycle loops. Also, the method is independent of computer hardware. That is, there will be significant savings in all computers. Also, the examples have shown that the results can be used for initial design and as initial estimates for final design. Current and future work will
Ind. Eng. Chem. Res., Vol. 36, No. 10, 1997 4297
address the computational aspects related to optimization of flowsheets with distillation trains and recycles.
n
ln Φi ) C0 +
∑ Ckpi,k
(A1)
k)1
Acknowledgment The authors gratefully acknowledge Elf Aquitaine Production for providing financial support for this project and the Engineering Research Center: Phase equilibria and Separation Processes (IVC-SEP) of the Chemical Engineering Department of the Technical University of Denmark. We are also thankful to Elf Atochem for providing petrochemical plant data.
provided that the mixture parameters [or any of their power transforms] can be expressed as linear combinations of the pure-component parameters, the weighting factors being the component mole fractions (∀k, pk ) ∑xipi,k). For example, with the two-parameter SRK equation of state (Soave, 1972) the values of pi,1 and pi,2 are (ai)1/2 and bi when kij ) 0 (conditions for the linear combination). At equilibrium the isofugacity criterion
fiL ) fiV
Nomenclature a, b ) mixture parameters in eq A4 ai, bi ) pure-component parameters in the equation of state B ) bottom products flowrate C1, C2, C3, C4 ) first four classes used for proposed clustering di ) flowrate of component i in distillate f ) objective function in Table 5 fij ) molar feed flowrate of component j on plate i Fi ) total molar feed flowrate to plate i FTRA ) feed tray location hi ) liquid molar enthalpy on plate i Hi ) vapor molar enthalpy on plate i Kij ) equilibrium constant of component j on plate i lij ) molar liquid flow of component j on plate i Li ) total liquid molar flowrate on plate i M ) number of pseudocomponents Mwi ) molecular weight of component i N ) number of components n ) number of independent parameters in the equation of state P ) number of plates Pi ) pressure on plate i Qi ) heat flowrate to stage i Qreb ) reboiler heat duty Qcond ) condenser heat duty R ) gas constant in eq A4 R ) reflux ratio Ri ) recovery ratio Si ) side products flow from plate i Ti ) temperature on plate i Tci, Pci ) critical properties of component i v ) molar volume in eq A4 vij ) molar vapor flow of component j on plate i Vi ) total vapor molar flowrate on plate i xiF ) molar composition of component i in the feed Z ) compressibility factor in eq A4 Greek Letters ωi ) acentric factor of component i ∆C0, ∆C1, ∆C2 ) regression coefficients in eq 5 Superscripts L ) liquid phase V ) vapor phase
Appendix Leibovici et al. (1996) developed a method to predict the detailed composition of the phases resulting from a flash calculation just performed on a lumped mixture. If one considers any equation of state P ) f(T,V) involving n parameters, the fugacity coefficient of component i in a mixture of N components can be expressed simply as
(A2)
has to be fulfilled for all components, and the logarithm of equilibrium constant can therefore be written as n
ln(Ki) ) ∆C0 +
∑ ∆Ckpi,k
(A3)
k)1
where ∆Ck ) CkL - CkV. This relation shows that, in the case where no binary interaction coefficients (kij) are used to represent the mixture, the coefficients ∆Ck are known analytically and the K-values of all components which were lumped together can be easily obtained from eq A3 assuming that the phase parameters have not changed. Binary interaction coefficients kij are small, typically zero for hydrocarbon-hydrocarbon interactions and in the range 0-0.15 otherwise. If they are used to describe the lumped system, then the coefficients ∆Ck can be approximated by any standard least-squares fit method and requested K-values can be predicted from eq A3. Note that in order to obtain the (n + 1)∆Ck values the number of pseudocomponents in the lumped system has to be at least n + 1. For the SRK equation of state with no binary interaction coefficients the analytical expressions for ∆C0, ∆C1, and ∆C2 are
(
∆C0 ) ln
∆C1 )
∆C2 )
(
)
vV - bV vL - bL
) ( ) ( ) ( )
2xaV bV 2xaL bL ln 1 + V - L ln 1 + L V b RT v b RT v
ZL - 1 aL bL ZV - 1 + 1 + ln bL (bL)2RT vL bV bV aV ln 1 + V (A4) V 2 (b ) RT v
In practice the procedure contains the following steps: Step 1. Lump a mixture of N components into M pseudocomponents. Step 2. Perform a flash calculation with the lumped system in order to obtain the equilibrium constants of the pseudocomponents. Step 3. Determine the ∆Ck values from the lumped system. In the case of all kij ) 0, use eq A4. In the case of some kij * 0, use regression (note that the condition n + 1 e M e N has to be fulfilled). Step 4. Calculate Ki from eq A3 for all components in the full system.
4298 Ind. Eng. Chem. Res., Vol. 36, No. 10, 1997
This procedure will work even better if the lumping does not affect, or only slightly affects, the phase parameters and the phase split. If this is the case, then the ∆Ck are the same for the lumped scheme and the original one. Supporting Information Available: Tables of output results and internal flowrates for 79-component mixture distillation (3 pages). Ordering information is given on any current masthead page. Literature Cited Andersen, P. M. Steady-State Simulation of Chemical Processes. Ph.D. Thesis, Technical University of Denmark, Lyngby, Denmark, 1985. Boston, J. F.; Sullivan, S. L. A New Class of Solution Methods for Multicomponent Multistage Separation Processes. Can. J. Chem. Eng. 1974, 52, 52-63. Buzzi Ferraris, G. A Powerful Improvement of the Global NewtonRaphson Method for Multistaged Multicomponent Separators. Comput. Chem. Eng. 1983, 7, 73-85. Christiansen, L. J.; Michelsen, M. L.; Fredenslund, Aa. NaphtaliSandholm Distillation Calculations for NGL Mixtures Near the Critical Region. Comput. Chem. Eng. 1979, 3, 535-542. Hofeling, B. S.; Seader, J. D. A Modified Naphtali-Sandholm Method for General Systems of Interlinked Multistaged Separators. AIChE J. 1978, 24, 1131-1134. Holland, C. D. Fundamentals and Modelling of Separation Processes; Prentice-Hill: Englewood Cliffs, NJ, 1975. La´ng, P.; Szalma´s, G.; Chika´ny, G.; Keme´ny, S. Modelling of a Crude Distillation Column. Comput. Chem. Eng. 1991, 15, 133-139. Leibovici, C. A Consistent Procedure for the Estimation of Properties Associated to Lumped Systems. Fluid Phase Equilib. 1993, 87, 189-197.
Leibovici, C. (Elf Aquitaine Production). Personal communication, 1996. Leibovici, C.; Stenby, E.; Knudsen, K. A Consistent Procedure for Pseudo-Component Delumping. Fluid Phase Equilib. 1996, 117, 225-232. Letourneau, J. J.; Joulia, X.; Koehret, B. A Heuristic Arrangement Method for the Simulation of Interlinked Separation Columns. Proc. PSE ’94 1995, 127-132. Michelsen, M. Simplified Flash Calculations for Cubic Equations of State. Ind. Eng. Chem. Process Des. Dev. 1986, 25, 184-188. Montel, F.; Gouel, P. L. A New Lumping Scheme of Analytical Data for Compositional Studies. Presented at the 59th Annual Technical Conference and Exhibition, Houston, Sept 16-19, 1984; Paper SPE 13119. Naphtali, L. M.; Sandholm, D. P. Multicomponent Separation Calculations by Linearization. AIChE J. 1971, 17, 148-153. PRO/II manual, Version 3.3, Oct 1992 (Simulation Sciences Inc., USA). Smith, B. D. Design of Equilibrium Stage Processes; McGraw-Hill Book Co.: New York, 1963. Soave, G. Equilibrium constants from a modified Redlich-Kwong equation of state. Chem. Eng. Sci. 1972, 27, 1197-1203. Wayburn, T. L.; Seader, J. D. Solution of Systems of Interlinked Distillation Columns by Differential Homotopy-Continuation Methods. Proc. Second Int. Conf. Found. Comput.-Aided Process Des. 1983, 765-862.
Received for review April 10, 1997 Revised manuscript received July 25, 1997 Accepted July 26, 1997X IE970270Z
X Abstract published in Advance ACS Abstracts, September 15, 1997.
