Phase-Plane Analysis of Homogeneous Azeotropic Distillation Columns

Changes in holdup are compared to the global mass balance to derive a vector field describing the trajectories of the product compositions. The vector...
0 downloads 0 Views 207KB Size
Ind. Eng. Chem. Res. 2002, 41, 3963-3969

3963

Phase-Plane Analysis of Homogeneous Azeotropic Distillation Columns Eleonora Bonanomi and Manfred Morari* Automatic Control Laboratory Swiss Federal Institute of Technology, ETHZ, CH-8092 Zu¨ rich, Switzerland

In this paper, it is shown that despite the large number of states and equations describing a distillation column the dynamics can be analyzed by lumping the states into a few variables characterizing the total holdup. Changes in holdup are compared to the global mass balance to derive a vector field describing the trajectories of the product compositions. The vector field contains detailed information about the dynamic behavior and is obtained without performing a complete dynamic simulation. Here, the procedure is applied to homogeneous azeotropic distillation, but with appropriate adjustments, it should be applicable to any column. 1. Introduction Among the classical unit operations, there are many processes that exhibit multiple steady states (MSS); the continuous stirred tank reactor (CSTR) is the best known example.1 In the classical CSTR analysis, MSSs are explained by the balance of heat fluxes. In particular, the heat produced by a chemical reaction is removed by a coolant. While the removal of heat depends linearly on the reactor temperature, heat generation depends nonlinearly on the reactor temperature. The two fluxes are equal at the operating points. Therefore, representing each heat flux as a function of the reactor temperature, the intersections of the two curves correspond to steady states. If there are multiple intersections, MSS are present for certain operating conditions and the socalled “slope condition” can be applied as a necessary condition for stability. The CSTR with the “slope condition” is an example of a simple representation of the physical aspects leading to MSS and instability. In this paper, an analogous approach for distillation columns is proposed. Distillation columns, in contrast to a CSTR, are characterized by a large number of states, and the dynamics is described by hundreds of equations. As a consequence, a description similar to the CSTR is not easily developed and a reduction of the number of variables is necessary. Significant progress in understanding the existence of multiplicity in distillation columns has been made with the ∞/∞ analysis.2 This technique needs only basic (experimental) vaporliquid equilibrium information about the mixture and gives necessary and sufficient conditions for the existence of MSS in a column with an infinite number of stages and operating at infinite reflux. The development of appropriate software allowed to perform of a complete bifurcation analysis of the system using a detailed model of the column; for a given column, the steady states were identified and their stability was investigated.2 Another step was made by Cornelius Dorn,3 who analyzed how the column holdup is changing to derive insight into the dynamic behavior of the column. He developed a methodology to qualitatively characterize the stability conditions of steady-state profiles. Starting

from a perturbation of the steady-state profile, he studied how its shape is affected by changes in the total holdup composition. This allowed him to determine if the profiles are locally stable or unstable. In this paper his idea is extended in the sense that the qualitative interpretation is now quantified. The column dynamics is described with a reduced number of variables. Using this simple description, the dynamics and phenomena leading to instability can be understood better. During the transient, the holdup of each component is changing. The accumulation or depletion of a component pushes the column profile in a direction that corresponds to profiles with a higher or lower quantity of that component in the column. The transient proceeds until the driving force is zero. The driving force for changes of the total column holdup is calculated as the difference between the actual product composition and its steady state. From the driving force, the vector field providing a global view of the behavior of the column for any distillate composition is available. In particular, the direction and the magnitude of movement of the column composition after a perturbation are calculated and represented in the distillate composition space without performing a dynamic simulation. Therefore, the analysis is not local but gives results for every physically possible distillate composition. As shown later, MSS and instability are explained in a way similar to that for the CSTR. The analysis is applied to a homogeneous azeotropic distillation column, but with minor adjustments, it is suitable for other kinds of columns. In the next section, general features regarding the behavior of columns during the transient are introduced. These observations constitute the basis of and define the limitations of the proposed analysis. Section 3 presents the procedure to determine the number of steady states and their stability for a given column. Section 4 introduces the vector field representing the direction of movement of the distillate composition after a perturbation. 2. Column Profiles during Transient

