## Calculational Methods

Rigorous Computer-Based Calculation Procedures It is obvious that a set of curves such as shown in Fig. 13-104 for a binary mixture is quite tedious to obtain by hand methods. The curves shown in Fig. 13-99 for a multicomponent batch distillation are extremely difficult to develop by hand methods. Therefore, since the early 1960s, when large digital computers became available, interest has been generated in developing rigorous calculation procedures for binary and multicomponent batch distillation. For binary mixtures of constant relative volatility, Huckaba and Danly [Am. Inst. Chem. Eng. J., 6, 335 (1960)] developed a computer program that assumed constant-mass tray holdups, adiabatic tray operation, and linear enthalpy relationships but did include energy balances around each tray and permitted use of nonequilibrium trays by means of specified tray efficiencies. Experimental data were provided to validate the simulation. Meadows [Chem. Eng. Prog. Symp. Ser. 46, 59, 48 (1963)] presented a multi-component-batch-distillation model that included equations for energy, material, and volume balances around theoretical trays. The only assumptions made were perfect mixing on each tray, negligible vapor holdup, adiabatic operation, and constant-volume tray holdup. Distefano [Am. Inst. Chem. Eng. J., 14, 190 (1968)] extended the model and developed a computer-based-solution procedure that was used to simulate successfully several commercial batch-distillation columns. Boston et al. [Foundations of Computer-Aided Chemical Process Design, vol. II, ed. by Mah and Seider, American Institute of Chemical Engineers, New York, 1981, p. 203) further extended the model, provided a variety of practical sets of specifications, and utilized modern numerical procedures and equation formulations to handle efficiently the nonlinear and often stiff nature of the multicomponent-batch-distillation problem. The simpler model of Distefano is used here to illustrate this nonlinear and stiff nature.

Consider the simple batch- or multicomponent-distillation operation in Fig. 13-105. The still consists of a pot or reboiler, a column with N theoretical trays or equivalent packing, and a condenser with an accompanying reflux drum. The mixture to be distilled is charged to the reboiler, to which heat is then supplied. Vapor leaving the top tray is totally condensed and drained into the reflux drum. Initially, no distillate is withdrawn from the system, but instead a total-reflux condition is established at a fixed overhead vapor rate. Then, starting at time t = 0, distillate is removed at a constant molal rate and sent to a receiver that is not shown in Fig. 13-105. Simultaneously, a fixed reflux ratio is established such that the overhead vapor rate is not changed from that at total reflux. Alternatively, heat input to the reboiler can be maintained constant and distillate rate allowed to vary accordingly. The equations of Distefano for a batch distillation operated in this manner are as follows (after minor rearrangement), where i,j refers to the ith of C components in the mixture and the fth of N theoretical plates.

Component mole balances for total-condenser-reflux drum, trays, and reboiler, respectively:

dM0 -i dxio lit

Mj dMj dt

—,j+1—+1 m, xij+1 i = 1 to C j = 1 to N (13-150)

Mn+1

Total mole balance for total-condenser-reflux drum and trays respectively:

Energy balance around ,th tray:

where HV and HL are molar vapor and liquid enthalpies respectively. Phase equilibriums:

Mole-fraction sum:

Molar holdups in condenser-reflux drum, on trays, and in reboiler: Mo = G0P0 (13-157)

where G is the constant-volume holdup, M%+1 is the initial molar charge to reboiler, and p is the liquid molar density.

Energy balances around condenser and reboiler respectively:

Equation (13-160) is replaced by the following overall energy-balance equation if QN+1 is to be specified rather than D:

dM0 lit

