I n d . Eng. Chem. Res. 1988,27,1466-1473
1466
Ohnishi, K.; Hori, Y.; Tomari, Y. Bunseki Kagaku 1977,26, 74. Sakai, A.; Kamaishi, T.; Kinoshita, K. J. Tech. Assoc. Pulp Paper Znd. 1979,62,45. Schwochau, K.; Astheimer, L.; Schenk, H. J.; Witte, E. G. 2.Naturforsch. 1982,37b, 214. Tabushi, I.; Kobuke, Y. Israel J. Chem. 1985,25,217. Tabushi, I.; Kobuke, Y.; Nishiya, T. Tetrahedron Lett. 1979a,3515. Tabushi, I.; Kobuke, Y.; Nishiya, T. Nature (London) 197913,280, 665.
Tabushi, I.; Kobuke, Y.; Ando, K.; Kishimoto, M.; Ohara, E. J . Am. Chem. SOC.1980,102,5947. Tabushi, I.; Kobuke, Y.; Nakayama, N.; Aoki, T.; Yoshizawa, A. Ind. Eng. Chem. Prod. Res. Dev. 1984,23,445. Tabushi, I.; Yoshizawa, A.; Mizuno, H. J. Am. Chem. SOC.1985,107, 4585.
Receiued for review November 3, 1987 Accepted March 11, 1988
Simulation of Continuous-Contact Separation Processes: Multicomponent Batch Distillation David M. Hitcht Department of Chemical Engineering, North Carolina State University, Raleigh, North Carolina 27695-7905
Ronald
W.Rousseau*
School of Chemical Engineering, Georgia Institute of Technology, Atlanta, Georgia 30332-0100
The present work describes the use of a numerical procedure for simulation of batch multicomponent distillation in a continuous-contact (packed-column) operation. The method uses a relaxation procedure to simulate both start-up and product withdrawal periods of column operation. Example calculations on a three-component mixture are used to show the effects of varying packing height, boilup rate, reflux ratio, and condenser holdup and to demonstrate that the computational procedure is rapid and stable. Additionally, use of the method in determining optimum operating strategy is illustrated. Although a substantial amount of work has been done recently on numerical methods for packed-column design and simulation, both the quantity and quality of available packed-column algorithms still rate far below the levels achieved in staged-columnsimulation. Deficiencies remain evident in several key areas. For example, one major problem that has yet to be addressed adequately is the inability of current numerical schemes to solve, without the threat of encountering convergence difficulties, the complicated systems of nonlinear model equations required for a rigorous treatment of simultaneous mass- and heat-transfer phenomena. Another significant deficiency is the lack of flexibility of most current packed-column algorithms; Le., they are typically capable of performing calculations on only one type of separation process under a relatively narrow range of operating conditions. Some other aspects of packed-column simulation that have not been considered in a rigorous way include (1)unsteadystate operations of continuous processes and (2) batch processes such as batch distillation. Considering the problem of nonconvergence to be of paramount importance, it was decided that the primary goal of the work presented here would be to develop a highly stable yet rigorous numerical procedure for packed-column simulation. On the basis of the experience gained from the development of staged column algorithms, a relaxation approach was chosen because of its demonstrated stability of convergence for any physically realizable process. This stability is obtained because a relaxation algorithm closely models the transient behavior of the process from startup to steady-state, i.e., if a steady-state condition is ever achieved. A procedure of this type has been used by von Stockar and Wilke (1977a,b) to converge on the steady-state operations of a simple, three-compo-
* To whom correspondence
should be addressed. Current address: Tennessee Eastman Company, P.O. Box 1972, Kingsport, T N 37662. 0888-5885/88/2627-1466$01.50/0
nent, packed-column absorption process. Their algorithm was claimed to have “superior convergence properties”. Along with the assurance of convergence stability, an added bonus derived from the use of a relaxation procedure is the capability of simulating unsteady-state as well as steady-state operations. Therefore, a numerical algorithm derived from a relaxation approach may be structured with built-in flexibility so that it could be applied to a variety of continuous-contact separation processes. Following these guidelines, an implicit numerical integration procedure, known as DYNCAL, was developed. This analysis-format program was described in a previous article (Hitch et al., 1987). DYNCAL was used to solve the unsteady-state differential-contactor model equations for packed-column absorption and stripping operations involving nonideal, wide-boiling mixtures. Excellent agreement was obtained between DYNCAL predictions and experimental data obtained from a pilot-scale acid gas removal system that utilized refrigerated methanol to absorb carbon dioxide from a coal gasifier off-gas stream. A similar algorithm, named DISTCAL, has now been developed to handle the inherently unsteady-state problem of batch distillation involving narrow-boiling mixtures. So far as we can determine, prior to this work, no rigorous numerical procedure for the simulation of batch distillation in a packed column has ever been presented in the scientific literature. The details of the formulation of the model equations and the mathematical structure of the numerical algorithm, DISTCAL, are discussed below. Simulation predictions are evaluated qualitatively through the use of a parametric study to test the validity of the mathematical model and the numerical solution procedure. In addition, a rudimentary product yield optimization is performed to demonstrate a potential application for the simulation.
Development of Algorithm Before going into the details of the development of the mathematical model and numerical integration scheme 0 1988 American Chemical Society
r BEGIN
Input feed composition and amount and fix boilup rate
J
Determine still temperature, heat duty, and vapor composition
Compute physical and transport property values at time N
1
Solve Equations (1) through (4) to yield new bulk and interfacial mole fractions
1
Ind. Eng. Chem. Res., Vol. 27, No. 8, 1988 1467 the composition. Therefore, it is preferable to use a total energy balance and a total material balance to compute bulk liquid and vapor flow rates, while determining column temperatures via bubble-point (equilibrium) calculations. This type of pairing of dependent variables with model equations is analogous to that of the BP convergence method for staged columns (King, 1980). An analysis-format program, DISTCAL, has been developed using a paired sequential approach to simulate batch distillation in a packed column. As with DYNCAL, the algorithm is built primarily around systems of finite-difference equations that are solved implicitly using matrix decomposition techniques. However, because of the inherent differences between distillation and absorption/ stripping operations, changes have been made in both the formulation of the mathematical model and the computational procedure followed by the program, DISTCAL. The revised model for the batch distillation simulation is given below: unsteady-state individual component material balances
Do bubblepoint calculation with Equation (7) to yield new interfacial temperature
L Solve Equations (5) and (6) simultaneously for new liquid and gas total flow rates
Jl
N = N + 1
mass-transfer relationships a. liquid side abLj/at = ~ ( L I ~ ) + / ~k:;a(xj* z - xj) b. gas side a b , / a t = -a(Gyj)/az - k&a(y; - y j * )
Operation Complete?
equilibrium relationships yj* = m .x .* I 1 unsteady-state, adiabatic energy balance Figure 1. Simplified flow diagram for simulation program DISTCAL.
utilized by DISTCAL, it may be beneficial to review some of the subtle but significant differences between the separation processes of distillation and absorption or stripping. Absorption/stripping operations typically involve a wide-boiling set of components. Indeed, the liquid solvent used should have a low volatility in comparison to the volatility of the gaseous solute or solutes and an even lower volatility with respect to that of the inert or slightly soluble constituents of the gas. For this reason, the total liquid and vapor flow rates are determined most strongly by the equilibrium and mass-transfer relationships. In other words, over a wide range of temperatures, the most volatile components will be present predominantly in the vapor phase and the liquid phase will be made up primarily of the solvent; thus, although changes in the total enthalpy of each phase will have a pronounced effect on the liquid and vapor temperatures, there will not be a dramatic variation in total phase flow rates. Therefore, energy balances and heat-transfer expressions represent the most effective means for predicting column temperatures, while tobd flow rates are best determined by the equilibrium and mass-transfer relationships for absorption/stripping operations. This particular scheme for pairing dependent variables with model equations is analogous to the sum rates (SR) convergence method for staged columns (King, 1980). Distillation, however, typically involves a narrow-boiling set of constituents. In this case, variations in total enthalpy will change the flow rates of each phase but will have little effect on column temperatures, which depend mostly on
dBLHL/dt = d ( L ' H ~ ) / d z d(GrHG)/dZ unsteady-state, total material balance aB,/at
+ aBG/at = aLyaz - aGr/az
summation equation
The only major differences between this unsteady-state model and the one presented earlier (Hitch et al., 1987) for absorption/stripping operations are (1)the use of bulk liquid- and vapor-phase mole fractions and total flow rates rather than component flow rates as model variables and (2) the deletion of the heat-transfer expressions. Removal of the heat-transfer relationships was based on the assumption that, at any vertical position in the column, the bulk liquid- and gas-phase temperatures will be very close to the interfacial temperature. This assumption is common in the case of staged distillation calculations and is likely to be valid for relatively ideal systems such as the simple hydrocarbon mixtures studied in the next section. In any event, the current methods available for the estimation of heat-transfer coefficients for packed columns do not justify the additional computational effort required to evaluate the bulk-phase temperatures in a distillation process. The complete set of 4nc + 3 variables included in this new, modified model are listed in Table I. As with the paired simultaneous approach used by the packed-absorber simulation program DYNCAL, the paired sequential method linearizes the system of partial differential model equations by matching each model variable
1468 Ind. Eng. Chem. Res., Vol. 27, No. 8, 1988 Table I. Model Variables 1. nc liquid component mole fractions 2. nc gas component mole fractions 3. nc liquid interfacial mole fractions 4. nc gas interfacial mole fractions 5. 1 liquid-phase total flow rate 6. 1 gas-phase total flow rate 7. 1 interfacial temperature
XI Y1 X j*
Yj: L
G’
T*
with the particular model equation that most directly determines the value of the variable. In this case, bulk and interfacial mole fractions are paired with the individual component mass balances (eq l),the mass-transfer relationships (eq 2 and 3), and the equilibrium expressions (eq 4); interfacial temperatures are determined through iterative bubble-point calculations by using eq 7 ; and, finally, eq 5 and 6 are solved simultaneously to yield total liquid and vapor flows. The iterative nature of the bubble-point calculations necessitates a sequential computational format. The steps involved in this approach are summarized in Figure 1 which is a simplified flow diagram for the computer program DISTCAL. In batch distillation, there are essentially two modes of operation which Holland and Liapis (1983) refer to as the start-up and product periods. During the start-up period, the column is operated at total reflux until the reflux reaches a desired purity. Because no product is withdrawn, the column may achieve steady state during the start-up period. However, once the product period begins and distillate is withdrawn, unsteady-state operation of the column is inherent because the composition of the liquid in the stillpot is changing continually. The computer program DISTCAL is capable of simulating both modes of operation for batch distillation in a packed column. The numerical procedure is closely related to the operation of the physical process in that it first simulates the start-up period from the point when heat is first added to the stillpot until a steady-state operating condition is achieved. The input variables for the simulation are the quantity and composition of the feed charge (assumed to be at its bubble point), the desired boilup rate from the still, and the height of packing in the column. On the basis of this information, a stillpot temperature and stillpot vapor composition are determined via a bubble-point calculation. This fixes the boundary conditions on the variables associated with the vapor phase at the bottom of the column. The entire column is then “filled” with vapor at the reboiler vapor composition,temperature, and prescribed flow rate (boilup rate). The simulation assumes that the column is equipped with a total condenser, thereby fixing the boundary conditions on the variables associated with the liquid phase, at the top of the column, based on the composition, flow rate, and bubble-point temperature of the (condensed) overhead vapor. Note: Rather than include the highly system-specific heat capacities of the column and internals, we have chosen to neglect these terms in the present analysis. Consequently, the length of the start-up period determined from this simulation may be much shorter than would be required in practice. An improved estimate of start-up behavior could be obtained in a rather straightforward way by including a lumped heat capacity term for the column and internals in the model. Next, physical and transport properties such as densities, viscosities, diffusivities, etc., are estimated (see Appendix) and the coefficients in the block-tridiagonal matrix formed from the finite-difference approximations to eq 1-4 are computed. This matrix is then decomposed to yield new values for the bulk and interfacial mole fractions
throughout the finite-difference grid at an arbitrary time step into the future. This part of the simulation procedure is entirely analogous to the simultaneous-convergence approach of the absorber simulation, DYNCAL, except that only compositions and not temperatures or flow rates are computed by solving the matrix. Instead, column interfacial temperatures are determined in the next step of the sequence via bubble-point calculations using eq 7 and the newly computed interfacial mole fractions. Finally, these new temperatures and compositions are used to compute liquid and vapor specific enthalpies, thus allowing the simultaneous solution of the total energy and material balances, eq 5 and 6, respectively, to yield new total liquid and vapor flow rates throughout the column. At this point, the time counter is incremented by one and the above sequence of calculations is repeated based on the new values for the dependent model variables. Calculations for the start-up period are continued until the values of the calculated variables are no longer changing from one time step to the next, indicating that a steadystate operating condition has been achieved. At this point, the product period is initiated by making a step change from a total reflux operation to one with a specified finite reflux ratio. The sequence of calculations for the product period is exactly the same as that given above for the start-up period. During either mode of operation, the simulation can provide complete concentration and temperature profiles throughout the column, as well as reboiler and condenser heat duties, all as a function of time. During simulation of the product period, distillate compositions and flow rates are also available. Product period computations may be continued for any specified length of time or until a desired stillpot liquid composition is achieved, assuming this composition is physically obtainable in a batch distillation operation. Illustrations from DISTCAL simulations of packed-column batch distillation using simple, saturated hydrocarbon mixtures, are presented in the next section.
Batch Distillation of Low Molecular Weight Alkanes: Parametric Study and Qualitative Evaluation of Algorithm Predictions Faced with the above-mentioned absence of useful data in the literature dealing with batch distillation in packed columns, we decided to attempt a qualitative verification of the DISTCAL algorithm through the application of a parametric study of program input specifications. The input specifications most suitable for variation in a parametric study of an analysis-format simulation are (1) height of packing, (2) boilup rate, and (3) reflux ratio. In addition to these variables, a study of the influence of condenser holdup on algorithm predictions was made. In all DISTCALsimulations, the column was specified to have a nominal diameter of 3.0 f t and to be packed with 0.5in.-diameter Berl saddles. The column operating pressure was assumed constant at 300 psig. The chemical system chosen for simulation in this particular study consisted of three low molecular weight alkanes: propane, n-butane, and n-hexane. The procedure followed for the parametric study was always initialized with the specification of a charge consisting of equal amounts (200 mol) of each of the three alkanes mentioned above. Next, specifications for the height of packing, boilup rate, and condenser holdup were made and the start-up portion of the calculation was begun, “operating” at total reflux until a steady-state operating condition was achieved as described in the section above. At this point, a finite reflux ratio was specified and calculations for the product period were begun. These
Ind. Eng. Chem. Res., Vol. 27, No. 8, 1988 1469 DISTCAL Simulation
DISTCAL S i m u l s t m n 1
0.2
0.
0
1
2
3
4
5
7
6
8
9 1 0 1 1 1 2 1 3 1 4 1 5
Time (houra)
-Height
Time (hours)
-Bsae
= 10 ft. ..... Height = 20 ft ..... Height = 40 it.
input specifications
Figure 2. Condenser concentrations during product period using base input specifications.
Figure 3. Effect of varying packing height on product concentrations. DISTCAL Simulation
Table 11. h D u t SDecifications for DISTCAL Simulations boilup height of rate, condenser R / V ratio holdup, mol simulation packing, ft mol/h 150.0 0.75 0.0 run la 10.0 0.75 run 2 20.0 0.0 150.0 0.0 150.0 0.75 40.0 run 3 100.0 0.0 0.75 run 4 10.0 200.0 0.0 0.75 run 5 10.0 0.50 10.0 0.0 150.0 run 6 0.90 0.0 150.0 run 7 10.0 0.75 10.0 25.0 150.0 run 8 0.75 50.0 150.0 run 9 10.0
1.
"";,...p-./
t
Concentration
Base input specifications. Feed composition = 200 mol of propane, 200 mol of n-butane, and 200 mol of n-hexane. Column pressure = 300 psig.
calculations continued until the composition of the liquid remaining in the still was essentially pure n-hexane. The base input specifications for this parametric study were a column height of 10 ft, a boilup rate of 150 mol/h, a ratio of liquid reflux to overhead vapor flow (R/V ratio) of 0.75, and a condenser holdup of 0.0 mol. Figure 2 is a plot of the composition of the material in the total condenser as a function of time as computed by DISTCAL for the product period and the base case values given above. At the start of the product period, the distillate is enriched in the most volatile component, propane. As the amount of propane present in the still is depleted, the intermediate component, n-butane, becomes the dominant constituent present in the condenser. Finally, once the removal of the more volatile components is complete, only n-hexane, the least volatile component, appears in the distillate. Because these plots are indicative of the qualitative behavior of the system as a function of operating conditions, the focus of this parametric study will be on how changes in input specifications affect the appearance of distillate compo-
-
Time (houra) _ .
Boilup Rats
_ _Boilup _ . . Rat.
1W lb-mole/hr
= 150 lb-mole/hr
..... Boilup Rate = 200 lb-mole/hr
Figure 4. Effect of varying boilup rate on product concentrations.
sition versus time curves such as these. In this parametric study, a total of nine simulations were performed by varying each of the four input specifications in turn while keeping the others constant at their particular base values. A complete list of the input specifications for
1470 Ind. Eng. Chem. Res., Vol. 27, No. 8, 1988 DISTCAL Simulation 1q
'
w/f 0
2
i
,
5
I,,
4
5
~,
t
'
1 8
10
12 14 16 18 20 22 2 4 2 8 2 8 30 32
34 36 38 40 Time (hours)
Time (hours) ~
R / V Ratio = 0 5 0 R / V Ratio = 0 75 R / V R a t i o = 0 90
Figure 5. Effect of varying reflux ratio on product concentrations.
these nine simulations is given in Table 11. Figures 3-6 show the effects of the different input specification variations on the predicted condenser mole fractions of each of the three chemical components in the system. These plots show the transient response of the column during the product period only. The initial value of the mole fraction (i.e., at time = 0 h) shown for each component is actually the steady-state concentration of that particular constituent achieved during the start-up period, operating at total reflux. Figure 3 demonstrates the effect on the mole fractions of propane, n-butane, and n-hexane, respectively, of a change in the specified height of packing. The three packing heights studied were 10, 20, and 40 ft. A quick examination of these plots reveals that, as expected, a greater degree of separation is obtainable as the packing height increases. This figure shows that increasing the height of packing results in an increase in the purity of propane in the distillate product during the initial part of the product period. Likewise, the purity of the intermediate-volatility component, n-butane, is enhanced during the middle of the product period by an increased packed height. Finally, it can be seen that the mole fraction of n-hexane, the least volatile component, asymptotically approaches 1.0 as time increases for all three column packing heights studied; this asymptote is approached most rapidly for the tallest packing height of 40 ft, indicating that the more volatile components are removed faster with a larger packing section. All of these observations conform to the known fact that, for a masstransfer-controlled system, the net rate of mass transfer and thus the ultimate degree of separation obtainable at any given point in time will increase with an increase in interfacial contact area. In Figure 4, the impact of varying the boilup rate on the time-dependent composition of the distillate is examined. Three different boilup rates were used: 100, 150, and 200
__ Holdup = ~ ~ ~ ....
0 . 0 lb-moles . . Holdup = 25.0 lb-moles Holdup = 50.0 lb-moles
Figure 6. Effect of varying condenser holdup on product concentrations.
mol/h. These plots reveal that increasing the boilup rate has only a minor effect on the ultimate degree of separation obtained, but it does significantly accelerate the removal of the higher volatility components from the still. Again, this result is consistent with the intuitive expectation that changes in vapor flow rate (at constant reflux ratio) will not drastically affect the separation achieved in a column in which only rectification occurs. What little impact the changes in boilup rate do have can be attributed to a slight enhancement of the gas- and liquid-phase mass-transfer coefficients as a result of the higher throughput rates. Figure 5 demonstrates how changes in the reflux to overhead vapor (Rl V ) ratio affect the distillate composition. As explained above, the simulation of the product period is initiated by making a step change from a total reflux condition (RIV = 1.0) to operating with partial reflux. This figure shows that the highest purity of propane in the condenser is obtained under conditions of total reflux. Once the reflux-to-vapor ratio is lowered, the propane concentration drops by an amount that is roughly proportional to the decrease in R / V . Thus, as expected, the higher the RIV ratio, the higher the maximum obtainable propane purity at the start of the product period. Similarly, one can see that the maximum purity of the intermediate component, n-butane, is achieved with the highest value of R / V. Naturally, for a given boilup rate, a higher R / V ratio will increase the time necessary to remove completely all of the higher volatility components; therefore, the "peak" indicating a maximum in n-butane concentration is made wider and occurs later as R l V increases. As before, all of these plots are consistent with the behavior expected based on our practical knowledge of the chemical system. Figure 6 illustrates what can happen to column performance as a result of increasing condenser holdup. All of the curves shown to this point were based on the as-
Ind. Eng. Chem. Res., Vol. 27, No. 8, 1988 1471 sumption of negligible holdup in the total condenser. This implies that the reflux being returned to the column must have the exact same composition as the overhead vapor being withdrawn at any given point in time. This assumption is probably not valid for a typical industrial facility. Instead, there will be some finite residence time associated with the condenser and the condensate receiving tank. If the condensate receiving vessel is considered to be well-mixed,then the composition of the liquid in every location within the tank is uniformly equal to the composition of the liquid exiting the tank. However, at any given point in time, the composition in the tank will, in general, be different from the composition of the overhead vapor flowing out of the column. The extent to which the composition of the liquid in the receiving tank differs from the composition of the overhead vapor depends on the size of the tank. From Figure 6, it can be seen that the greater the condenser holdup, the greater the concentration of propane in the distillate at any given point in time. The reason for this is that accumulation of liquid in the condenser receiving vessel effectively “averages”the composition of the vapor currently leaving the column with the compositions of all the overhead vapor which preceded it. Since the concentration of propane in the overhead vapor decreases monotonically with time, the condensate that has been held up in the receiving tank must have a higher concentration of propane than that of the overhead vapor currently being produced. By the same token, as the condenser holdup increases, the peak representing the maximum obtainable purity of n-butane occurs later and is reduced in height as a result of the “averaging” effect mentioned above. As before, these results are all in line with intuitive expectations. Aside from the predictions of condenser mole fractions, a statement should be made about the validity of the calculations of total vapor and liquid flow rates. Figure 7 reveals how the total overhead vapor flow rate varies with time during the product period of a simulation using the base input specifications. Notice the sudden increase in the overhead vapor flow rate (depicted by the uppermost curve) at the start of the product period. This is due to the decrease in the amount of cold liquid being returned as reflux to the column. After this initial rapid increase in the vapor flow, there is a long period of steady but slow decline followed by a major drop. The vapor flow rate finally “bottoms out” at about 15 mol/(h.ft2) and then begins a rapid recovery almost to the original rate of flow. The three lower curves have been added to Figure 7 to help explain the large dip in the vapor flow curve. It readily can be seen that the onset of the rapid drop in the overhead vapor flow roughly coincides with the beginning of the decline in the reboiler concentration of the intermediate component, n-butane. As the concentration of n-butane approaches zero, the overhead vapor flow rate reaches a minimum and then begins to rise. In order to explain this behavior, recall that the simulation is being run with a constant molar boilup rate from the reboiler. However, the overhead vapor flow depends on the heat of vaporization/condensation of the material at the top of the column which in turn is a function of composition. As long as a significant amount of higher volatility (lower heat of vaporization) material remains in the still, the vapor flowing up the column is composed mainly of the lower boiling components and the heat supplied to the reboiler can support a relatively constant vapor flow through ut the column. However, once this low-boiling materia is depleted, the average heat of vaporization of the material
P
DISTCIIL Simulation 251
,10
-.;
,”
D
-\
l
d
- _, 2
4
8
6
10
12
14
0
Time (hours)
__
Overhead Flow
Propane Conc ..... h-Butane Conc Y-Hexane Conc
Figure 7. Variation of overhead vapor flow rate during product period using base input specifications.
flowing up the column increases dramatically and the overhead vapor flow rate drops until the rate of heat addition to the reboiler increases enough to compensate. Thus, the vapor flow predictions of DISTCAL are once again consistent with physical intuition.
Simulation Applications: Optimization of Product Yield Having demonstrated that the batch distillation simulation program DISTCAL is qualitatively sound, there is every reason to expect that, given a set of valid physical parameters, the algorithm is capable of simulating an actual batch distillation process with quantitative accuracy. Naturally, if this were the case, then the algorithm would be a useful predictive tool for analyzing and optimizing various process operating conditions before they are attempted in the plant, and the algorithm could be used in devising and driving a control strategy for the operation. On the basis of the assumption that the current set of model parameter correlations available to DISTCAL are valid for the conditions of the simulation, a rudimentary attempt at optimizing the yield of propane in the distillate product has been made. The results of this optimization effort are depicted in Figures 8-10. The goal of the optimization was to maximize the yield of propane while maintaining a concentration of propane in the continuously accumulated distillate product of no less than 78 mol %. This optimization was accomplished by increasing the R l V ratio and boilup rate as the cumulative product composition began to dip below the prescribed value of 78 mol % propane. In order to establish a basis for comparison, the simulation was first run at the standard conditions listed above with a constant R f V ratio of 0.75 and a constant boilup rate of 150 mol/h until the minimum acceptable distillate product concentration of 78% propane was reached. As can be seen from Figure 8, this composition
1472 Ind. Eng. Chem. Res., Vol. 27, No. 8, 1988 DISTCAL Simulation
DISTCAL Simulatioi
low
I O
tl 0
10-
,/ “
I
P n
50-
07
3 0
0 5
1 5
1 0
Time (houra)
2 0
2.5
3 0
3 5
e o
40
Time (hours)
__ Total .....
Amount of Proprns Yole F T a c t x o n Proproe
~
Figure 8. Yield of propane in a 78 mol % product obtained using a constant R / V ratio of 0.75 and a constant boilup rate of 150 mol/h.
Total Amount of Propane Mole Fraction Propane
Figure 10. Yield of propane in a 78 mol % product obtained using a variable R/ V ratio (0.75-0.995) and a variable boilup rate (15C-750 mol/h).
DISTCAL Simulation l0cq
1 0
!
+OQ
d
I’
d
0.0
~
~
0.5 1.0
~
1.5
:
2.0
:
~
0.5
~
3.0 3.5
:
4.0
~! 4.5
~,
5.0
5.5
bo
6 . 0 6.5
Time (hours)
-T o t s 1 ..~.. Mole
Amount of Propane
Fraction Propane
Figure 9. Yield of propane in a 78 mol % product obtained using a variable R / V ratio (0.75-0.995) and a constant boilup rate of 150 mol/h.
was obtained after approximately 1 h of operating time with a total yield of propane in the distillate product of about 30 mol.
Figure 9 depicts a DISTCAL simulation in which the product period was initiated with the standard operating conditions, i.e., R / V = 0.75 and boilup rate = 150 mol/h. The R / V ratio was then increased incrementally whenever the cumulative propane product concentration began to drop below the minimum until near total reflux was obtained and no further product was being withdrawn. During this time, the boilup rate maintained constant at 150 mol/h. Figure 8 shows that this process took just under 7 h and that a total propane yield of about 70 mol was achieved. Thus, although the product yield was over twice that of the previous run,the operating time increased by a factor of 7 . In the final simulation, illustrated by Figure 10, both the R / V ratio and the boilup rate were increased incrementally as the propane product concentration dropped. The increases were made so as to maintain a constant product flow rate as demonstrated by the constant slope of the curve, indicating propane yield as a function of time. This procedure was followed until an R / V ratio of 0.995 was achieved, a t which point the boilup rate had reached a relatively high value of 750 mol/h. The time required for this operation was about 4 h and the resulting propane yield was approximately 90 mol, thus representing a 3-fold increase in propane yield at the expense of a 4-fold increase in operating time as compared with the original “nonoptimized” base case operation. As mentioned above, this optimization procedure is crude at best and is meant only to illustrate the capabilities of the algorithm. A more sophisticated optimization approach would no doubt produce a better operating scheme. Conclusion From the onset of this work, it was known that many significant gaps existed in the literature dealing with
Ind. Eng. Chem. Res., Vol. 27, No. 8, 1988 1473 packed-column simulation. One almost universal area of deficiency is the problem of poor convergence characteristics when dealing with anything but the simplest of physical systems. The source of such inadequate convergence properties can almost always be traced to the difficulties involved in numerically handling the large systems of highly nonlinear differential equations that describe the phenomena of simultaneous heat and mass transfer in continuous contactors. In addition to unsatisfactory convergence attributes, some other areas of a problematical nature include the unavailability of simulation procedures capable of addressing unsteady-state operating conditions and, consequently, the complete lack of algorithms dealing with batch processes such as batch distillation. The objective of this research effort, then, was to develop a numerical approach that would address the problem of unsteady-state calculations in packed columns, paying particular attention to the need for an assured stability of convergence. Toward this end, the development of the DISTCALsimulation program represents a considerable improvement in the software available for transient simulation of separation processes in packed columns. The algorithm accommodates a rigorous set of model equations that may include both nonideal VLE relationships and nonlinear mass-transfer expressions. In addition, the relaxation procedure used is very stable, and a very good approximation to the transient behavior of the system may be obtained even when employing relatively large time steps in the implicit integration procedure. The accuracy of the DISTCAL simulation predictions has been tested qualitatively, and in every instance, the algorithm was found to conform to the expectations of its behavior derived from practical knowledge of the processes being simulated as well as from intuitive reasoning. The ultimate accuracy with which any algorithm can simulate a particular situation, of course, depends on the validity of the correlations used to compute the physical parameters demanded by the model. These correlations as well as the sources for all physical property data used in this work are given in the Appendix. At this point, no system-specific adjustable parameters are used in the simulation; however, adjustable parameters may be readily included should a higher level of accuracy be required than can be realized through the use of general empirical correlations.
Nomenclature a = gas-liquid interfacial area per unit volume of packing, area/volume b , = molar gas-phase holdup of component j , moles/volume b , = molar liquid-phase holdup of component j, moles/ volume BG = total molar gas-phase holdup, moles/volume BL = total molar liquid-phase holdup, moles/volume G' = total superficial molar flow rate of gas phase, moles/ (timesarea)
k 'y, = gas-phase mass-transfer coefficient for component j , moles/ (time-area) k '.,= liquid-phase mass-transfer coefficient for component J , moles/(time.area) L' = total superficial molar flow rate of liquid phase, moles/ (timearea) mi= equilibrium distribution coefficient for component j , mole fraction/mole fraction nc = number of chemical components t = time x j = bulk liquid mole fracton of component j x,* = interfacial liquid mole fraction of component j yj = bulk gas mole fraction of component j y,* = interfacial gas mole fraction of component j z = axial coordinate, length
Appendix: Source Material for Physical Property Data Interfacial area: Treybal (1980, pp 204-205). Liquid-phase specific enthalpy: Holland and Liapis (1983, pp 172-173). Gas-phase specific enthalpy: Holland and Liapis (1983, pp 172-173). Liquid-phase mass-transfer coefficients: Shulman et al. (1955). Gas-phase mass-transfer coefficients: Treybal (1980, p 203). Vapor/liquid equilibrium distribution coefficients: Holland and Liapis (1983, pp 172-173). Liquid density correlation: Yen and Woods from Reid et al. (1977). Liquid surface tension: Reid et al. (1977). Liquid-phase viscosity: Reid et al. (1977, Appendix A). Gas-phase viscosity correlation: Dean and Stiel from Reid et al. (1977, p 419). Literature Cited Hitch, D. M.; Rousseau, R. W.; Ferrell, J. K. Ind. Eng. Chem. Res. 1987,26,1092. Holland, C. D.; Liapis, A. I. Computer Methods for Solving Dynamic Separation Problems, 1st ed.; McGraw-Hill: New York, 1980. Reid, R. C.; Prausnitz, J. M.; Sherwood, T. K. The Properties of Gases and Liquids, 3rd ed.; McGraw-Hill: New York, 1977. Shulman, H. L.; Ullrich, C. F.; Proulx, A. A.; Zimmerman, J. 0. "Wetted and Effective Interfacial Areas, Gas and Liquid Phase Mass Transfer Rates". AIChE J. 1955,I, 253. Treybal, R. E. Mass-Transfer Operations, 3rd ed.; McGraw-Hill: New York, 1980. Von Stockar, U.;Wilke, C. R. "Rigorous and Short-Cut Design Calculations for Gas Absorption Involving Large Heat Effects. 1. A New Computational Method for Packed Gas Absorbers". Ind. Eng. Chem. Fundam. 19774 16,88. Von Stockar, U.;Wilke, C. R. "Rigorous and Short-Cut Design Calculations for Gas Absorption Involving Large Heat Effects. 2. Rapid Short-Cut Design Procedure for Packed Gas Absorbers". Ind. Eng. Chem. Fundam. 1977b,16,94. Received f o r review June 29, 1987 Revised manuscript received February 29, 1988 Accepted April 15, 1988