* Author to whom correspondence should be addressed. Phone: +41 1 632-2271. Fax: +41 1 632-1211. E-mail: [email protected].

As pointed out by Skogestad and Morari4 and later confirmed by Kumar and Daoutidis,5 the transient of

10.1021/ie010727b CCC: $22.00 © 2002 American Chemical Society Published on Web 07/16/2002

3964

Ind. Eng. Chem. Res., Vol. 41, No. 16, 2002

columns occurs in two different time scales. The fast dynamics is related to the internal distribution of the holdup in the column, while the slow dynamics corresponds to the movement of the distillate and bottoms composition to reach positions that satisfy the material balance. Moreover, starting from the material balance equation describing the dynamics of the column and letting the reflux go to infinity, one obtains the equation for the distillation line.6 Therefore, in a column operating at infinite reflux, the profile coincides with a distillation line at any time during transient. The consequence of this observation is the following. Assume a column with an initial profile and output compositions that do not correspond to steady state. For any initial column profile, the holdup in the column first quickly readjusts to coincide with a distillation line and then the distillate and bottoms move to eventually satisfy the material balance. During the movement of the bottoms and distillate, the profile maintains the property of coinciding with a part of a distillation line. If, for any reason, there is a displacement from a distillation line, the time necessary for the profile to coincide again with a distillation line is negligible. In the qualitative analysis, Dorn3 uses the distillation lines. Considering the similarity of these lines with residue curves, one can use residue curve maps and obtain a geometrical tool to investigate the stability of a column as a complement of the ∞/∞ analysis. In the following, the column profile is calculated using a stage by stage procedure: the distillation line equation yk+1 ) xk is substituted with the mass balance equation around the top of the column, yk+1 ) (1/V)(Lxki + DxD i ), and the profile is calculated going down stage by stage starting from the top of the column. Notice that, in the limit for infinite reflux, the mass balance equation is again the distillation line equation. Using the stage by stage procedure, the contribution of the external fluxes is not neglected and the analysis is suitable also for a column operating at high but finite reflux. 3. Stability Analysis 3.1. Nullclines. Consider a distillation column model based on the constant molar overflow assumption. The global mass balance for each component of the mixture is

d dt

n

(

Mkxki ) ) FxFi - Bx˜ Bi - DxD ∑ i k)1

(1)

Here x˜ Bi is the bottoms composition calculated for a specified xD i with the procedure presented in appendix A. This procedure determines the profile and bottoms composition based on the discussions of the previous section: at any time the profile along the column coincides with a part of a distillation line. Therefore, for a given distillate and proceeding down stage by stage, the composition on each stage is calculated. When ˆ Bi is the steady-state mass balance FxFi - DxD i ) Bx substituted in eq 1, one obtains

d dt

n

(

d

ˆ Bi - x˜ Bi ) ∑ Mkxki ) ) Mtot dt(xavg i ) ) B(x k)1

(2)

If the holdup on each stage is the same, a new variable xavg is introduced and calculated as follows: i

n

xavg ) i

∑ xki k)1 n

(3)

