Ind. Eng. Chem. Res. 1988,27, 1900-1910
1900
Registry No, Methanol, 67-56-1;ethylene glycol, 107-21-1.
Conclusion Permeation rates and selectivities of pervaporation of methanol-ethylene glycol mixtures through a cellophane membrane were found to be functions of operating temperature and feed composition. The results indicated that the separation mechanism was governed by the sorptiondiffusion behavior of the permeants, the liquid-polymer interactions in the membrane, and the intrachain hydrogen bonding of ethylene glycol, the effect of liquid-liquid interactions being negligible. Above 60% methanol content in the feed, the effect of diffusion, however, was much less pronounced, and selective separation of methanol was attributed to the preferential sorption of methanol at the upstream side of the membrane. A more detailed analytical study of the liquid-polymer interactions is likely to elucidate further the exact mechanism of pervaporation of such binary mixtures.
Literature Cited Baker, T. H.; Fisher, G. T.; Roth, J. A. “Vapour-Liquid Equilibrium and Refractive Indices of the Methanol-Ethylene Glycol System”. J . Chem. Eng. Data 1964,9, 11-12. Binning, C. R.; Lee, R. J.; Jennings, J. R.; Martin, E. C. “Separation of Liquid Mixtures by Permeation”. Znd. Eng. Chem. 1961,53, 45-50. Brun, J. P.; Larchet, C.; Melet, R.; Bulvestre, G. “Modelling of the Pervaporation of Binary Mixtures Through Moderately Swelling, Non-Reacting Membranes”. J. Membrane Sci. 1985,23,257-283. Fels, M.; Huang, R. Y. M. “Theoretical Interpretation of the Effect of Mixture Composition on Separation of Liquids in Polymers”. J . Macromol. Sci. Phys. 1971,B5(1), 89-110. Greenlaw, F. W.; Shelden, R. A,; Thompson, E. V. “Dependence of Diffusive Permeation Rates on Upstream and Downstream Pressures 11. Two Component Permeant”. J. Membrane Sci. 1977, 2,333-348. Itoh, T.; Toya, H.; Ishihara, K.; Shinohara, I. “Design of Polymer Membrane with Permselectivity for Water-Ethanol Mixutre. I1 Preparation of Crosslinked Poly(methy1Acrylate) Membrane an Diethylene Triamine and its Permselectivity”. J . Appl. Polym. Sci. 1985,30, 179-188. Larchet, C.; Brun, J. P.; Guillou, M. “Separation of Benzene-nHeptane Mixture by Pervaporation with Elastomeric Membrane. Performance of Membranes”. J. Membrane Sci. 1983,15,81-96. Mark, H. F., Gaylord, N. G., Bikales, N. M., Eds. “Cellulose”. In Encyclopedia of Polymer Science and Technology; Interscience: New York, 1965;Vol. 111. Rautenbach, R.; Albrecht, R. “Separation of Organic Binary Mixtures by Pervaporation”. J. Membrane Sci. 1980, 7, 203-223. Rogers, C. E.; Fels, M.; Li, N. N. “Separation by Permeation through Polymeric Membranes”. In Recent Developments in Separation Science; Li, N. N., Ed.; CRC: Cleveland, OH, 1976;Vol. 11. Suzuki, F.; Onozato K. “Pervaporation of CH30H-H20 Mixture by Poly(methy1 L-glutamate) Membrane and Synergetic Effect of their Mixture on Diffusion Rate”. J. App. Polym. Sci. 1983,28, 1949-1956. Suzuki, F.; Onozato, K.; Takahashi, N. “Pervaporation of Athermal Mixture of Benzene-Toluene by Poly(ethy1ene terephthalate) Membrane and Synergetic Effect on Concentration Dependence of Diffusion Rate”. J . Appl. Polym. Sci. 1982,27, 2179-2188. Suzuki, F.; Onozato, K.; Yaegashi, H.; Masuko, T. “Pervaporation of Organic Solvents by Poly[bis( 2,2,2-trifluoroethoxy)phosphazene] Membrane”. J . Appl. Polym. Sci. 1987,34,2197-2204.
Acknowledgment We acknowledge the assistance extended to us by the Advanced Cryogencis Centre of Research, Calcutta, for kindly providing us with liquid nitrogen. We are also grateful to M/S. Kesoram Rayon Ltd., India, for supplying the cellophane films for our studies. Nomenclature
J, = total permeation rate, kg/(h m2) Ji = individual permeation rate, g/(cm2 s) J, = permeation ratio T = absolute temperature, K E, = apparent activation energy, kcal/mol Do= diffusivity at zero concentration, cm2/s D = integral diffusivity, cm2/s C, = g sorbed/cm3 dry membrane L = thickness of membrane, cm W , = weight of dry membrane, mg W , = weight of solvent saturated membrane, mg W , = weight of membrane after time t, mg Subscripts
Received for reuiew December 2, 1987 Revised manuscript receiued May 17,1988 Accepted June 4,1988
M = methanol G = ethylene glycol
Simulation and Optimization of Three-phase Distillation Processes Jeffrey P. Kingsley and Angelo Lucia* Department of Chemical Engineering, Clarkson University, Potsdam, New York 13676
The relationship between two-phase and three-phase solutions to heterogeneous distillation simulation and optimization problems is studied. Simulation results illustrate that bifurcation points can occur in the two-phase solutions and that, in regions of steady-state two-phase solution multiplicity, not all two-phase solutions can be used to compute a three-phase solution. Numerical difficulties are attributed to marked differences in the separation performed by the two- and three-phase solutions and to ill-conditioning when a second liquid phase is present in small amount. The optimization study shows that multiple local optima with different numbers of three-phase trays can occur. These local optima need not be spurious and indicate that some global optimization method is required to guarantee the correct solution is located. Accordingly, the tunneling algorithm for global optimization is used to solve a three-phase distillation optimization problem. Results show that tunneling can be used to successfully find the global solution. In this manuscript, we are primarily concerned with building an understanding of the various solutions to heterogeneous distillation processes. In particular, we are
* Author to whom correspondence should be addressed. 0888-5885/88/2627-1900$01.50/0
interested in the relationships between homogeneous (two-phase) and heterogeneous (three-phase) solutions to the same problem and in the associated implications of these relationships. Recently there have been a number of simulation papers 0 1988 American Chemical Society
Ind. Eng. Chem. Res., Vol. 27, No. 10, 1988 1901 initialize 8 solve homogeneous distillation
compute LLE on each homogeneous stage
I
2. The homogeneous solution(s) is often the starting point for heterogeneous distillation calculations, regardless of how the model equations are solved (see, for example, Ferraris and Morbidelli (1982)). No research, to our knowledge, addresses the relationship between these two facts, A second and perhaps more pragmatic objective was to use, if possible, any of this insight to strengthen existing algorithms for three-phase distillation calculations. The equations that model the jth stage of heterogeneous distillation column are the component material balances
perform VP VLLE flash on stages with false LLE behavior
solve heterogeneous distillation removing extraneous liquid phases
Figure 1. Heterogeneous distillation algorithm.
that report multiple two-phase solutions to distillations involving potentially heterogeneous mixtures (Magnussen et al., 1979; Prokopakis and Seider, 1983; Kovach and Seider, 1987) and many recent papers on three-phase distillation (Prokopakis et al., 1981; Ferraris and Morbidelli, 1981,1982; Kovach and Seider, 1987). However, no research has been done to connect the two. Furthermore, no studies concerned with the optimization of three-phase distillation processes have been conducted. Because mathematical analysis of these problems alone presents a formidable task, it seems reasonable to begin understanding the solution structure for heterogeneous distillation processes using a computational approach. The main results of this research are the following: 1. By use of the azeotropic distillation of ethanol and water with a benzene entrainer, parameterization in both Murphree stage efficiency and bottoms flow shows that bifurcation points in the two-phase solutions can occur. Solution curves from two-phase solutions to a three-phase solution to the same problem were also studied. Results show that some homogeneous solutions cannot be used to find the three-phase solution to a given problem because of large qualitative differences in separation and that this is particularly problematic in regions of steady-state multiplicity. 2. The level sets (or contours) of a representative objective function for the ethanol-water separation were used to study the optimization of heterogeneous distillations. Local optima were shown to occur on different steady-state solution surfaces. These optima are not spurious and suggest that some global analysis is required to ensure that the correct solution is found. 3. The tunneling method of global optimization (see, e.g., Levy and Montalvo (1985)) was applied to this three-phase distillation problem and was successful in finding the global solution from several reasonable starting points.
Simulation One objective of the simulation portion of this work was to build an understanding of two- and three-phase solutions to heterogeneous distillations. This study was motivated by two facts: 1. Multiple two-phase (or homogeneous) solutions to distillations involving mixtures that are potentially heterogeneous have been reported before (see Magnussen et al. (1979), Prokopakis and Seider (1983),and Kovach and Seider (1987)).
equilibrium equations
energy balance (1
+ Sf)LfHf+ (1 + Sf')Lf'Hf' + (1 + S)')VjH)' -
FjHT - Lf-lHj-1 - LfflHffi- Vj+iHs, + Qj = 0 (5)
summation equations
E(xj;i - yj,i) = 0
(7)
The stages in the column are numbered from top to bottom, and specifications are handled by including them directly in the equation set. The distillation algorithm that we use is given in flow-chart form in Figure 1. It is a modification of the algorithm given in Ferraris and Morbidelli (1981, 1982). All subproblems are solved by using a direct prediction Newton-like method (i.e., no line searching, trust regions, etc.) because we have found them to be reliable. However, this approach sometimes requires a careful initialization of the homogeneous distillation subproblem (see Kingsley (1986) for details). Also,unlike that of Ferraris and Morbidelli (1981,1982),our algorithm does not require the number of three-phase trays to be specified beforehand. Heterogeneous behavior on each tray is determined by a VP flash. Tangent plane analysis (see Baker et al. (1981) and Michelsen (1982a,b))could also be used; however, a VP flash was selected because temperature is a variable, but the energy balance need not be included in the equation set. There are many interesting questions that can be asked concerning the relationship between homogeneous and heterogeneous solutions to the same distillation. These include such questions as the following: 1. Do multiple homogeneous solutions exist for distillations that exhibit three-phase behavior on some of the stages in the column? 2. If so, are all homogeneous solutions connected to a single heterogeneous solution or are they connected to different three-phase solutions? Are some homogeneous solutions not connected to heterogeneous solutions? 3. What is the nature of solution multiplicity for multiparameterized homogeneous solutions to these dis-
1902 Ind. Eng. Chem. Res., Vol. 27, No. 10, 1988 Table I. Problem Statement for Ethanol Dehydration Column 42 stages partial condenser, condenser duty = 0 kJ/h bottom column pressure = 1.2156 X lo5 Pa top column pressure = 1.013 X lo5 Pa two liquid feeds: (1) stage 1 flow rate, 373.0 kmol/h; temp, 298.0 K (2) stage 4 flow rate, 100.0 kmol/h; temp, 354.0 K feed pressure equals pressure of feed stage feed mole fractions component feed 1 feed 2 ethanol 0.2627 0.8700 benzene 0.5737 0.0000 water 0.1635 0.1300 n
Table 11. Two-Phase Solutions to Ethanol Dehydration Column for Efficiency of 1.00 xB,EtOH
B, kmol/h 93.00 90.00 88.00 87.73" 87.50 87.00 86.89 86.868 86.863" 86.00 85.00 82.00
high purity 0.962 78 0.980 33 0.992 69 0.994 40 0.995 87 0.999 08 0.999 78 0.99991 0.999 92
intermediate
low purity
0.888 32 0.930 85 0.983 33 0.996 45 0.999 20 0.999 92
0.888 32 0.866 92 0.85704 0.848 26 0.847 79 0.847 68 0.833 31 0.821 74 0.798 48
" Approximate values at turning point. Table 111. Two-Phase Solutions to Ethanol Dehydration Column for Efficiency of 0.20
A I
V low purity SOlUtlon8
J
Figure 2. Homogeneous solutions.
tillations? Do multiple solutions occur as a consequence of regular turning points (as in the case of singly parameterized solutions (see Magnussen et al. (1979)))or do the solutions bifurcate? How are these special points of interest (turning points and bifurcations) in the homogeneous solution manifold related to the heterogeneous solution(s)? For the purpose of illustration, we have chosen the azeotropic distillation of ethanol and water with benzene as the entrainer. The definition of the column is given in Table I. In this work, the UNIQUAC equation was used to model liquid-phase activity coefficients, and the vapor phase was modeled by using the B-truncated virial equation of state (see Prausnitz et al. (1980)). Furthermore, the interaction parameters for the UNIQUAC equation (Table 4-4, p 65) and the vapor pressure equation and constants (Appendix C-2, p 150) given in Prausnitz et al. (1980) were used for all studies in this paper. Figure 2 shows the ethanol composition in the bottoms stream as a function of bottoms flow rate and efficiency (for the top four stages) for the two-phase solutions. Note that the geometry of the surface is very complicated. It contains regular turning points, steep fronts, and a bi-
B , kmol/h 89.00 87.00 85.75 83.75 83.45 83.25 83.00 82.75
xB,EtOH
0.968 17 0.979 71 0.987 19 0.999 54 0.981 25 0.964 90 0.942 28 0.921 50
B, kmol/h 82.72 82.69 82.66 82.62 82.59 82.50 82.00
xB,EtOH
0.919 50 0.917 62 0.915 84 0.91363 0.91208 0.907 91 0.892 50
furcation point. Also unless explicitly stated, a value of 1.00 was used for the Murphree tray efficiency on the remaining stages for all illustrations in this paper. The front edge (or boundary) of the surface shown in Figure 2 is the homogeneous solution curve for an efficiency of 1.00. Selected computational results used to generate this curve are presented in Table 11. This solution curve is labeled ABCD and contains two regular turning points, one each at points B and C. The approximate values of the bottoms flow rate at those turning points are 86.863 and 87.730 kmol/h. The corresponding ethanol compositions in the bottoms stream are 0.999 92 and 0.888 32 respectively. Any homogeneous solution on the top part of this curve (i.e., between points C and D) gives relatively high-purity ethanol in the bottom of the column and can be easily calculated by a direct prediction Newton-like method from a reasonable starting point. We have used the thermodynamically constrained hybrid method of Venkataraman and Lucia (1986) and a starting strategy based on the use of the short-cut techniques (see Venkataraman and Lucia (1988)) and have had no difficulty calculating these high-purity ethanol solutions. The calculations typically converged in 10 iterations. Furthermore, any low-purity ethanol solution between points A and B is also easily calculated (in approximately 10-15 iterations) by direct prediction Newton-like methods using an ideal solution starting point. The intermediate solutions between points B and C must be computed by using a continuation method and often require small steps in the continuation parameter and 30-50 Newton-like iterations to go from one solution to the next. The back edge of the surface is the solution curve for an efficiency of 0.20 and is labeled EFGHI. Selected two-phase solutions for this set of tray efficiencies are shown in Table 111. This solution curve has no turning points, but it does have an extremely steep front between points F and G. The high-purity ethanol solutions between points G and I are again easy to compute by using direct prediction Newton-like methods even though the ethanol composition goes through a maximum value at point H.
Ind. Eng. Chem. Res., Vol. 27, No. 10, 1988 1903 Table IV. Two- and Three-phase Solutions to Ethanol Dehydration Column for Efficiency of 0.50
B, kmol/ h 110.00 100.00 90.70 90.00 87.00 86.50 86.347 86.00 85.629" 85.60 85.572" 80.00 78.495 75.00 70.00
high purity 0.876 12 0.920 13 0.968 66 0.972 72 0.990 88 0.994 02 0.994 98 0.997 20 0.999 63 0.999 75 0.999 75
M(0)b intermediate
0.954 51 0.976 11 0.999 75
low purity
M(1) 0.870 76 0.926 51 0.975 14 0.979 18 0.997 19 0.992 21 0.969 66 0.92949 0.917 57 0.91690 0.916 29 0.874 86 0.869 06
0.95451 0.943 27 0.939 11 0.877 06 0.870 11 0.856 10 0.838 47
M(2) 0.863 14 0.922 32 0.977 75 0.981 71 0.999 30 0.977 67 0.967 78
M(3) 0.859 67 0.919 26 0.977 90
"Approximate values a t turning point. bM(J) is the solution curve (manifold) for which there are J three-phase trays in the column.
The low-purity ethanol solutions are easily calculated in a similar manner. To calculate solutions along the steep front, a continuation method must be used. In fact, to span the bottoms flow rate range from 82.62 to 82.75 kmol/h, a step size of 0.002 kmol/h and as many as 46 Newton-like iterations were required. The collection of these steep fronts gives the appearance of a steep ledge in the surface. There is a bifurcation point at an efficiency of 0.447 and a bottoms flow rate of 85.52 kmol/h. It is shown as point J in Figure 2. In the neighborhood of the bifurcation point, it is quite easy to jump from the high-purity solutions to the low-purity solutions without noticing any apparent discontinuity in ethanol composition. We now turn our attention to the heterogeneous solutions. Here we are interested in studying the geometry of two-phase and three-phase solutions to the same problem, particularly in regions where there are multiple homogeneous solutions. We are also interested in the numerical implications of the relationship between these solutions. Figure 3 shows projected plane curves for the ethanol, benzene, and water distillation defined in Table I. In this figure, the efficiency of the top four trays was fixed at 0.50, and the ethanol composition in the bottoms stream for various solutions (both homogeneous and heterogeneous) was parameterized in bottoms flow rate over the range 70.0-110.0 kmol/h. These solution curves are labeled M(O), M(1), M(2), etc., where the number in parentheses indicates the number of three-phase trays in the column for a given curve (manifold). This particular efficiency was selected because it is reasonably close to the value of efficiency that gives the singular bifurcation in the homogeneous case &e., t = 0.447). Numerical results for some of these homogeneous and heterogeneous solutions are given in Table IV. The left most curve labeled M(0) is the plane curve for the two-phase solutions. I t contains two regular turning points at A and B. The bottoms flow rate has values of 85.63 and 85.57 kmol/h at points A and B, respectively, and the corresponding ethanol compositions are 0.954 51 and 0.999 75. Note that the low-purity and high-purity solutions are beginning to approach each other because we are in the neighborhood of the bifurcation point. The next curve to the right, M(1), is the solution curve for the heterogeneous solutions in which the top stage is a three-phase tray. It splits off from M(0) at point C, exhibits no turning points, goes through a maximum ethanol composition value Gust like the homogeneous curve), and then follows the qualitative behavior of the high-purity
8
-ti
1.00
-
0.95
-
0.QO-
x
0.80 70.
80.
9 0.
100.
6 (kgmdlhr)
Figure 3. Plane curves for homogeneous and heterogeneous solutions.
homogeneous solutions at high values of bottoms flow rate. Similarly, the M(1) curve gives rise to M(2), the solutions in which the top two stages are three-phase trays, at point D. The M(2) curve has the same characteristics as MU). Finally M(2), in turn, gives birth at E to solutions in which the top three stages are three-phase trays. It is denoted by M(3). However, unlike the branching in the previous curves, this solution appears after M(2) has passed through its maximum ethanol composition. We stopped this part of the numerical investigation here; however, there is no reason to assume that this pattern of bifurcation (i.e., propogation of solutions with one more three-phase tray) stops. Convergence to the three-phase solution with the correct number of three-phase trays by Newton-like methods can depend on several factors. First, the separation performed by a heterogeneous solution must be qualitatively similar to the one performed by the homogeneous solution from which the computations are initiated. Second, no heterogeneous solution (intermediate or final) should have a tray in which the second liquid phase on that tray is present
1904 Ind. Eng. Chem. Res., Vol. 27, No. 10,1988
L from the corresponding homogeneous solution at point multiplicity, only the lowlpurity solution can be used to find the heterogeneous solution. There is a straightforward explanation for this. In cases where the calculations successfully transform a homogeneous solution into a heterogeneous one, the separations for all solutions are qualitatively the same. When at least two solutions give qualitatively different separations, numerical difficulties can occur. For example, the separation denoted by point K is a high-purity two-phase solution. The one denoted by point L is the corresponding heterogeneous solution for the same specifications and is a relatively low-purity solution. This difference in separation causes numerical difficulties in going from the homogeneous to the heterogeneous solution in this situation because the two-phase distillation solution is a poor starting point for the three-phase problem. Ill-conditioning can also cause numerical difficulties when one of the solutions along the path from a homogeneous solution to the fiial heterogeneous solution is near a point where a new three-phase tray appears. Newton-like methods will normally break down in this situation. Poor Newton-like directions result in either the phase splitting calculations or the subsequent distillation calculations because of this, and the algorithm usually terminates with an incorrect number of three-phase trays. Points C, D, and E are examples of ill-conditioning difficulties. The conditions under which convergence difficulties occur are independent of the number of two-phase solutions or the number of three-phase trays in the final heterogeneous solution. However, the fact that qualitatively similar separations are needed for convergence suggests that, in regions of two-phase solution multiplicity, some homogeneous solutions may be good starting points and others may be very poor starting points, and this becomes more pronounced as the homogeneous solutions depart from a bifurcation point. Circumventing computational difficulties in the simulation of heterogeneous distillation processes is a challenging task. The complex geometry of the homogeneous and heterogeneous solutions makes it difficult to handle large qualitative differences in separation and/or ill-conditioning in the distillations. Figure 4, which represents the composite homogeneous and heterogeneous solutions parameterized in bottoms flow rate and efficiency of the top four stages, shows just how complex this geometry can be. It also suggests that, at this time, parameterization (see Kovach and Seider (1987)) may be the only way to routinely solve these problems. However, parameterization is not without its share of difficulties. Locating the correct
M(O)
t
x'EtoH'
Figure 4. Composite homogeneous and heterogeneous solutions.
path to follow can present problems and formidable computational efforts are sometimes required. We would like to end by partially answering the questions posed at the beginning of this section and with some comments on the utility of the algorithm presented in this section. 1. Clearly multiple homogeneous solutions exist for distillations that exhibit three-phase behavior for identical sets of specifications. 2. 0ur.limited experience seems to indicate that all homogeneous solutions ultimately lead to the same heterogeneous solution (see Figure 3). Of course, there can be spurious heterogeneous solutions with an incorrect number of three-phase trays along the way. Spurious or not, we have been unable to produce multiplicity in the solutions of any heterogeneous distillation. Thus, we believe that, in regions of homogeneous solution multiplicity, all two-phase solutions are connected to the same heterogeneous solution. However, not all homogeneous solutions can be used to compute a heterogeneous solution. 3. Multiparameterizedhomogeneous solutions can exhibit both regular turning points and bifurcations. In the particular distillation studied in this paper, the bifurcation point is a singular bifurcation. It is unclear at this time what the relationship between these so-called singular points (i.e., turning and bifurcation points) in the homogeneous solutions and heterogeneous solutions to the same problem is, since three-phase solutions exist outside the range of specifications that give two-phase solution multiplicity. 4. In general, our Newton-like algorithm for heterogeneous distillation process simulation is quite useful. We have used it to successfully simulate the relatively few examples that appear in the literature. These examples include the 1-propanol, 1-butanol, water column originally studied by Block and Hegner (1976); the 2-propanol, cyclohexane, water column in Ross and Seider (1980); and the acetonitrile, acrylonitrile, water separation presented by Ferraris and Morbidelli (1982). Furthermore, the
Ind. Eng. Chem. Res., Vol. 27, No. 10, 1988 1905 Table V. Profit Function Values for Two-Phase Solutions to Ethanol Dehydration Column for Efficiency of 0.80 xB.EtOH
B, kmol/h
high purity
93.00 90.00 88.00 87.1745 87.10 87.00 86.81 86.785O 86.00 85.00
0.962 39 0.979 87 0.992 19 0.997 45 0.997 92 0.998 56 0.999 78 0.999 91
a Approximate
intermediate
0.940 08 0.957 96 0.971 19 0.996 13 0.999 91
low purity
0.94008 0.925 38 0.918 55 0.910 29 0.909 41 0.890 84 0.876 47
intermediate
low purity
1323.40 1457.01 1553.07 1594.36 1598.10 1603.14 1612.45 1612.83
108.85 561.14 891.86 1517.29 1612.83
108.85 -275.99 -460.95 -690.89 -715.70 -1264.52 -1722.01
values at turning point.
modifications that we have made to the approach of Ferraris and Morbidelli (1981,1982) show that the number of three-phase trays does not need to be specified a priori but can be correctly determined during the course of solution. There is also no need for line searching or switching to equation decoupling techniques as suggested by Ferraris and Morbidelli. In contrast, parametric continuation methods are useful for computing intermediate solutions and solutions along extremely steep fronts but should be used judiciously. Certainly, a large amount of analysis and numerical experimentation is still needed to find ways of improving both our fundamental understanding of heterogeneous distillation and the reliability and efficiency of associated computer tools. We are, in our opinion, a long way from establishing global results for distillation similar to those that exist for the single-stage flash process (e.g., tangent plane analysis).
Optimization This section of the manuscript is concerned with the optimization of three-phase distillation processes, an area that, to ow knowledge, has not been addressed in the open literature before. For the purposes of continuity and illustration, we have again used the ethanol, benzene, water column defined in Table I. The objective function in this investigation is given by fp
high purity
lurninp Pt
.
solutions
C *
IOW
purity solutions
Y
I
I(Et OH)
tv
Figure 5. Contours of profit function on homogeneous solutions.
= c ~ B [ x ( E ~ O H-) c] ~~ .Q~ R / ~ ( E ~ O H ) (8)
where f p is a profit function and c1 and c2 are the selling price for 200 proof ethanol and an operating cost coefficient, respectively. Although c1 is largely dependent on state taxes, we have used a value of $100/kmol. Furthermore, the operating costs have been correlated to the reboiler duty and ethanol composition in the bottoms stream, and we note that the term QR/x(EtOH)in eq 8 increases as r(Et0H) increases. The coefficient c2reflects, among other things, the cost of low-pressure steam and labor. The value of c2 was $4.5 X 104/kJ. Finally, B , x(EtOH), and QR are the bottoms flow rate, ethanol mole fraction in the bottoms stream, and reboiler duty, respectively. The equality constraints in the optimization problem are eq 1-5 rewritten in terms of component flow rates plus the specification equation Q, = 0 kJ/h, where Q, is the condenser duty. They are highly nonlinear. The inequality constraints are bounds on variables, particularly bounds on component flow rates in product streams. For this problem, the most important bounds are those on the ethanol, benzene, and water component flow rates in the bottoms stream. We have used an upper bound of 88.9 kmol/h for ethanol and values of 1 X 10-lo
Table VI. Profit Function Values for Three-phase Solutions to Ethanol Dehydration Column for Efficiency of 0.80
B, kmol/ h
XB,EtOH
fp, $/h
95.00 93.00 90.00 88.00 87.70 87.60 87.54 87.50 87.40 87.10
0.958 23 0.968 72 0.985 35 0.997 09 0.998 89 0.999 49 0.999 83 0.998 67 0.993 20 0.974 09
1427.46 1492.94 1601.17 1680.24 1692.40 1696.31 1698.20 1663.23 1508.77 978.15
and 1 X 10-lokmol/h for the lower bounds on benzene and water, respectively. Also, the efficiency of the top four trays in the column was allowed to vary between 1.0 and 0.20. To study this optimization problem, we studied the behavior of the level sets (or contours) of the objective function on various constraint surfaces. Figures 5 and 6 show the contours of the objective function on surfaces M(0) and M(1). Tables V and VI contain some of the numerical solutions used to generate these figures. Similar
1906 Ind. Eng. Chem. Res., Vol. 27, No. 10, 1988
uppor turnlng polnt manifold
lowaf turning point manifold
Figure 6. Contours of profit function on heterogeneous solutions.
figures can be constructed for M(2), M(3), etc. Note that there is an optimum point on each of the surfaces M(0) and M(1). The optimum on M(0) is very close to the collection of turning points for the high-purity solutions. It is actually located at the point A = (x(EtOH), B, 7) = (0.99991, 86.888, 0.923) where the value of the objective function is $1623.59/h. The maximum value of fp on M(l) is $1700.5/h and is very close to the ridge of maximum ethanol composition on this surface. It is located at A = (0.99983,87.60,0.770) in Figure 6. At this solution, the amounts of the two liquid phases on tray one are 65.27 and 360.14 kmol/h, respectively. The maximum of fp on M(2) lies on the boundary between M(l) and M(2) at the point (x(EtOH), B, 7) = (0.9997, 87.29, 0.628). Figure 7 shows the level sets of fp on the composite surface formed from the intersection of M(0) and M(1). The important thing to note is that neither optimal solution is spurious. Phase splitting calculations or tangent plane analysis applied to all stages shows that the local maximum at point A in Figure 7 truly operates in the homogeneous regime. Point B in that figure, in which the top stage in the column is a three-phase tray, is the global maximum.
Preliminary Optimization Calculations We tried solving this optimization problem using a local optimization method. The algorithm used for these preliminary optimization calculations is a modification of the successive quadratic programming approach of Wilson (1963), Han (1976), and Powell (1978). It is a Newton-like algorithm that has been tailored to meet specific application needs in the area of separation process optimization. It uses a feasible starting point and solves a linear leastsquares problem for initial values for the Lagrange multipliers. Hessian matrix approximations are built using some analytical second derivative information and thermodynamically constrained quasi-Newton updates. Finally Powell's line search function and a quadratic programming method based on an active set strategy are used. Details associated with automatic starting point generation, initial multiplier estimation, Hessian matrix approximation,and line searching strategy can be found in the recent paper by Kumar and Lucia (1987). For these studies, the cal-
1
x ( Et OH)
Figure 7. Contours of profit function on composite solutions.
culations were not restricted to any particular surface. Rather, the Newton-like algorithm was allowed to add and delete three-phase trays during the course of iteration. We felt this was a reasonable strategy because the line search function that is used forces the constraints to be satisfied fairly well even though the path is in principle an infeasible one. The results for these calculations were mixed. If the algorithm was initialized with a reasonable two-phase distillation solution, then it usually converged to the local maximum at point A in Figure 7 without much difficulty provided there was no change in the number of phases on any tray in the column during the course of iteration. For example, starting at the high-purity solution (x(EtOH),B, 7) = (0.96278, 93.0, 1.00), the optimization algorithm converged to point A in Figure 7 in 21 iterations and 25 function and gradient calls. The objective function value improved from $1331.71/h to $1623.59/h. In contrast, line searching difficulties or loss of descent direction and subsequent failure often occurred if the starting point was a low-purity solution and forced the optimization algorithm to follow the steep part of the constraint surface. The same behavior was observed if a high-purity three-phase solution in which the top tray is a three-phase tray was used as the starting point and provided there was no change in the number of three-phase trays during the course of iteration. Convergence to the global solution was observed, but line searching and loss of descent direction difficulties and failure occurred if the steep side of the constraint surface had to be followed. By far the most serious difficulties occurred when the optimization algorithm attempted to add or delete a phase to a tray in the column. For example, the addition of a second liquid phase to the top stage resulted in a change in the number of constraints in the active set in the form of n, phase equilibrium equations. The mass balance and energy balance equations for the top stage also had to be
,
Ind. Eng. Chem. Res., Vol. 27, No. 10, 1988 1907 modified to accommodate the new liquid phase. The existing Lagrange multipliers for this stage had to be modified, and some new ones had to be estimated. This is where serious computational difficulties arose. The new iterate and corresponding updated multipliers were often a poor solution approximation and resulted in very poor Newton-like steps. Furthermore, even if the estimation of the new multipliers was reasonable, the optimization calculations frequently encountered difficulties due to the ill-conditioning of the constraint Jacobian, large constraint violations, or discontinuities in the line search function caused by the addition of the second liquid phase. Many of these difficulties manifested themselves as line search or loss of descent direction failures. In general, the algorithm frequently failed in this situation. The occurrence of multiple optima in heterogeneous distillation and the serious algorithmic difficulties that it presents for local methods suggests that some global optimization method is required to help guarantee that the correct solution is found. In the next section, we briefly describe the tunneling method for global optimization.
The tunneling method of global optimization is a general approach that calculates a constrained global minimum by alternating between a minimization phase and a tunneling phase until a better solution cannot be found. It is usually convenient to implement the tunneling phase using the same constrained minimization algorithm that is used for the minimization phase, although this is not strictly necessary. The basic steps of the tunneling method are iterative and given in the following algorithm:
The Tunneling Method of Global Optimization In comparison to local methods such as the generalized reduced gradient (GRG) method (see Wolfe (1967), Abadie and Carpentier (1969)) and successive quadratic programming (SQP) algorithms (Wilson, 1963; Han, 1976; Powell, 1978), there has only been a modest effort in developing global methods for optimization. Consequently, there are relatively few general techniques available for solving such problems, and many applications have relied on specialpurpose methods for finding global solutions. Good introductions as well as various surveys of the field can be found in Dixon and Szego (1975) and Boggs et al. (1985). Methods for global optimization are usually either deterministic or stochastic, although some features of both are typically used in problem solving. Within the class of deterministic methods, there are trajectory (relaxation) methods (see Branin (1972),Branin and Hoo (1972),Hardy (1975)), methods based on regions of attraction (see Treccani et al. (1972) and Corles (1975)), and function modification techniques (see Goldstein and Price (1971), Levy and Montalvo (1985), and Ge and &in (1987)). Stochastic methods, which involve some random or probabilistic aspects, include pure random search techniques, clustering methods (Kan and Timmer, 1985), and methods based on stochastic gradient differential equations (see Pentini et al. (1985)). Some of these techniques have found application in chemical engineering. For example, Van Dongen and Doherty (1983) have used trajectory methods to find various local minima in the Gibbs free-energy function for vapor-liquid equilibria involving heterogeneous mixtures. Moreover, relaxation methods have been used for solving steady-state simulation problems since the late 1950s (see Rose et al. (1958), Jelinek and Hlavacek (1976), and Ketchum (1979)). The tangent plane analysis of Gibbs (see Michelsen (1982a,b)) is a function modification method, and random search techniques have also been used to solve chemical engineering optimization problems (see Luus (1975),Gaines and Gaddy (1976),and Ballman and Gaddy (1977)). Steady-state process optimization problems are frequently large, nonlinearly constrained problems with few degrees of freedom. In light of this, the tunneling method (see Levy and Montalvo (1985)) appears to be the most promising because it has been tested on small unconstrained and nonlinearly constrained mathematical benchmark problems and given reasonable results.
where c ( x ) is a vector of constraints and m’is the number of equality constraints. Denote this solution by x i and its corresponding objective function value by f ( x i ) .
tunneling algorithm 1. Set k = 1, and choose x i 2. Starting with xg, calculate
min f ( x )
(9)
subject to j = 1,...,m’
(10)
j = m’+l, ...,m
(11)
cj(x) = 0
cj(x) I 0
3. Choose
and solve min T [ x , x L f ( x ; ) l
(12)
subject to j = 1,...,m’
(13)
j = m’+l, ...,m
(14)
cj(x) = 0 cj(x) I 0
for x k such that f ( X k ) I f ( x ; ) 4. If x k = x i and f ( X k ) I f ( x ; ) , set = x k , set k = k + 1,and go to step 2. Otherwise, stop. The user must supply an initial estimate of the solution for the first minimization phase. Step 2 of the algorithm is the minimization phase. It calculates a local constrained minimum and defies the tunneling level to be used in step 3 of the algorithm. Note that the form of the tunneling function forces the objective function value to improve and prevents x i from being located again because the denominator goes to zero as x approaches x i . This technique is called deflation. The pole strength in step 3 is determined using the heuristic rules given in Levy and Montalvo (1985, p 26), and the unknown variables in the tunneling phase are initialized by perturbing the current local minimizer x i . It may be necessary to try a variety of starting points here, and the number required and way in which they are computed may be strongly problem dependent. Also for problems in which several minima exist at the same tunnelin# level (i.e., x ; , x i , ..., $ such that f ( x i ) = f ( x 3 = ... = f ( x i ) ) , there are provisions in the algorithm for deflating each of these minima to avoid cycling between them. This is not common in our applications and has been omitted here for simplicity. Those interested in this aspect of the algorithm are referred again to Levy and Montalvo (1985). Also we caution the reader to note that although the tunneling phase can be terminated whenever f ( x k ) I f ( x i ) , it must be understood that for constrained problems both the points and z k must, in principle, be feasible points. Thus, some care must be exercised if the minimization
1908 Ind. Eng. Chem. Res., Vol. 27, No. 10, 1988
algorithm being used is an infeasible path altorithm like successive quadratic programming. If the tunnel exit (i.e., the point x k ) is a better estimate of the global minimizer, then it is used as the starting point for the next minimization phase and the algorithm proceeds iteratively. If, after trying various starting points for the tunneling phase, a better estimate cannot be found, then the point x i is assumed to be the global minimizer with corresponding objective function value f ( x i ) and the algorithm is terminated. Again, this is a subjective decision that is often best made by the user because of its problem dependence.
An Application of Tunneling to Three-phase Distillation Because of the rather mixed results obtained using just the modified successive quadratic programming algorithm, we decided to try the tunneling method on the ethanol dehydration column. We were particularly interested in seeing if tunneling could be used to help predict solutions with different numbers of three-phase trays in the column, since this is a situation that presents the most serious challenge to the successive quadratic programming algorithm. Because the profit function for the ethanol-water separation exhibits a unique maximum for a given number of three-phase trays and because the tunneling algorithm is not designed to handle constrained problems in which the number of constraints changes from one minimization phase to the next, some modification in the basic algorithm is required. The algorithm must include some provision for changing the number of phases on any given tray in the column. The starting point in step 1of the algorithm is chosen as any high-purity two-phase solution to the column calculated by solving a simulation problem with selected values of certain decision variables. (For this particular problem, the total bottoms flow rate and the efficiency in the recityfing section of the column were designated as decision variables.) The optimal solution to the two-phase column is computed by using the SQP algorithm with the number of phases on each tray in the column held fixed, and the tunneling phase is entered. If other local optima for the two-phase column existed, the tunneling algorithm would proceed according to the steps outlined in the previous section. However, for this application, the exit to the tunneling phase cannot be chosen according to step 4 of the algorithm because the condition x k = x i such that f ( x k ) I f ( x i ) will never occur because for a fixed number of three-phase trays the optimal solution is unique. Instead we check for liquid-phase instability on each tray in the column during each iteration of the tunneling phase. If liquid-phase instability is detected, a simulation problem is solved by using the algorithm outlined in Figure 1 with the current iterative values of the decision variables held fixed. If the number of three-phase trays at the solution of the simulation problem does not change, then tunneling is continued. If the number of three-phase trays does change, then the objective function (-fp in this case) is evaluated at each feasible solution. Let m be the current number of threephase trays from the previous minimization phase and denote the feasible solutions to the simulation problem by x i , x i , ...,x i , where the superscrpt denotes the number of three-phase trays in the column. Select the first index j = m+l, ..., p such that f ( x h ) < f(x!). If this condition cannot be met, continue the tunneling phase. If f(x6) < f ( x r ) for some j = m+l, ...,p , then the next minimization phase is begun with xi+l = x i and the number of threephase trays is held fixed at j . This procedure is repeated until sequential minimization phases fail to produce any
improvement in the objective function.
Numerical Results In this section, we summarize the numerical results obtained by applying the tunneling algorithm to the ethanol dehydration column. The results reported here are for a particular case; however, they are typical of many of the optimization calculations performed on this column. The initial minimization (-fp) phase was started with a high-purity two-phase solution calculated by solving a simulation problem with the specifications Q, = 0 kJ/h, B = 90 kmol/h, and q = 1.00. These calculations converged in nine iterations using the thermodynamically constrained hybrid method and starting strategy given in Venkataraman and Lucia (1988). The corresponding ethanol purity in the bottom of the column and reboiler duty were x (EtOH) = 0.98033 and Q R = 1.55729 X lo7 kJ/h, respectively. The value of the profit function at this starting point was $1466.63/h. The minimization calculations required 14 iterations and 16 function and gradient calls using the algorithm of Kumar and Lucia (1987) and a tolerance on the two-norm of the Kuhn-Tucker conditions of to reach point A in Figure 7 where x ; = (x(EtOH), B, 7) = (0.99991, 86.888, 0.923) and f p ( x i ) = $1623.59/h. The initial tunneling phase was started by finding XI = 0.6 using the rules in Levy and Montalvo and perturbing the decision variables. The total bottoms flow rate (with the composition fixed) and the efficiency were changed by 0.2 kmol/h and -0.10, respectively. All other variables were held fixed at their values given by x i . The modified SQP algorithm of Kumar and Lucia was also used for the tunneling calculations, and on the second iteration of tunneling, a liquid-phase split was detected on tray one in the column. The values of the total bottoms flow rate and efficiency on this iteration were B = 88.73 kmol/h and 7 = 0.85. This iterate is an infeasible point where the 2-norm of the constraints is satisfied to an accuracy of 7.6 X Next, a three-phase distillation simulation problem was solved (to an accuracy of lo+) with the specifications Q, = 0. kJ/h, B = 88.73 kmol/h, and q = 0.85. These calculations required three iterations to solve the two-phase problem and five iterations to solve the simulation problem in which the top tray in the column was a three-phase tray. The corresponding profit function values were $1526.76/h and $1631.88/h, respectively. Since the addition of a three-phase tray resulted in an improvement in the profit function, the next minimization phase was entered. The three-phase solution to the simulation problem was used as the starting point for the second minimization with the number of three-phase trays held fixed at one. These calculations converged in 27 iterations and 31 function and gradient calls to point B in Figure 7. This is the global solution and given by xd = (x(EtOH), B, 7) = (0.99983, 87.560, 0.770) and f p ( x d ) = $1700.5/h. A second tunneling phase was initialized with a tunneling level of -f,(x;) = 1700.5 and X2 = 0.6. The unknown variables for the tunneling calculations were initialized by perturbing the total bottoms flow rate and efficiency by 0.2 kmol/h and -0.10, respectively, and by keeping the remaining unknown variables fixed at their values at XI. After four iterations of tunneling, a liquid-phase split was detected on tray 2 of the column. Here the values of the bottoms flow rate and efficiency were 91.55 kmol/h and 0.67, and the constraints were satisfied to an accuracy of 8.6 x A three-phase simulation problem was solved with the specifications Q, = 0. kJ/h, B = 91.55 kmol/h, and 9 = 0.67. Three and four iterations were required to solve the
Ind. Eng. Chem. Res., Vol. 27, No. 10, 1988 1909 simulation problems for a solution with one and two three-phase trays, respectively. The corresponding values for these feasible solutions were $1508.62/h and $1516.74/h, respectively. Since the addition of a second liquid phase on tray 2 resulted in an improvement in the profit function, a third minimization phase was entered. The solution to the simulation problem was used as the starting point for these calculations, and the number of three-phase trays was fixed at two. Ten iterations and 10 function and gradient calls were required before the SQP algorithm terminated because of line searching difficulties. Further inspection of the tenth iterate from the optimization calculations showed that the second liquid phase on tray 2 in the column had a total flow rate of 0.067 kmol/h and that the next iterate wanted to remove the second liquid phase on that tray. This discontinuity in the number of phases resulted in line searching failure, and it was concluded that the local optimum here was at the point where there is only one liquid phase on tray 2. Furthermore, since the corresponding value of the profit function at this point was $1665.20/h and was less than f & f ) = $1700.5/h, x i was assumed to be the global optimal solution. Finally, we note that a number of other starting points for the third minimization gave the same result that xf was the global solution.
Acknowledgment This work was supported by the Office of Basic Energy Science, U.S.Department of Energy, under Grant DEFG02-86ER13552. The authors also thank T. Zukas and Dr. S. Venkataraman for their help with some of the numerical calculations.
Nomenclature B = total bottoms flow rate c ( x ) = constraint functions f ( x ) , f p = objective function, profit function F = total feed flow rate H = molar enthalpy K = equilibrium ratio I = vector of liquid molar flow rates L = total liquid flow rate n, = number of components p = pressure Q = heat duty S = side-stream total flow rate fraction T(x, z*, A) = tunneling function V = total molar vapor flow rate x , xk, x* = vector of liquid mole fraction, iterate, optimal value x(Et0H) = ethanol mole fraction in bottom stream y = vapor mole fraction z = feed mole fraction Subscripts
B = bottoms i = component index j = stage index k = iteration number Superscripts F = feed L = liquid phase V = vapor phase
I = first liquid phase I1 = second liquid phase Greek Symbols X = pole strength r)
= Murphree stage efficiency
Registry NO. EtOH, 64-17-5;C&,
71-43-2.
Literature Cited Abadie, J.; Carpentier, J. “Generalization of the Wolfe Reduced Gradient Method to the Case of Nonlinear Constraints”. In Optimization; Fletcher, R., Ed.; Academic: London, 1969. Baker, L. E.; Pierce, A. C.; Luks, K. D. “Gibbs Energy Analysis of Phase Equilibria”. SPE/DOE Second Joint Symposium on Enhanced Oil Recovery, Tulsa, OK, Society of Petroleum Engineers, Dallas, TX, 1981;paper SPE/DOE 9806. Ballman, S. H.; Gaddy, J. D. “Optimization of Methanol Process by Flowsheet Simulation”. Ind. Eng. Chem. Process Des. Dev. 1977, 16,337-341. Block, U.; Hegner, B. “Development and Application of a Simulation Model for Three-phase Distillation”. AIChE J. 1976,22,582-589. Boggs, P. T.; Byrd, R. H.; Schnabel, R. B. Numerical Optimization 1984;SIAM: Philadelphia, PA, 1985. Branin, F. H. “Widely Convergent Method for Finding Multiple Solutions of Simultaneous Nonlinear Equations”. IBM J. Res. Dev. 1972,504-522. Branin, F. H.; Hoo, S. K. “A Method for Finding Multiple Extrema of a Function of n Variables”. In Numerical Methods of Nonlinear Optimization;Fletcher, R., Ed.; Academic: London, 1972. Corles, C. R. “The Use of Regions of Attraction to Identify Global Minima”. In Towards Global optimization; Dixon, L. C. W., Szego, G. P., Eds.; North-Holland: Amsterdam, 1975. Dixon, L. C. W.; Szego, G. P. Towards Global Optimization; North-Holland: Amsterdam, 1975. Ferraris, G. B.; Morbidelli, M. “Distillation Models for Two Partially Immiscible Liquids”. AZChE J. 1981,27,881-888. Ferraris, G. B.; Morbidelli, M. “An Approximate Mathematical Model For Three-phase Multistaged Separators”. AIChE J. 1982, 28,49-55. Gaines, L. D.; Gaddy, J. L. “Process Optimization by Flow Sheet Simulation”. Znd. Eng. Chem. Process Des. Deu. 1976, 15, 206-21 1. Ge, R. P.; Bin, Y. F. “A Class of Filled Functions for Finding Global Minimzers of a Function of Several Variables”. JOTA 1987,54, 241-252. Goldstein, A. A.; Price, J. F. “On Descent from Local Minima”. Math. Comp. 1971,25,569-574. Han, S.P. “Superlinearly Convergent Variable Metric Algorithms for Nonlinear Programming Problems”. Math. Prog. 1976,11, 263-282. Hardy, J. W. “An Implemented Extension of Branin’s Method”. In Towards Global Optimization;Dixon, L. C. W., Szego, G. P., Eds.; North-Holland Amsterdam, 1975. Jelinek, J.; Hlavacek, V. “Steady-State Counter Current Equilibrium Stage Separation with Chemical Reaction by Relaxation Method”. Chem. Eng. Commun. 1976,2,74-85. Kan, A. H. G. R.; Timmer, G. T. “A Stochastic Approach to Global Optimization”. In Numerical Optimization 1984;Boggs, P. T., Byrd, R. H., Schnabel, R. B., Eds.; SIAM Philadelphia, PA, 1985. Ketchum, R. G. “A Combined Relaxation Newton Method as a Global Approach to the Computation of Thermal Separation Processes”. Chem. Eng. Sci. 1979,34,387-395. Kingsley, J. P. “Heterogeneous Flash and Distillation Calculations Using Newton’s Method”. Ph.D. Dissertation, Clarkson University, New York, 1986. Kovach, J. W.; Seider, W. D. “Heterogeneous Azeotropic Distillation-Homotopy Continuation Methods”. Comput. Chem. Eng. 1987,11,593-605. Kumar, A.; Lucia, A. “Separation Process Optimization Calculations”. IMACS J. Appl. Num. Math. 1987,3,409-425. Levy, A. V.; Montalvo, A. “The Tunneling Algorithm for the Global Minimization of Functions”. SZAM J. Sci. Stat. Comput. 1985, 6,15-29. Lucia, A.; Kumar, A. “Distillation Optimization”. Comput. Chem. Eng. 1988,in press. Luus, R. “Optimization of Multistage Recycle Systems by Direct Search”. Can. J . Chem. Eng. 1975,53,217-220. Magnussen, T.; Michelsen, M. L.; Fredunslund, A. “Azeotropic Distillation Using UNIFAC”. Inst. Chem. Eng. Symp. Ser. No. 56,Thrid Intl. Symp. on Distillation; I C E Rugby, England, 1979; 4.2/1-19. Michelsen,M. L. “The Isothermal Flash Problem, Part I. Stability”. Fluid Phase Equilib. 1982a,9,1-19. Michelsen, M. L. “The Isothermal Flash Problem. Part 11. Phase Split Calculation”. Fluid Phase Equilib. 1982b,9,21-40.
I n d . Eng. Chem. Res. 1988, 27, 1910-1916
1910
Pentini, F. A.; Parisi, V.; Zirilli, F. “Global Optimization and Stochastic Differential Equations”. JOTA 1985, 47, 1-16. Powell, M. J. D. “A Fast Algorithm for Nonlinearly Constrained Optimization Problems”. In Lecture Notes in Mathematics; Watson, G. A,, Ed.; Springer-Verlag: Berlin, 1978. Prausnitz, J. M.; Anderson, T. F.; Grens, E. A.; Eckert, C. A.; Hsieh, R.; O’Connell, J. P. Computer Calculations for Multicomponent Vapor-Liquid and Liquid-Liquid Equilibria; Prentice-Hall: Englewood Cliffs, NJ, 1980. Prokopakis, G. J.; Seider, W. D. “Feasible Specification in Azeotropic Distillation”. AZChE J . 1983, 29, 49-60. Prokopakis, G. J.; Seider, W. D.; Ross, B. A. “Azeotropic Distillation Towers With Two Liquid Phases”. Foundations of Computer Aided Chemical Process Design Vol. II.; Mah, R. S. H., Seider, W. D., Eds.; AIChE: New York, 1981. Rose, A.; Sweeny, R. F.; Schrodt, V. N. “Continuous Distillation Calculations by Relaxation Method”. Ind. Eng. Chem. 1958,50, 737-740. Ross, B. A.; Seider, W. D. “Simulation of Three-phase Distillation Towers”. Comput. Chem. Eng. 1980,5, 7-20. Stadtherr, M. A.; Chen, H. S. “Strategies for Simultaneous Modular Flowsheeting and Optimization”. In Foundations of Computer
Aided Process Design; Westerberg, A. W., Chien, H. H., Eds.; CACHE: Ann Arbor, MI, 1984. Treccani, G.; Trabattoni, L.; Szego, G. P. “A Method for Isolation of Minima”. In Minimization Algorithms, Mathematical Theories and Computer Results; Szego, G. P., Ed.; Academic: London, 1972. Van Dongen, D. B.; Doherty, M. F. ’Calculation and Stability of Multiple Equilibrium Points for Nonideal Mixtures”. Fluid Phase Equilib. 1983, 14, 335-343. Venkataraman, S.; Lucia, A. “Exploiting The Gibbs-Duhem Equation in Separation Calculations”. AZChE J . 1986,32, 1057-1066. Venkataraman, S.; Lucia, A. “Solving Distillation Problems by Newton-like Methods”. Comput. Chem. Eng. 1988, 12, 55-69. Wilson, R. B. “A Simplified Algorithm for Concave Programming”. Ph.D. Dissertation, Harvard University, Cambridge, MA, 1963. Wolfe, P. “Methods in Nonlinear Programming”. In Nonlinear Programming; Abadie, J., Ed.; Wiley: New York, 1967. Received for review November 19, 1987 Revised manuscript received June 16, 1988 Accepted July 6, 1988
Studies on Gas Holdup in a Bubble Column Operated at Elevated Temperatures Zou Renjun*+ Hebei Academy of Sciences, Shijiazhuang, The People’s Republic of China
Jiang Xinzhen, Li Baozhang, Zu Yong, and Zhang Laiqi Chemical Engineering Department, Northwestern University, Xi”.,
The PeopleS Republic of China
Gas holdup was studied in a bubble column operated a t elevated temperatures for the systems air-water, air-alcohol, and air-5% NaCl solution. The influence of temperature on the gas holdup was studied focally nearing the boiling points. The bubble column is 100 mm in diameter and 1.05 in height, and the gas distributor has a single nozzle that is 10 mm in diameter. The gas and liquid flowed upward concurrently through the column. The gas and liquid superficial velocities were 1-16 and 0.7 cm-s-l, respectively. The temperature range was 25-96.56 “C. A gas holdup correlation implicating the effect of the operating temperature was developed with an average deviation of 3.1%. The result of this paper can be used for the scale-up and the design of bubble columns. Bubble columns are widely used in the fields of the chemical, biochemical, and food industries as well as environmental engineering as absorbers, desorbers, or reactors; see Schugerl et al. (1977), Jiang and Zhang (1981), Jiang (1983), Shah et al. (1982), and Zou (1985). Gas holdup is the basic parameter indicating the hydrodynamical characteristics of bubble columns. It affects directly the geometric sizes of bubble columns, and the gas-liquid interfacial areas thus affects the mass-transfer rates of bubble columns. So it is one of the necessary and important parameters for the design of bubble columns. Hence, many investigators used measuring methods to study the gas holdup in bubble columns to correlate the gas holdup, the physical properties of the system, and the sizes of the apparatus, as well as the operating conditions. On the basis of the experiments, many gas holdup correlations were suggested by different investigators. However, most of the previous correlations were developed based on the experimental data in bubble columns operated at ordinary temperatures, less than 40 “C, but the commercial bubble columns are commonly operated at elevated temperatures, often nearing the boiling points. It is wondered In accordance with the authors’ wishes, their family names are listed first. Zou Renjun is concurrently a t Hebei Institute of Technology and Northwestern University.
whether the previous correlations can be used for the scale-up of the commercial bubble columns. Although Quick and Deckwer (1981) measured the gas holdup in a bubble column operated at the temperature range 60-170 “ C for some liquid hydrocarbons, they didn’t developed a new correlation to correlate quantitatively their data. Jiang (1983) has suggested a gas holdup correlation. In the present paper, a complementary gas holdup correlation implicating the effect of the temperature in a bubble column operated at elevated temperatures is developed. Experimental Section The flow diagram of the bubble column is shown in Figure 1. The bubble column is made of stainless steel. It is 100 mm in inside diameter and 1.05m in height. Four pressure taps were drilled in the wall at 250-mm intervals, the lowest, tap 1, being 150 mm above the bottom plate of the column. The gas and liquid spargers are of single nozzle type, and the nozzle’s inside diameter is 10 mm. The gas and liquid nozzles are located 50 mm and 25 mm above the bottom plate of the column, respectively. Taps 1 and 3 are also used as the insertion points of the thermocouples. The column was operated continuously,and the gas and liquid flowed upward concurrently. The feed superficial velocities were 1-16 cms-’ with respect to the gas flow and
088~-58S5/S8/2627-1910$01.50/0 0 1988 American Chemical Society
