Ind. Eng. Chem. Res. 2008, 47, 9983–9995
9983
Analysis of Extractive Distillation with Mathematical Programming Abdulfatah M. Emhamed, Barbara Czuczai,* Endre Rev, and Zolta´n Lelkes Department of Chemical and EnVironmental Process Engineering, Budapest UniVersity of Technology and Economics, H-1111, Budapest, Budafoki ut 8, Budapest, Hungary
Extractive distillation processes are investigated by optimization. The proposed MINLP model and optimization tool are able to analyze complex distillation processes. Application of the method is demonstrated on two different extractive processes: acetone/methanol mixture using water as heavy solvent and ethanol/water mixture using methanol as light solvent. The two systems are analyzed separately; the analysis is performed both from technological and economical aspects. The optimal structures are found widely independent of the weights of different cost parts. Thus, small changes in the economic environment will not alter the optimal configuration. 1. Introduction Distillation is one of the most widely used techniques to separate liquid mixtures. The separation is based on the difference in boiling temperature of the components. By successively heating up and condensing the mixture, the lower boiling component is enriched in the vapor phase while the higher boiling component is enriched in the liquid phase. Not all liquid mixtures can be separated by ordinary fractional distillation, however. When the components to be separated constitute an azeotropic mixture, separation becomes impossible with conventional distillation. This work deals with extractive distillation, a distillation technique where another component, the so-called entrainer or solvent, is added to the system. Extractive distillation is widely used in the chemical industry and is commonly applied in a continuous multicolumn processes.1 A typical extractive distillation process configuration is shown in Figure 1. Design of a homogeneous azeotropic sequence performing the desired separation is usually carried out in two steps:2 (1) screening for a potential solvent, (2) synthesizing a separation sequence for the selected solvent. Separation of a binary mixture with extractive distillation is usually carried out in two columns. One of the azeotropic constituents is produced as a pure product in the first, extractive, column. The other product of that column contains the solvent and the other azeotropic constituent. The other azeotropic constituent is separated from the solvent in the second, solvent recovery, column. The solvent is usually recycled to the extractive column. One of the most important steps in developing a successful (economical) extractive distillation sequence is solvent selection. The solvent must be carefully chosen.3–5 Langston et al.,6 using HYSYS software, investigated the effect of solvents for separating several binary mixtures with extractive distillation. In another work Munoz et al.7 used HYSYS to compare extractive distillation and pressure swing distillation for the separation of isobutyl alcohol and isobutyl acetate using water as solvent. Hilal et al.8 studied the possibility of reducing the specific consumption of extractive solvent, to separate some binary systems, by increasing the distance between the mixture feed stage and extractive agent feed stage, and by splitting the extractive agent feed. Studying these distillation processes is important for deeply understanding them and successfully designing economically * To whom correspondence should be addressed. E-mail: ufo@ mail.bme.hu. Tel.: +36-1-463-1189. Fax: +36-1-463-3197.
Figure 1. Technological scheme of a typical extractive distillation process with solvent recovery.
optimal distillation systems. Such research is often performed using process simulators, but application of optimization tools is more convenient. When the problem is to find the optimal or minimum reflux ratio for a fixed configuration and product composition, the user of the simulator software has to make a lot of trials to determine the optimal solution. The same results can be achieved in one step using an NLP model, without making a lot of useless attempts. Similarly, the optimal structure can be calculated in one step with some iteration cycles by solving an MINLP problem. Many interesting and effective approaches have been developed in the past decade for optimizing distillation columns. These models simultaneously determine the configuration of the column and the operation parameters of the distillation process. As the number of stages is also to be determined, logic (or binary) variables are needed to model their existence. Some researchers apply GDP-models.9 Here logic variables are mapped to the stages, representing their existence. Some other authors apply MINLP models.10 In those cases binary variables are used for assigning the stage onto which the reflux is fed. When optimizing a complex distillation system, usually characterized with large problem size, the huge number of strongly nonlinear equations involves a serious difficulty in solving the model, and good initial values are needed for solving the NLP subproblems during the outer approximation algorithm. In a recent work by Farkas et al.11 the original outer approximation algorithm is modified in order to provide good initial values in each iteration. In the present work, an analysis method based on applying an optimization tool is developed for complex distillation
10.1021/ie071556z CCC: $40.75 2008 American Chemical Society Published on Web 11/17/2008
9984 Ind. Eng. Chem. Res., Vol. 47, No. 24, 2008 Table 2. Values of Factors at Levels
bpay Cs Cw UF
-2
-
center
+
+2
2.7874 1.6645 × 10-6 4.2867 × 10-8 1.4715
4 1.7667 × 10-6 4.3923 × 10-8 1.5657
7 2.01335 × 10-6 4.6472 × 10-8 1.7934
10 2.26 × 10-6 4.902 × 10-8 2.021
11.2426 2.3622 × 10-6 5.0076 × 10-8 2.1153
design limits and collecting preliminary data for optimal design and analysis from an economical point of view. The feasible region of the process is explored in order to find good bounds and initial estimates of the design and operation parameters. Using the initial data, MINLP is applied for determining the optimal structure and operational parameters in one step, in order to study the effect of cost factors on the optimal configuration and analyze the process from an economical point of view. The extractive distillation (in general) is found to be a stable process in the sense that the optimal structure is rather insensitive to the cost factors. Two different types of extractive distillation systems, namely acetone/methanol azeotropic mixture with water as heavy entrainer, and ethanol/water azeotropic mixture with methanol as light entrainer,12 are used as examples to demonstrate the application of the method. 2. Proposed Method
Figure 2. Superstructure of a single distillation column.
Figure 3. Visualization of the values of factors in the orthonormal matrix. Table 1. The Orthonormal Matrix test
bpay
Cs
Cw
UF
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25
-1 -1 -1 -1 -1 -1 -1 -1 +1 +1 +1 +1 +1 +1 +1 +1 -2 +2 0 0 0 0 0 0 0
-1 -1 -1 -1 +1 +1 +1 +1 -1 -1 -1 -1 +1 +1 +1 +1 0 0 -2 +2 0 0 0 0 0
-1 -1 +1 +1 -1 -1 +1 +1 -1 -1 +1 +1 -1 -1 +1 +1 0 0 0 0 -2 +2 0 0 0
-1 +1 -1 +1 -1 +1 -1 +1 -1 +1 -1 +1 -1 +1 -1 +1 0 0 0 0 0 0 -2 +2 0
processes. The analysis method is used for deeply understanding the studied processes through exploring their operational and
The main objective of this paper is to present an optimizationbased method for studying complex distillation processes from an economical (financial) point of view. The optimization is carried out using mixed-integer nonlinear programming. Because of the size and complexity of the problem, this requires a good (feasible) initial configuration and tight bounds for the number of stages. Therefore, preparatory to solving the MINLP model, the first step of the analysis is exploring the feasible region. After completing the analysis, the next phase is the optimization itself, where the optimal configuration is determined at different cost-factors combinations. The feasible region of the extractive column is determined by four parameters jointly: the three column section stage numbers (rectifying, extractive, and stripping sections) and the reflux ratio. This implies that the first step is to find the minimum number of stages in each column section with which the separation may be feasible. These minima can be used as lower bounds in subsequent optimization steps. The relationship between actual stage numbers in the column sections and minimum and maximum reflux ratios constitue the boundary of the feasible region. These data can be used for providing tight bounds and good, near-optimal, initial values for the reflux ratio in subsequent optimization steps. Additionally, they can be utilized to understand the behavior of these processes from a chemical engineering point of view.
Figure 4. Process scheme of acetone-methanol separation.
Ind. Eng. Chem. Res., Vol. 47, No. 24, 2008 9985
Figure 5. The total cost as function of stage numbers section by section. Table 4
Table 3. Run Results of the Experimental Design for Acetone/ Methanol System test
bpay
Cs
Cw
UF
total cost
∑QR
∑N
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25
-1 -1 -1 -1 -1 -1 -1 -1 +1 +1 +1 +1 +1 +1 +1 +1 -2 +2 0 0 0 0 0 0 0
-1 -1 -1 -1 +1 +1 +1 +1 -1 -1 -1 -1 +1 +1 +1 +1 0 0 -2 +2 0 0 0 0 0
-1 -1 +1 +1 -1 -1 +1 +1 -1 -1 +1 +1 -1 -1 +1 +1 0 0 0 0 -2 +2 0 0 0
-1 +1 -1 +1 -1 +1 -1 +1 -1 +1 -1 +1 -1 +1 -1 +1 0 0 0 0 0 0 -2 +2 0
249.823 281.8288 212.0692 277.4572 253.6883 285.3449 238.3446 295.5067 129.9131 155.1512 130.0188 152.6864 151.5152 172.2935 148.3437 166.6504 346.767 147.346 162.036 171.702 175.06 188.364 165.273 200.174 187.691
12.04 12.044 9.613 10.793 12.036 11.066 10.933 10.723 10.665 10.672 10.721 11.596 10.663 10.786 10.886 10.79 10.668 12.029 9.202 8.835 10.845 10.791 11.596 10.791 12.032
63 63 70 69 63 64 64 70 71 79 70 64 79 77 69 69 71 63 81 74 67 78 64 69 63
Errors in estimating possible future cost factors can affect the efficiency of optimization. If the price of the energy is estimated inaccurately, the calculated optimum may be far from the real future optimum. How the different cost-factors will affect the optimal configuration is investigated with MINLP in the second, final, step using the previously determined bounds and initial data. As the effect of cost factors is investigated through solving mathematical programming problems, it seems to be reasonable to apply mathematical programming also in the phase of exploration of the feasible region, in the preparatory first step. Therefore, instead of using the usual process simulators, the minimal number of stages and minimal and maximal reflux ratios are also determined through optimization. 2.1. Model and Software. The MINLP model developed by Farkas et al.11 is used for analyzing extractive distillation. This model works with minimum number of binary variables. This minimum number is achieved by applying different types of units containing different equilibrium stage numbers in the superstructure of a single column (see Figure 2). If the binary variable of a unit equals one then the unit is included in the
bpay
Cs
Cw
(a) Coefficients of the Cost Factors to Total Cost b′0 ) 201.802 -58.476 6.821 -1.983 bj bjj 31.603 -8.491 -1.069 (b) Coefficients of Interactions to Total Cost bpay 1.208 2.258 Cs 1.906 Cw UF (c) Coefficients of the Cost Factors to ∑N b′0 ) 69.36 2.0343 -0.195 0.0778 bj bjj -1.9 3.35 0.85 (d) Coefficients of Interactions to ∑N bpay 0.87 -3.375 Cs -0.5 Cw UF (e) Coefficients of the Cost Factors to ∑QR b′0 ) 10.91 -0.0272 -0.039 -0.205 bj bjj 0.403 -0.762 0.12 (f) Coefficients of Interactions to ∑QR bpay -0.049 0.396 Cs 0.0924 Cw UF
UF
16.128 -0.564 -6.201 -1.087 3.365
0.654 -2.15 -0.375 0.25 -0.5
-0.0113 0.326 0.0566 -0.2012 0.162
structure, all the equilibrium stages contained in that unit work, and their input and output streams may take positive value. If the binary variable of a unit equals zero then the unit is not included in the structure, none of the equilibrium stages contained in that unit works, and their input and output streams must equal zero. In the latter case, the liquid and vapor streams flow through transmitter units instead of the units containing equilibrium stages. The superstructures of the extractive distillation processes detailed below are built up using the above summarized concept. As the number of nonlinear equations increase with the numbers of components, columns, and stages, it becomes more and more difficult to solve the model with the optimization tools. To facilitate the solution, a modification of the original outer approximation algorithm is suggested and applied by Farkas et al.11 In each iteration during the solution process, initial values are calculated from stage to stage in the column for all the variables in the model. The initial values are in accordance with and calculated from the array of binary variables. By modifying the outer approximation algorithm in this way, it has become
9986 Ind. Eng. Chem. Res., Vol. 47, No. 24, 2008
Figure 6. (a) The total cost as function of bpay and Cs; (b) total cost as function of bpay and UF; (c) total cost as function of Cs and UF; (d) ∑N as function of bpay and Cs; (e) ∑N as function of bpay and UF; (f) ∑N as function of Cs and UF; (g) ∑QR as function of bpay and Cs; (h) ∑QR as function of bpay and UF; (i) ∑QR as function of Cs and UF.
possible to reliably solve the model and to use the method as an analyzing tool for studying complex distillation processes. The step of calculating initial values is also performed when only NLP problems are solved, namely, when the minimum reflux ratio is computed for fixed stage numbers. AIMMS 3.713 is used for implementing the modified algorithm. 2.2. Exploration of the Feasible Region. As it has already been shown,2,14 the feasible region of extractive distillation is different from that of conventional distillation. Whereas separation of the components and purity of the products increase monotonously with both the stage number and the reflux ratio in conventional distillation processes, extractive distillation can behave in a different way. Both the stage numbers2 and reflux ratio2,14 of the extractive column may have both lower and upper bounds.
Determining these values can be difficult, and can be more easily accomplished by optimization than with a trial and error procedure using process simulators. A lot of trials may be necessary in the latter, see, for example, the iterative method given by Chadda, Malone, and Doherty.15 Minimum and maximum reflux ratios for a certain number of stages can be found in one step using NLP. The minimum and maximum stage numbers can be found by systematic search in less time and with less trouble than with use of simulators. For this aim, each section of the extractive column is allowed to have maximum 31 trays in the present paper. A 3-dimensional grid covering the [0, 31] × [0, 31] × [0, 31] domain of the stage number combinations is created. When determining the borders of the feasible region, the section stage numbers are fixed as parameters, the binary variables are eliminated from the model, and
Ind. Eng. Chem. Res., Vol. 47, No. 24, 2008 9987
Figure 7. (a) NR as function of bpay and UF; (b) NR as function of Cs and UF; (c) NR as function of bpay and Cs; (d) NE as function of Cs and UF; (e) NE as function of bpay and UF; (f) NE as function of bpay and Cs; (g) NS as function of Cs and UF; (h) NS as function of bpay and UF; (i) NS as function of bpay and Cs.
NLP problems are solved. Instead of developing partial models for determining the feasible region of the extractive column, the model of the full process can be applied with fixed stage numbers and the reflux ratio in the conventional column. The resulting NLP problems are solved in AIMSS 3.7, using CONOPT 3.14G as NLP solver. After finding the minimum stage number in each section, the relationship between the sections stage numbers (as independent variables) and minimum and maximum reflux ratios (as dependent variables) is explored. These five data altogether determine the feasible region. How these data depend on each other, that is, interaction between sections, is also studied.
2.3. Studying the Effect of Cost Factors on the Optimal Configuration. As a result of the above detailed study, good initial values and bounds are available that can be utilized to investigate the effect of changing cost factors. Here the total process is optimized, that is, the operation parameters (reflux ratios, reboiler duties, etc.) as well as the stage numbers in the columns are considered unknown, and thus none of the binary variables are fixed. This implies that an MINLP is solved because the model contains both binary and continuous variables and nonlinear equations. Unlike the exploration of the feasible region, investigation of the effect of cost factors cannot be done using process simulators like
9988 Ind. Eng. Chem. Res., Vol. 47, No. 24, 2008
Figure 8. Ethanol-water separation applying methanol solvent.
ASPEN. As the optimum structure needs to be determined at different parameter settings, this step requires the use of optimization tools, and solving MINLP-s. Solving large MINLP problems with huge number of nonlinear (and nonconvex) constraints may require good bounds and initial values. The results of the preparatory step exploring the feasible region can provide us with tight bounds for both the continuous and the discrete variables. The 3-dimensional grid described in section 2.2 results in a number of NLP problems, some of which are infeasible because of the insufficient number of stages. These infeasibilities show the lower bound of the number of stages of each column section. The minimum of the reflux ratios of all the feasible runs provides the lower bound, the maximum provides the upper bound of the reflux ratio. Those tight bounds help the solver to find the solution over a reduced search space. Guaranteeing tight bounds described in section 2.2 might have been carried out using engineering criteria or heuristics. Using the same optimization tool for both task, however, is an automated and more convenient way. The cost function given by ref16 for calculating the total annual cost of distillation columns is used as objective function to be minimized. The form of the cost function is the following: total cost ) btaxCs∆HRV˙ + btaxCw∆HCV˙ + UF
f(N, D) (1) bpay
f(N, D) ) 12.3[615 + 324D2 + 486(6 + 0.76N)D] + 245N(0.7 + 1.5D2) (2) where Cs is unit price of steam, ∆HR˙V is reboiler duty, Cw is unit price of cooling water, UF is update factor of the price of column installation, bpay is payback period, N is number of stages, and D is column diameter. The four cost factors (Cs, Cw, UF, bpay) are chosen to investigate their effect. A steady increase of the objective function with increasing Cs, Cw, UF, and decreasing bpay is predictable from its mathematical form. Our goal here is to explore which factors have significant effect on the total annual cost and the optimal structure and how the optimal structure changes with these cost factors. An orthonormal matrix containing the set values of the cost factors is designed for this aim, based on the composite design used for multiple regression in the statistical literature.17 Two basic, reasonably estimated, levels (two different values) are set to each cost factor. The process is optimized at each possible combination of these values; this involves 24 experiments. To determine the nonlinear parts in the regression polynomial, an additional experiment is made at the center point of the design and at excessive points (star-points) outside the core of the design, as is visualized in Figure 3.
The orthonormal matrix is presented in Table 1. The columns of the table represent the four cost factors; the “-1” and “+1” signs denote the lower and upper levels of the factors; “0” marks the center point; “-2” and “+2” are used for the star-points. The cost factor values at the two levels and in the center point are shown in Table 2. Table 2 assigns 25 MINLP problems. These are solved in AIMMS 3.7 with modified outer approximation algorithm. CPLEX 10.0 is applied as MILP solver; CONOPT 3.14A is used for NLP subproblems. The total cost calculated according to eq 1 is used as objective function. The value of total cost, the sum of the stage numbers, and the total reboiler duty, are the recorded measurement values. After determining the three measurement values at each design point, multivariate nonlinear regression is applied to calculate the parameters of the following polynomial: Yˆi ) b0 +
∑ b x + ∑ b x x + ∑ b (x j ji
j
jk ji ki
j