It reflects the total holdup of component i in the column.6 The dynamics of n states, xki where k ) 1, ..., n, related to component i, are lumped into a new variable. Moreover, the right-hand side shows clearly the driving force during the transient: xˆ Bi - x˜ Bi . This term can be interpreted as the displacement of the output composition from steady state, and it indicates if a component is accumulating or depleting in the ) 0 holds for each column. At steady state, (d/dt)xavg i component, x˜ Bi ) xˆ Bi , and the driving force is zero. The ) 0 for each component can be nullclines (d/dt)xavg i calculated and plotted in the distillate composition space. Their intersections correspond to steady states. In addition, they divide the space around the steadystate distillate composition into six regions corresponding to the regions used for the qualitative analysis.6 When a disturbance moves the distillate to one of these regions, then each component accumulates or depletes in the column and a sector of directions of possible movements of the distillate composition can be calculated. To simplify the presentation, each step of the stability analysis is explained in detail in the next section using an example. Before that, it is worth noticing the analogy with the CSTR. For the column, as for the CSTR, the multiplicity is explained using the balance between two quantities. The energy balance of the CSTR is here the material balance. The exchange of energy with the external source, the linear term, is here the global mass balance, and the internal production of energy, the nonlinear term, is substituted by the constraint on the profile given by the distillation line and representing the internal distribution in the column. In the CSTR problem, the temperature is used as an independent variable, while the distillate composition is chosen for the column. Finally, the dependent variable is the energy for the CSTR and the bottoms composition for the column. 3.2. Direction of Movement. The mixture considered (acetone, benzene, and n-heptane) is of the class 001.7 The minimum-boiling binary azeotrope between acetone and n-heptane is at 93.56% (molar basis) of acetone. The column has 21 stages. The feed of 100 kmol/h consists of 92% acetone, 7.93% benzene, and 0.07% n-heptane and enters on stage 15. A high reflux (reflux/feed ) 500) is used (the method is valid for any “sufficiently high” reflux). The distillate flow is 0.6 kmol/ h. In Figure 1a, the surfaces xˆ Bi and x˜ Bi for acetone are plotted as a function of the distillate composition xD. For each specific distillate xD, the distance between the two surfaces is proportional to the driving force determining the transient of the system. The intersection of the two surfaces corresponds to a zero driving force for the accumulation of acetone and is the nullcline for this component. In the same way (Figure 1b,c), the nullclines for the other two components can be calculated. The three nullclines for acetone, benzene, and heptane are projected in the distillate composition space in Figure 1d; they intersect at the steady state. It is worth noticing here the meaning of the two surfaces. The xˆ Bi surface is a plane corresponding to the linear equation given by

Ind. Eng. Chem. Res., Vol. 41, No. 16, 2002 3965

Figure 1. (a-c) Bottoms molar fraction (x˜ Bi and xˆ Bi ) as a function of the distillate composition. (d) The projection of the nullclines in the distillate composition space identifying six regions.

the steady-state global mass balance. The x˜ Bi surface is not a plane and reflects the nonlinear connection between distillate and bottoms compositions. As for the CSTR, the interaction of these two surfaces determines the presence of MSS and their stability as explained in the following. As Figure 1d shows, each nullcline divides the distillate composition space into two regions. The region where x˜ Bi > xˆ Bi corresponds to an accumulation of component i in the total column holdup when xD is in that region, and the other region corresponds to depletion. This comparison is repeated for each component in each region yielding a complete picture of the change in holdup (cf. Table 1). As the holdup is changing, the shape of the column profile and the distillate and bottoms compositions will adjust according to this change. The dynamics of the holdup describes the dynamics of the whole column and is represented by the average composition in the column. For each distinct distillate, there are three average mole fractions, one for each component, given by ∑kxki k ) nxavg i , where xi is calculated by means of the stage by

Table 1. Comparison of x˜ B ˆB i and x i for Each Component in Each Regiona region

acetone

benzene

heptane

1 2 3 4 5 6

accumulation depletion depletion depletion accumulation accumulation

accumulation accumulation accumulation depletion depletion depletion

depletion depletion accumulation accumulation accumulation depletion

a Depletion in column holdup: x ˜ Bi > xˆ Bi . Accumulation in column holdup: x˜ Bi < xˆ Bi .

stage procedure. In other words, the complete dynamics is lumped into the dynamics of the three average molar fractions. Thus, the analysis is reduced to a lower dimension system. Using avg i , the interpretation of what is happening in the column is much easier because there is no need to follow the evolution of the composition on all stages. Consider, for example, a perturbation that brings the distillate composition to region 5. Then there is an accumulation of heptane and acetone and a depletion of benzene. The accumulation of heptane implies move-

3966

Ind. Eng. Chem. Res., Vol. 41, No. 16, 2002

Figure 2. (a-c) Contour lines of constant holdup. The semicircles represent the possible direction of movement of the distillate. (d) Combination of conditions in parts a-c yielding the shown sector. The distillate moves in a direction lying within this sector.