ld(MjH,J

With D and R specified, Eqs. (13-149) to (13-161) represent a coupled set of (2CN + 3C + 4N + 7) equations constituting an initial-value problem in an equal number of time-dependent unknown variables, namely, (CN + 2C)xi,j, (CN + C)yi,j, (N)Lp (N + 1)V;, (N + 2)Tj, (N + 2)Mj, Q0, and QN+x, where initial conditions at t = 0 for all unknown variables are obtained by determining the total-reflux steady-state condition for specifications on the number of theoretical stages, amount and composition of initial charge, volume holdups, and molar vapor rate leaving the top stage and entering the condenser.

Various procedures for solving Eqs. (13-149) to (13-161), ranging from a complete tearing method to solve the equations one at a time, as shown by Distefano, to a complete simultaneous method, have been studied. Regardless of the method used, the following considerations generally apply:

1. Derivatives or rates of change of tray and condenser-reflux drum liquid holdup with respect to time are sufficiently small compared with total flow rates that these derivatives can be approximated by incremental changes over the previous time step. Derivatives of liquid enthalpy with respect to time everywhere can be approximated in the same way. The derivative of the liquid holdup in the reboiler can likewise be approximated in the same way except when reflux ratios are low.

2. Ordinary differential Eqs. (13-149) to (13-151) for rates of change of liquid-phase mole fractions are nonlinear because the coefficients of xi j change with time. Therefore, numerical methods of integration with respect to time must be employed. Furthermore, the equations may be difficult to integrate rapidly and accurately because they may constitute a so-called stiff system as considered by Gear (Numerical Initial Value Problems in Ordinary Differential Equations, Prentice Hall, Englewood Cliffs, N.J., 1971). The choice of time step for simple explicit numerical procedures (such as the Euler and Runge-Kutta methods) of integrating sets of ordinary differential equations in initial-value problems may be governed by either stability or truncation-error considerations. Truncation errors in the dependent variables may be scarcely noticeable and generally accumulate gradually with time. Instability generally causes sudden and severe errors that are very noticeable. When the equations are stiff, stability controls and extremely small time steps may be necessary to prevent instability. A common measure of the severity of stiffness is the stiffness ratio iXimax/!Ximin, where X is an eigenvalue for the jacobian matrix of the set of ordinary differential equations. For Eqs. (13-149) to (13151), the jacobian matrix is tridiagonal if the equations and variables are arranged by stage (top down) for each component in order.

For a general jacobian matrix pertaining to C components and N theoretical trays, as shown by Distefano [Am. Inst. Chem. Eng. J., 14, 946 (1968)]. Gerschgorin's circle theorem (Varga, Matrix Iterative Analysis, Prentice Hall, Englewood Cliffs, N.J., 1962) may be employed to obtain bounds on the maximum and minimum absolute eigenvalues. Accordingly, iXimax ^ max

The maximum absolute eigenvalue corresponds to the component with the largest K value (KLj) and the tray with the smallest holdup. Therefore, if the derivative term and any variation in Lj, Vj, and Ki j are neglected, i^ima, = 2

Lj + kijvj

In a similar development, the minimum upper limit on the eigenvalue corresponds to the component with the largest K value and to the largest holdup, which occurs in the reboiler. Thus i^imin =

Therefore, the lower bound on the stiffness ratio at the beginning of batch distillation is given approximately by iXimax = 2 (MN + t iXimin ( Mn where MN+1 and MN are the molar holdups in the reboiler initially and on the bottom tray respectively. In the sample problem presented by Distefano (ibid.) for the smallest charge, the approximate initial-stiffness ratio is of the order of 250, which is not considered to be a particularly large value. Using an explicit integration method, almost 600 time increments, which were controlled by stability criteria, were required to distill 98 percent of the charge.

At the other extreme of Distefano's sample problems, for the largest initial charge, the maximum-stiffness ratio is of the order of 1500, which is considered to be a relatively large value. In this case, more than 10,000 time steps are required to distill 90 percent of the initial change, and the problem is better handled by a stiff integrator.

In Distefano's method, Eqs. (13-149) to (13-161) are solved with an initial condition of total reflux at L0 equal to D(R + 1) from the specifications. At t = 0, L0 is reduced so as to begin distillate withdrawal. The computational procedure is then as follows:

1. Replace Lj1 by Lj1 - D, but retain Vj1 and all other initial values from the total-reflux calculation.

2. Replace the holdup derivatives in Eqs. (13-149) to (13-151) by total-stage material-balance equations (e.g., dMj /dt = Vj+1 + Lj _ 1 -Vj - Lj) and solve the resulting equations one at a time by the predictor step of an explicit integration method for a time increment that is determined by stability and truncation considerations. If the mole fractions for a particular stage do not sum to 1, normalize them.

3. Compute a new set of stage temperatures from Eq. (13-156). Calculate a corresponding set of vapor-phase mole fractions from Eq. (13-155).

## Post a comment