ment of the distillate composition toward higher holdup for heptane. As mentioned before, the holdup is represented by xavg. The contour lines of constant holdup composition as a function of the distillate are plotted in Figure 2. The holdup of each component is constant along the corresponding contour lines, and it increases or decreases when moving from one contour line to another. Based on these contour lines, it is possible to identify the cone of the feasible directions of movement. All feasible directions for increasing the holdup of heptane lie in the semicircle “H” plotted in Figure 2c. The accumulation of acetone implies the movement of the distillate in a direction that must lie in the semicircle “L”. Finally, the depletion of benzene brings the column in the direction of lower benzene holdup; the direction is in the semicircle “I”. These three conditions together restrict the trajectories for the distillate to fall within the cone highlighted in Figure 2d. The same analysis is repeated for all of the six regions. In Figure 3a, the six sectors and the steadystate distillate composition are illustrated. The arrows show the direction of the distillate in each sector: the

Figure 4. Projection of the nullclines in the distillate composition space for a column with three steady states corresponding to points 1-3.

column would move back to the steady state, which is therefore stable. The result is confirmed in Figure 3, representing distillate trajectories calculated with a dynamic simulation. Coming back again to the analogy with the CSTR, these results can be interpreted in the following way. If the distillate is perturbed, the profile and the bottoms will readjust quickly in order to belong to the same distillation line. The mass balance is not satisfied anymore, and the different components accumulate or deplete in the column. For a stable steady state, the accumulation and depletion of each component direct the column toward the original equilibrium point. For an unstable steady state, the directions point away from the original point toward a stable point. In the same way, the analysis can be applied to other cases where MSS exist. One case with MSS is shown in Figure 4. The column has 46 stages; the feed of 100 kmol/h consists of 89.37% acetone, 0.7% benzene, and 9.93% n-heptane and enters on stage 41. The reflux is 4800 kmol/h, and the distillate flow rate is 7.5 kmol/h. The nullclines intersect in three points corresponding to three steady states. The stability of each equilibrium point can be studied in the same way as has been done

Figure 3. (a) Predicted directions of movement of the distillate after perturbation. (b) Dynamic simulation for different perturbations of the column in the neighborhood of a stable steady state.

Ind. Eng. Chem. Res., Vol. 41, No. 16, 2002 3967

Figure 5. Vector field for two different columns: (a) one steady state; (b) limit cycles. The arrows represent the direction and the magnitude of the movement.

for the previous column. For a stable steady state, the feasible directions point toward the steady state; for unstable points, the feasible direction points away from the steady state. Figure 4 shows the nullclines for this system, and the six regions around each steady state can be identified. For point 1, they correspond exactly to the qualitative regions,6 but differences are present for points 2 and 3. For these two points, the nullclines are not parallel to the three edges of the composition space, but they assume a more complex shape. This is because in the quantitative analysis the bottoms is not fixed as it was in the qualitative study. The movement of the bottoms resulting from the distillate movement cannot be neglected in a finite column. In fact, it is the movement of the bottoms that determines the nonlinear shape of the nullclines with the connected multiple intersections. Identifying the nullclines and their intersections allows us to calculate the position of the steady states. A first impression about the stability can also be derived. However, in some cases the cone of the possible directions is very large. Therefore, exact conclusions about stability cannot be extracted. In the next section, we will go one step further: the vector field representing the trajectories of the distillate during the transient is calculated. 4. Vector Field In this section, the results obtained previously are used to construct the vector field: for every distillate composition in the distillate space, the vector representing the direction and magnitude of movement is calculated. Assume a distillate composition xD that does not correspond to a steady state. In the previous section, two values of the bottoms composition have been calculated: x˜ B given by the stage by stage procedure and xˆ B that satisfies the steady-state global mass balance. Then, when the two values were compared, the nullclines and the steady states were identified. The equation representing the dynamics of the total column holdup of each component is again considered:

B d avg x ) tot(xˆ Bi - x˜ Bi ) dt i M

(4)

Because the right-hand side is known, the incremental change of the total column holdup in a unit of time, ∆xavg, can be computed for each distillate composition xD. Among the possible states to represent the column, and therefore the vector field, the distillate is the most significant, also considering that it has been chosen as an independent variable. The vector fields in Figure 6 are plotted in the distillate space; therefore, for every ∆xavg/∆t, the corresponding ∆xD/∆t ) (x′D - xD)/∆t must are be computed. The initial conditions xD and xavg i known. As explained before, using eq 4, the changes of the total column holdup, ∆xavg i , are evaluated in a unit of time for each of the three components. The column will move to a point x′D that corresponds to the new total column holdup. To determine x′D, the contour line of the is calculated for each component by interpolanew xavg i tion. The intersection of the three contour lines gives the point x′D, and the vector x′D - xD identifies the direction and magnitude of the vector field for the initial point xD. Applying the procedure described above for a set of xD forming an appropriate grid in the distillate composition space, a complete vector field is constructed (Figure 5). (In complete analogy, a vector field for the bottoms composition can be formulated; see appendix C.) Because the vector field at a point is tangent to the trajectory through that point, trajectories for any given distillate composition can be sketched. When the vector field is studied, equilibria and their stability properties as well as limit cycles can easily be determined. Figure 6 represents the same columns of Figure 5 considering the direction but not the magnitude of the movement. The continuous lines are trajectories evaluated with a dynamical simulator. Comparing the arrows with the lines shows a good agreement between the vector field and the evolution of the distillation unit. 5. Conclusion A new approach for studying the dynamic behavior and stability of distillation columns was presented. Following the qualitative interpretation of the column dynamics introduced by Dorn,3 a quantitative approach was proposed here. The only assumption underlying the analysis is that the reflux flow rate must be sufficiently high compared to the feed flow rate. Under this assumption the column profile during the transient must

3968

Ind. Eng. Chem. Res., Vol. 41, No. 16, 2002

Figure 6. Comparison between the vector field and dynamic simulation for two different columns: (a) one steady state; (b) limit cycles.

Figure 7. Flows entering and leaving different sections of a distillation column: condenser, second stage, and reboiler.

coincide with a part of a distillation line. The interaction of this nonlinear constraint and the steady-state material balance in determining the column holdup is the basic phenomenon behind MSS. The mismatch between the transient output compositions and the steady-state output compositions is the driving force of the system. Considering this and the fact that the column transient can be represented by the average composition in the column, the vector field was evaluated. It represents the direction of movement of the column after a perturbation, and it gives a complete picture of the column behavior without the need for numerous dynamic simulations. Especially in the presence of MSS and oscillations, this is a simple tool to follow the transient of the column. In this paper, the method was applied to a homogeneous azeotropic distillation column, but it can be extended to heterogeneous and reactive columns. Appendix A: “Stage by Stage” Procedure Given is a column with n stages and feed (F, xF). The steady-state column profile can be determined through the classical McCabe-Thiele method8 for any feasible distillate and reflux flow rate. Starting with a first guess for distillate composition xD, the complete composition profile is calculated. Then the procedure is repeated with a new xD until the global mass balance is satisfied. If the procedure is not repeated, for each xD a bottom x˜ B is obtained based on the following assumption. Considering a column with infinite reflux and the two dynamics discussed in section 2, at the end of the fast

dynamics the column profile can be approximated by a distillation line. If the reflux is finite, the external fluxes cannot be neglected anymore, and instead of the distillation line equation, the material balance, eq 6, and the vapor-liquid equilibrium are used. In this condition, the column does not fulfill the global mass balance but will move, during the slow dynamics, to fulfill it as specified in section 2. In particular, x˜ B is calculated using the following procedure. First, the vapor and liquid flows in each section of the column are determined using mass balances. Then, the compositions on each stage are calculated. For a total condenser (Figure 7a), the compositions y2 and x1 and the distillate xD are the same. Figure 7b shows the flows for the first tray, stage 2. Here, the liquid composition x2 and the temperature on the tray are evaluated through a dew-point calculation for the known vapor composition y2. The composition y3 of the vapor rising from the tray below is then determined by the mass balances around the top of the column:

y3i )

1 (Lrectx2i + DxD i ) Vrect

(5)

The procedure is repeated for each stage of the rectifying section. For the feed tray (stage t) and all stages in the stripping section, the mass balance changes to

) yt+1 i

1 F (Lstripxti + DxD i - Fxi ) Vstrip

(6)

Ind. Eng. Chem. Res., Vol. 41, No. 16, 2002 3969

Concerning the partial reboiler (Figure 7c), yn and xn-1 are known from the previous calculations. Therefore, the bottoms composition x˜ B is computed with a dewpoint calculation for the known composition yn. Appendix B: Bottoms-Distillate Composition Space The described procedure is based on fixing the distillate composition and then evaluating the bottoms composition. The nullclines, the six regions, and the steady states are represented in the distillate composition space. This concept can be applied in a similar way starting with a fixed bottoms composition and evaluating the distillate composition. In this case, the three nullclines divide the neighborhood of the bottoms steady state into six regions. These lines can be directly obtained from the one previously calculated by applying the mass balance equations. For instance, the acetone nullcline is projected onto the bottoms space using the acetone material balance for each point of the line:

xBac )

FxFac - DxD ac B

(7)

The same is valid also for the other two components of the mixture. In the bottoms space, there are six regions and they correspond directly to the regions in the distillate space. The related regions are characterized by accumulation and depletion of the same components. The contour lines representing the average composition in the column are analogous. Consequently, the vector field for the bottoms composition space can be represented based on the same arguments as those applied for the distillate composition space. Appendix C: Simulation Details The thermodynamics parameters for the acetonebenzene-n-heptane mixture are reported in Gu¨ttinger

and Morari.9 In all of the simulations, the Antoine equation is applied to calculate the vapor pressure and the Wilson model is used for the liquid activity coefficient. The dynamic simulations are based on the CMO model integrated with DDASAC,10 an index-1 DAE solver. The column layout is specified in section 3.2, and the liquid holdup for the reboiler and condenser and on each tray is fixed to 3 kmol. Literature Cited (1) Van Heerden, C. Authothermic Processes: Properties and Reactor Design. Ind. Eng. Chem. 1953, 45 1242-1247. (2) Bekiaris, N.; Meski, G.; Radu, C.; Morari, M. Multiple Steady States in Homogeneous Azeotropic Distillation. Ind. Eng. Chem. Res. 1993, 32 (9), 2023-2038. (3) Dorn, C. Nonlinear Dynamics in Homogeneous Azeotropic Distillation. Ph.D. Thesis, ETH No. 13947, ETH, Zu¨rich, Switzerland, 2000; available via http://www.control.ethz.ch/∼disti. (4) Skogestad, S.; Morari, M. Understanding the Dynamic Behavior of Distillation Columns. Ind. Eng. Chem. Res. 1988, 27 (10), 1848-1862. (5) Kumar, A.; Daoutidis, P. Nonlinear Model Reduction and Control of High-Purity Distillation Columns. In 1999 American Control Conference, IEEE Service Center: Piscataway, NJ, 1999; Vol. 3, pp 2057-2061. (6) Dorn, C.; Morari, M. Qualitative Analysis for Homogeneous Azeotropic Distillation. 1. Local Stability. Ind. Eng. Chem. Res. 2002, 41, 3930-3942. (7) Matsuyama, H.; Nishimura, H. Topological and Thermodynamic Classification of Ternary Vapor-Liquid Equilibria. J. Chem. Eng. Jpn. 1977, 10 (3), 181-187. (8) Stichlmair, J.; Fair, J. Distillation: Principles and Practice; Wiley: New York, 1998. (9) Gu¨ttinger, T.; Morari, M. Comments on “Multiple Steady States in Homogeneous Azeotropic Distillation”. Ind. Eng. Chem. Res. 1996, 35 (8), 2816. (10) Stewart, W.; Caracotsios, M.; Sorensen, J. DDASAC Software Package Documentation, 1995.

Received for review September 4, 2001 Revised manuscript received April 23, 2002 Accepted May 15, 2002 IE010727B