Nunge, R. ., and R. Ad am , “Revers “Revers Osmosis in Laminar Flow through through Curved Tubes,” Desalination, 13,17 (1973). Patankar, S. V., V. S. Pratap, and D. B. Spalding, “Prediction “Prediction of L aminar Flow and Hea t Transfer in Helically Coiled Pipes,” Pipes,” 1.Fluid Mech., 2, 539 (1974). Peaceman, D. ., and H. H. Rachford, “Th e Numerical Solution f ParJ.Soc. ndust. Awl. Math., abolic and Elliptic Differential Equations,” J.Soc. 3,28 (1955). Raju, K. K., an S. L. Rathn a, “H eat Tran sfer for the Flow of of a Power-Law Fluid in a Curved Tube, Indian Imt. Sci., 52,34 (1970). Ranade, V. R., an J. J. Ulbrecht, “Th e Residence Tim e Distribution for Lam inar Flow of Non-New tonian Liquids through through H elical Tubes,” Chem. Eng. Cummun., 8,165 (1981). Ruthven, D. M., “T he Residence Residence Time Distribut Distribution ion for Ideal Laminar Flow in Helical Tubes,” Chem. Eng. S c i . , 26,1,113 (1971).
Seader, J. D., and L. M. outhwick, “Saponification Ethyl A cetate in Curved T ub e Reactors, Reactors, Chem. Eng. Cornm., 9,175 (1981). Srinivasan, P. S., S. S. Nand apurk ar, and F. A. Holland, “Pressure Drop an Hea t Trans fer in Coils, Coils, The Chem. Engr., 113, London (1968). Trivedi, N. nd K. Vasudeva, “Axial “Axial Dispersion in Laminar Flow in Helical Coils,” Chem. Eng. Sci., 30,317 (1975). Weissman, M. ., and L. F. Mockros, “Gas Transfer to Blood Flowing in Coiled Circular Tubes,” Proc. Am. Soc. Engrs. Eng. Mwh. Diu.), 94,857 (1968). (1968). Yao, L. S. an S. A. Berger, Berger, ‘‘Enhy Flow in a Curv ed Pipe,” Fluid Meck., 67,177 (1975). ManusMipt received Aug. 24,1982; reoLpion received Oct. 20,1983, and accepted accepted Oct. 20
Nonequilibrium Nonequilibrium Stage Stage Model Multicomponent Separation Processes Part
M etho of Solution Model Description and Metho
nonequilibrium stage model is developed for the simulation countercurrent multicomponent separation processes feature of the model is that th e component materi al and energy balance relations for each phase together with mass and energy transfer rate equations and equilibrium equati ons for the phase interface a re solved solved to fi nd the actua l separation directly. Computations of of stage efficiencies are entirely avoided procedure for solving the model equations simultaneously simultaneously using Newton’s method is outlined.
KRISHNAMURTHY and TAYLOR Department of Chemical Engineering
Clarkson University Potsdam, NY 13676
SCOPE Rigorous simulation of multistage processes processes such as distillation or absorption is, more often than not, based upon the equilibrium stage model. This model is well enough known not to need a det aile d description here (see, for example, example, the textbooks by King, 1980, p. 446; Henley and Seader, 1981, pp. 24, 557; Holland, 1975, 47; 1981, p. 6). Briefly, of of course, the model includes the assumption that the streams leaving leaving any particular stage are in equilibrium with each other. other. Component Component Ma teri al balances, balances, the equations of of phase Equil ibrium, Summation equations, and He at balance for each stage (the so-c so-call alled ed MESH equations)are equations)are solved using one of of the very many ingenious algorithms presently available to give product distributions, flow rates, temperatures, and so on. In ac tual operation, stages rarely, rarely, ever, ever, operat e at equilib rium despite attempts to approach this condition by proper design and choice of operating conditions. The usual way of dealing with departures from from equilibrium equilibrium is by incorporating a stage efficiency into the equilibrium relations. I t is with th introduction of this quantity that the problems begin. The fi rst problem is that there ar e several different definitions stage efficiency: Murphree (1925), Hausen (1953), generalized R. Krishnarnurthy presently with The BOC Group, Murray Hill, NJ 07974
AlChE Journal (Vol.
No
Hausen (Standart, (Standart, 1965), vaporization (Holland, 1975, p. 268), and others. There is by no means a consensus consensus on which definition is best. Arguments for and against various possibilities are presented by, amon g others, Standar (1965, 1971), Holland (1975, pp. 268,327), Holland a nd McMahon (1970), King (1980, p. 637), and Medina et al. (1978,1979). Possibly the most soundly based (in a thermodynamic sense), sense), the generalized Hausen efficiencies are ridiculously complicated to calculate; the least soundly based, based, t he Murphree efficiency, is the one most widely used because it is easily combined with the equilibrium equations. Therma l efficiencies efficiencies may also be defined but almost always ar e taken to have a value f unit y (except, (except, perhaps, when inert species are involved). involved). Whichever definition f stage efficiency is adopted, it must either be specified in advance or calculated from an equa tion derived by dividing by some reference separation t he actual separation obtained from a solution the component material balance equations for each phase. tions of the comp onent phase balances) have been proposed for binary systems (King, 1980, p. 612). The extension of these models to mul ticomponent systems provides further complications. In component system there are independent com-
March, 1985
Page 449
ponent efficiencies, which, for lack f anything bet tter to do, have usually been taken to be the same for all components (King, 1980, 636). In fact, th e individual component efficiencies are rarely eq ual, simply because different species exhibit varying facilities for mass transfer. Mass transfer in multicomponent mixtures is more complicated than in binary systems because of t he possible coupling between the individual concentration gradients. Phenomena such as reverse diffusion (diffusion of a species again st its own concentration gradient) or osmotic diffusion (diffusion of species even though no concentration gradi ent for that species exists) are possible in multicomponent systems but not in binaries (Toor, 1957). One of the intere sting consequences f these interaction effects is that the individual point efficiencies of different species are not constrained to lie between zero and one. Instead, they may be found anywhere in the range from --oo to Toor, 1964). Binary as well as multicomponent tray efficiencies may, for other reasons (e.g., weeping or entrainment), be greater than one. Experimental confirmation of Toor’s (1964) prediction has recently been published by Krishna et al. (1977). Vogelpohl (1979) provides further evidence that the component point efficiencies are unequal and unbounded. The subst anti al body of other dat a affirm ing that individual effici enci es are unequal is reviewed by Krishna et al. (1977). Models f mass transfer t hat are ab le to account for interaction effects now are available (see Krishna and St andar t, 1979, for a recent review). In general, they lead t o rate equati ons of th e form of Eq 13 below. These models have been used as a basis for developing methods for calculating efficiencies in multicomponent systems (Toor and Burchard, 1960; Toor, 1964; Diener and Gerster, 1968; Krishna et al., 1977; Medina et al.,
1979; Vogerphl, 1979). However, these models have been tested against a rather limited quantity of data , and, with the exception f the ea rly work of Toor a nd Burchard (1960), there has not yet been an attempt to combine these efficiency models with a method for simul ating multistage separat ion processes. In fact, some care in doing this would have to be exercised in order to avoid violating the material balances (King, 1980, p. 637). All things considered we believe that the someti mes arbitrary and ambiguous multicomponent stage efficiency adds unnec essary complexity to separation process modeling. It is our contention that multistage separation processes are more effective ly modeled by a sequenc e of nonequili brium stages for which t he soluti on of the conservat ion equations for each phase is used directly (rathe r tha n indirectly i n the form of an efficiency factor). It is the objective of this paper to present a fairly general nonequilibrium stage model of countercurrent multicomponent separation processes and to describe a method of solving the resulting set of equations. Other nonequilibrium models of staged equi pmen t in the sam e class as ours are due to Waggoner and Loud (1977) (distillation), to Waggoner and Burkhart (1978)(extraction),and to Ricker et al. (1981) extraction with backmixing). The simplifications made by them (but no by us) will be mentioned at appropriate points in the text. Probably the most sophisticated twophase model is that of Billingsley and Chirachavalla (1981) for continuous contact equi pmen t. Some f t he properties of component efficiencies in staged equipment (unequal and unbounded) are shared by component NTUs in continuous contact equipment, properties that Billingsley and Chirachavalla included in their model. The model proposed here may be (and is in Part II) applied to both stagewise and continuous contact equipment.
CONCLUSIONS AND SIGNIFICANCE nonequilibrium stage model of multicomponen t separation processes has been developed. The features of the model are: 1. The mass and energy conservation equations ar e split into two parts, one for each phase. The equations for each phase ar connected by mass and energy balances around the interface and by the assumption that the interface be a t thermodynamic equilibrium. 2. The process of simultaneous mass and energy transfer through the interface is modeled by means of ra te equations and transfer coefficients. Abstract and somewhat arbitrary concepts like efficiency or number of transfer units are avoided entirely.
NONEOUlLlBRlUM STAGE MODEL schematic representation of a nonequilibrium stage or section f a packed or wetted-wall column is shown in Figure 1. Packed towers and multistage columns consist of a sequence such stages. Vapor and liquid streams from adjacent stages are brought into contact on the stage and allowed to exchange mass and energy across their common interface represented in the diagram by the vertical wavy line. The stage is assumed to be at mechanical equilibrium; Provision is made for vapor and liquid feed streams, side stream diawoffs vapor and liquid, and for the addition or removal of heat. Steady state operation is assumed. The equations used to model the behavior this stage are now presented. Conservation Relations with all models of chemical processes, the analysis f the nonequilibrium st age starts with the construction f material and
Page 450
March, 1985
3. Resistances to mass and energy transfer offered by both flui d phases can be accounted for by using separ ate rate equations for each phase. Interfacial effects are ignored in the present development. 4. he model is formu late d in such a way as to make the describ ing equat ions and th e method of solving them largely indepen dent of the methods used to predict th e transfer rates as well as physical and thermodynamic properties of the system. 5. The e quat ions are best solved simultaneously using Newton’s met hod or one of its relatives.
energy balances. Th e mass balance for component
on stage is
where
In equilibrium stage calculations, Eqs 1- are solved subject to the requirement that the vapor and liquid streams leaving stage are in com plete thermal, mechanical, and chemical equilibrium. In our development of the nonequi librium stage model, Eqs 1are not useful as they stand. Instead, the stagebalances are split into two parts, one for each phase. For the vapor phase, the component mass balance is
AlChE Journal (Vol. 31, No.
_yl_m/-+ Interface
Ll-I
Liquid film
opor film
1.1
Vapor sidesiream
Bulk vapor
stage
Bulk liguid
iquid
Ni.j
Vapor feed
Liquid feed fk
Vj*l
HF
interface
41
Figure 2. Typical composition and temperature profiles in the region of the
Yi.p
TL
sL Figure
Schematic represen!atlon of
ry,.tj
nonequllibrlumstage.
(4
i,j+l-f;
and for th e liquid phase rf)lij
v$
,j-1-
The last terms in Eqs. an represent the net loss or gain of species du e to interphase transport. Formally, we may write
JV;
SN;da,;
SNkdaj
(6
where Nt is the molar flux of species at a particu lar point in the two-phase dispersion, and da represents the small amou nt of interfacial area through which that flux passes. We ado pt the convention that transfers from the vapor phase to the liquid phase are positive. It follows directly from Eqs. 1 , 4 , an that
Nh=O
M:
(7)
a result that may als be derived straightforwardlyby con structing a material balanc e around the en tire interface (hence the notation in Eq. 7). Equation is a statem ent of t he assump tion that there is no accum ulation of mass at the interface Th e energy balance for the vapor phase is
and for the liquid phase the energy balance reads rf)L,Hf -Lj-1Hf-,
Qf
FfHfF
(9
Here, represents the n et loss or gain of energy d ue to interphase transport. This term m ay be defined by
&r
JeYdal;
Jefdaj
and to liquid weeping through the tray floor (in these last two cases the m aterial balance relations would have to be modified). However, a comp rehensive model for estimating mass transfer rates in each of the possible flow regimes does not exist at present, and simpler approaches are needed. It is our purpose in this section to discuss he determ ination of the interphase ransport rates in a fairly general way. We consider it more appropriate to describe the precise form the rate relations may take in Part where the model is applied to some specific separations. Mass transfer is a rate process driven by gradients in concentration (m ore precisely, by g radients in chemical potential). Thus the rates of m ass transfer in th e vapor phase will depe nd on the differen ce between the bulk vapor compositions y; and th e vapor compositions at th e interface, yil (Figure 2) Similarly, the rates of mass transfer in the liquid phase will depend on the liquid composition at the interface, .z:~,nd on the bulk liquid composition, x$ (Figure 2). Many other factors (such as temperature, concentration , th e physical prop erties of th e fluids, system geometry an d h ydrodynamics, and the mass transfer rates themselves) influence th e actual rates that are achiev ed (Sherwood et al. 1975, p. 148). Current practice is to lump the effects of most of these variables into mass transfer coefficients and write NY
(11) Transport Relations
It is worth emphasizing that Eqs. 1-11 hold regardless of the models used to calculate the interphase transport rates J V ~ , N $ ~ & ~ ith , &a~mech . anistic model of sufficient complexity, it is possible, at least in principle, to account for mass transfer from bubbles in the froth as well as to entrained d roplets
k,"(yY
kf(x:
y!);
1,2.
f)
(12)
Rate equation s of this form a re used in th e two-phase models o Wag goner an d cowo rkers, f R icker et al., and of Billingsley and Chirachav alla. Rate equations of this form a re consistent with the material balance equations for the interface region only when all th k," have the sam e value and all the are equal (but not necessarily equal to k:). Under these circumstances the component efficiencies would all be equa l on any stage (they may , f course, vary from o ne stage to another) and there would be no point in writing this paper. more rigorous treatment of s transfer in multicomponent systems will lead us to the rate equations (Krishna and Standart, 1979):
c-
NY
(10)
where s the en ergy flux at som e particular point in the dispersion. An energy balance around the interface yields
AlChE Journal (Vol. 31, No 3)
interface.
Liquid sidestream
c-
kz(yr
kk(xL
zE)
yYNP
1,2..
(13)
rFNf.
1,Z..
(14)
where k & , k $ are multicomponent mass transfer coefficients. We shall assume that methods for obtaining these coefficients are available (see Part 11). If, as is the case with the m ore rigorous mass transfer models, the coefficients an depend on the mass transfer rates themselves, then we might represent the flux equations (1 an 13 by implicit an d non linear relations f th e general functional form NY
N:(k;,yf,y{,Nr,
1,2.
NF
Nt(k$,xf,xE,NE,
1,2.
C)
March, 1985
1,Z..
(15)
1 , ~ .
(16)
Page
It should be recognized that only of Eq. and of Eq. 16 can be written down. The flux f component is determined from the first fluxes and from the energy transfer rate equation as shown below The local energy flux, is made up of a con ductive heat flux and a convective contribution du e to the transport of enthalpy by interphase transport (Bird et al., 1960 , p. 566):
where are the partial molar enthalpies of species The conductive heat fluxes, q, are driven by temperature gradients in the fluid
h v ( ~ v Tr);
~ L ( T I
(18)
where h v and are heat transfer coefficients, which depend, among other things, on the mass transfer rates themselves. f we assume that these coefficients are available (from a suitable correlation or theoretical expression perhaps), then the local energy fluxes may be represented by eV(hV,TV,T‘,N[)
(19)
Interface Model
eL(hL,TL,Tr,Nf)
(20)
We adopt the conventional model of a phase interface: a singular surface offering no resistance to transport and where equilibrium prevails. The usual equations of phase equilibrium relate the mole fractions on each side of the interface:
The calculation of the total mass and energy transfer ra tes, ij and requires t he integration of the point flux relations (14-17) over some model flow path. (This, of course, is what’s done in the efficiency models, but with a different result in mind.) Many different models could be discussed here, but we shall confine ourselves to some of the simplest. In order to simplify the integratio ns required by Eqs. and 10, we assume that the interface state is the same throughout the dispersion on any stage1. (All other two-phas models make the same assumption, usually without explicitly stating it. ) We furt her assume that the mass transfer coefficients can be considered constant on any stage (again, this is the usual assumption). Then, by imposing a particular shape on the bulk phase composition profiles we find that the integrated total transport rates are equal to t he average fluxes multiplied by the total interfacial area,
where the mole fractions q$ and ?ikand temperatures, nd represent the integrated average bulk phase conditions I/ the bulk composition is assumed to be constant throughout the dispersion, then the average mole fractions simply are equal to the mole fractions of the streams leaving the stage: uij/Vj,?6 Zij/Lj. On the other hand, if the bulk composition varies linearly between the e ntering and leaving values, then th e average mole fractions are th e arithmetic averages:ij$ l/Z(uij/Vj u , , j + ,/V,+ I),?$ l/z(ltj/Lj Z i , j - l / L j - l ) . Finally, if we impose an exponential profile on the bulk mole fractions, the average composition is the logarithmic average of the entering and leaving values. Average temperatures can be calculated in a n analogous fashion. In our simulations f laboratory and pilot-plant-scale tray columns we use the arithmetic average composition for the vapor phase (the logarithmic average gives virtually identical results) and the outlet composition for the liquid phase. This correspondsalmost exactly with th e conventional model of plug flow of vapor throug a well-mixed liquid (see, for example, King, 1980, p.615). In our simulations of small-scale wetted-wall columns we use the ar ithmetic average for both phases. These simple models give good agreement with the experimental dat a (Part II), better even than some more complicated models. In writing Eqs. 21-24 have combined the interfacial area
Page
term directly with the heat and mass transfer coefficients. This is because many correlations (e.g., he AIChE method for bubble-cap trays) give the coefficient-area product. It is, of course, possible to use separate sources for these quantities, as would be done in a simulation of an operation in a wetted-wall column for which the interfacial area would be known fro m the geometry of the syste and the transfer coefficients only obtained from a correlation. Zuiderweg (1982) presents separate correlations for and for distillation on sieve trays. Notice also that Eqs. 21-24 are implicit in the mass transfer rates but not in the energy trans fer rates; there ar of Eq. 21, f Eq. 22, and one each of Eqs. 23 and 24. All this is done so that the rate equations (21-24) need be solved only once for each stage The assumptions that were made above may be relaxed if desired by furt her dividing each phase into a number of regions and writing separate balance and rate equations fo each region. Th e penalty is a large increase in the number o equations to be solved with a corresponding increase in the cost o obtaining a solution.
March, 1985
&j
(25)
:j
-cx:j-l=o
(26)
Kij(x$j&,T;,Pj) are the equilibrium ratios defined where in the usual way. We have already mentioned the reasoning leading to our considering the interfacial state to be uniform throughout the dispersion.
Variables and Functions
a Single Nonequilibrium Stage
Given the s tate of all feed streams, the flow rate of all side streams, heat loads, and pressures on the stage there is a total o 6c unknown quantities for each stage These are the component vapor flow rates ( v i j : in number), the component liquid flow rates (Zjj: c), the vapor temperature ( T y ) , he liquid temperature ( T f ) , the interface temperature ( T i ) , he vapor composition at the interface (y:,: c) he liquid compositionat the interface ( x t : c) he mass transfer rates ( N ; : ), and the energy transfer rates (Cv,Cf.).Th 6c independent equations that permit the calcufation of these unknowns are as follows: component material balances for the vapor ( M G : component material balances for the liquid ), component material balances around the interface ( M i j : ), he vapor phase energy balanc ( E Y ) , he liquid phase energy balance ( E L ) , he interface energy balance ( E i ) , the inte rfac e equilibrium refations(Q:j: ), the summation equations ( S y ; S,”), he vapor phase mass transfer rate eq uations (Eq. 21)(c he liquid phase mass transfer rate equations (Eq. 22) (c I) and t he energy transfer rate equations (23 and 24). The nonlinearity in these equations stems from th e presence o the values and enthalpies, as wel as the mass and energy transfer rate terms. reduction in this rather large nu mber of equations and unknowns can be obtained by eliminating certain equations that are simple linear combinations of variables. Recognizing that there is strictly independent transfer rates, the really on1 one set N,j(=N{= ay, we can eliminate the from Eq 21 and th from E q. 22 using the balance equations (Eq. 7) .Ntj
N b ( k ~ ~ j , y f j , j j r , , T Y , T ; , ~ k1,2. ~, 1,2..
C)
(27)
. c)~ j , N k ( k ~ ~ j , x b , x ~ ~ , ~ ~ , T1 ,;2, .N 1,2..
(28)
AlChE Journal (Vol. 31, No.
The energy transfer rate equations (Eqs. 23 an 24) are substituted into the interface energy balance to give &Y(hYuj,~TY,Tf,yK,.Nkj) j.
i;L
TI ,xkj,.Nkj)
(29)
Note that only mass transfer rates N,, appear in Eqs. 27-29; th transfer rates are eliminated from the lis (o Nz), an f independen t variables. Th e balance equations ( M f , :C) and the energy transfer rate equ ations (Eqs. 23 an 24) are removed from the set f ind epen dent function s. Two m ore variables xi an can easily be elim inated from the set of variables. They can be co mpu ted directly from th summation equations ( S y , S i ) , which, therefore, are dropped from the list f ind epen dent equation s. Th e final set of 5c indepen dent variables per stage is ordered into a vector (Xj) as follows
(olj,ozj
(XjIT
Ucj7T7,Tf,l1j,hj
.
&,>Y;j>dj,
kj,.NljA”j, I Yc-l,j,TjJljJzj~ E-1.f)
independen t equations corresponding to this set variables are ordered into a vector 5c
(sj)T
(M;,M;.
.
M ; , E ~ , E ~ , M ! ~ , MM~ $~, R. ; , R & .
Rc-i,j,Ef,Q!j,QLj. QEj,R!j,Rkj.
. Rk-1.1)
This concludes the form al developm ent of the mod el. It remains to show how what we refer to as the ME RQ equations (Material balances, Energy balances, Rate equations and Eq uilibrium relations) are solved and, in Part to present som e f our results.
SOLVING THE MERQ EQUATIONS
Methods of solving the eq uilibrium stage MESH equations fall into two groups; tearing methods (w here subsets the complete set f equatio ns are solved in seque nce) and sim ultaneou s correction (SC) metho ds (whe re all f th e equations ar e solved simultaneously). The sam e classification applies to the methods used to solve the various nonequilibrium models that have been devised. Waggoner and coworkers, for example, adop t a tearing strategy. Actually, they, as d o Ricker et al. (1981),consider their stages to be in thermal and mechanical equilibrium, and so they need only be concerned with solving what we might call the MRQ equations. Waggoner’s algorithm cannot be used whenever th e mass transfer rates are coupled, since the p hase balances for component depe nd o n the com position f all compo nents. Ricker et al. use Newton’s method to linearize all of the equ ations describing their stages. Billingsley and Chirachavalla (1981) use a collocation method to solve the differential equations that model their system. We have tried several different approach es to solving the MER equations and fou nd solving all of the eq uations simultaneousl using Newton’s method (or one its relatives) to be most effective. Newton’s Method and Related Fixed-Point Algorithms
(3(X)) enote an y collection f eq uations, nonlinear in ( X ) , he co llectio n f variables, w ith solution (0).The direct pred iction New ton correction is given by a solution of the equa tions linearized abou t the current estimate (Xk) of [Jkl(Xk+l
where
[ J k ] is
Xk)
-(3(Xk))
(30)
the Jacobian matrix with elements Jfj
d3f/dXj
(31)
Newton’s method is very v ersatile and qu ite robust; if it has a weakness, it is the com putation of [I]. or most engineering apAlChE Journal (Vol. 31,
o. 3)
plications f reason able com plexity, comp lete derivativ e information is rarely available in analytical form. In the present case, “unavailable” derivatives include derivatives of thermodynamic and transport properties with respect to quantities like temp erature, composition, and mass transfer rates (e.g., bKlj/dzij, dHf/dxf;i, dk$/dNkj, etc.). Finite-d ifference ap proxim ations of these derivatives can greatly increase the cost of s olution because many more physical property evaluations are required. Further, neglect f these derivatives is not generally jus tified; it can in crease the num ber of iterations or even cause failure. One way to avoid the repeated calculation of [J is to use a “quasi-Newton” m ethod (see Dennis and Mor6,1977).Broyden’s (1965) method , for example, computes an approximate Jacobian [I; ]the formula rom
where (2’) ( 3 ( X k + l ) Xk). Schubert’s (Xk+ (Xk)) an (S (1970) method is extension of the B royden meth od that preserves any known sparseness [J]. The quasi-Newton corre5tio is given by a solution the linear system (Eq. 30) with [ J k ] replacing [J
Unfortunately, quasi-Newton method s are not w ithout problems. Th e first is that they are not scale invariant (Newton’s method is and may, as a result, perform poorly on problems that are ill-conditioned. This is an impo rtant consideration in the present case since material and energy balance equations as well as transfer rate equations ar e part of the mod el. This particular set of eq uations is qu ite b adly scaled: I SI units are used then the MERQ equations include terms typically ranging over 10 orde rs f m agn itude [0(105) O(10-5)]. second prob lem so metim es arises because of the ne ed to supply an initial approximation to [J]. poor first approximation can increase the num ber iterations or even cause failure.
Hybrid Method
judicious comb ination of New ton’s metho d a nd a quasiNewton u pdate can, as shown recently b y Lucia a nd cow orkers (Lucia and Macchietto, 1983;Lucia and W estman, 1983),be a very efficien t mean s of solving the kinds of e qua tions foun d in chem ical process models. In wha t they call the hybrid metho d, the Jacob ian matrix is divided into two parts, a compu ted p art [Ck]and an approximated part [Ak]:
Ck contains any partial derivatives for which analytical expres sions can easily be derived. [Ak] is ad e up of a ny terms that contain p artial derivatives that ar e difficult to derive analytically and expensive to compute by numerical differentiation (e.g., dKfj/dxkj). In their algorithm, [Ck] s computed in each iteration while [Ak] only is updated from a n initial approximation (the null matrix in our application s) using the Broyden or Schu bert updates. 30 with [ J k ] given The correction to (Xk) again is provided by by Eq. 33. stepwise procedure for implementing the hybrid method can be found in the papers by Lucia and cow orkers cited above. In an early stage f this work, we compared the performance of th ree m ethods (Newton’s method, with finite-difference computations the “unav ailable” derivatives;Schubert’s metho d; and the hybrid metho d) of solving the MER Q equations. Almost always, Newton’s method took fewest iterations (-5) with the hybrid method requiring around ten iterations. However, the hybrid method typically takes 50% ess time (som etimes much less than 50%, ever m uch more). Schubert’s method generally worked wel in these exam ples because we ca n usually provide very good initial guesses f the indep end ent variables. Occasionally, this is not possible, and then Schubert’s metho d has failed. Newton’s method and the hyb rid alg orithm were used w ith com plete success to solve the problems described in Part
March, 1985
Page
Flgure 3. Block irldiagonalstructure
the Jacoblan for an absoi
?r.
Figure
Block trldlagonal structure of the Jacoblan for a dlsilllatlonCOlU with an equlllbrium reboiler and condenser.
APPLICATIONS OF THE MODEL A Multistage Process: An Absorber Perhaps
Absorbers frequently are modeled as simple sequences of equilibrium stages (King, 1980, chapter 10; Henley and Seader, 1981, chapter 15; Holland, 1975, chapter 4, 1981, chapter ). If, instead, we choose to model th e column by a sequence o nonequilib rium stag es, the vecto rs f variables and functions corresponding to the entire column are
where the variables and functions have been grouped by stage, as done by Naphtali and Sandholm (1971) n their formu lation f t he MESH equations. W hen group ed in this way, the vector of stage functions, 5f1),depend s only on the variablesfor the three adjacent stages, an 1. Thus the Jacobian matrix has the familiar ic k tridiagonal structure shown in Figure 3. The linear system (Eq. 30) is easily and efficiently solved when [J s bloek tridiagon al using the m atrix generalization of th e well-know n Thom as algorithm (Henley and Seader 1981, p. 600, for example, give the steps necessary to implem ent wh at is just Gaus sian elimination of submatrices). Th e submatrices V j ] , [ j ] , [ j] are extremely sparse; nd [ W j ] re almost empty. T he structure of these submatrices can be derived by straightforward differentiation of the ME RQ equations . Th e results f this exercise and diagram s of the sp arsity patterns of the submatrices are given in the app endix. The terms that are assigned to the computed and approximated parts of th submatrices are indicated there as well. It is interesting to com pare our to tal o 5c variables per stage with the 2E variables per stage in th e formulatio n f th e MESH equations due to Naphtali and Sandholm (1971)(see also He nley and Seader, 1981, p. 595; Holland, 1981, p. 581). For example, for a ten-stage column separating a ternary mixture, the equilibrium stage model comprises a mere 70 equations compared to 16 equations for the nonequilibrium stage model. For a 50-stage, five-component separation process, the equilibrium stage model boasts 550 equations. This compares to our nonequilibrium stage model with no less than 1,300equations. Clearly, there is a considerable differen ce in the size f the prob lem bein g solved. How ever, we consider large size only a mino r drawback f this noneq uilibrium stage model. W ith suitable sparse matrix-handling com puter codes, the large size of the system can easily be accommodated. Multlcomponent Distillation In a Multi-Nonequlllbrlum Stage Column
Distillation is the single most im portan t process presently employed by the chem ical process industries.In this section we discuss Page
March, 1985
the a ppl ication to this process of the ge neral stage model deve loped above. The formulationspresented h ere are used in our simulations of some actual multic om pon ent distillations described in Part Consider, first, a distillation column equipped with a partial condenser an d a p artial reboiler. Nonequilibrium models of multicomp onent co ndense rs/reboilers have been dev ised (for a review of recent d evelopm ents in these areas see, for example, chapter of volume of the recently published Heat Exchanger Design Handbook, Schlunder, 1983).However, these m odels involve the num erical solution f differential rather than algebraic equations Althoug h it certainly is possible to use the results of these models to com pute the MERQ functions, we have not attempted to do this yet. In stead, we p refer, for sim plicity, to solve algebraic equations only (at this stage anyway) and, therefore, we mod el the partia condenser/reboiler as equilibrium stages, as do Waggon er and Loud (1977) n their m uch simp ler nonequilibrium stage model. 2c variables correspond ing to this set of equations are the comp onent flow rates 1, an uij, and the temperature, Tj (see, e.g., Henley and Seader 1981, p. 594).Th equations that m odel these stages include the component material balances and equilibrium relationships. Th e energy balance is replaced by one specification equ ation (e.g., distillate rate specified). total condenser/reboiler may also be described by a set o 2c equations: component material balances, a specification function (to replace the energy balance as for a partial condenser/reboiler), equations equating the mole fractions i the product and reflux streams, and an equation specifying the Tbub, the bubble temperature the prod uct streams (e.g., point temp eratu re f th e product). For a column mod eled in either way, the Jacobian matrix retains its block band structure, but the subm atrices [Wl], [ U z ] , W,-I], an are no longer square. Rather, they a re rectangular with 5c 2c 1, as indicated in dimensions 2c or 5c Figure 4. It should be recognized that the steps f t he matrix generalization of the Thomas algorithm apply unchanged to the triangularization f th e matrix in Figu re 4. Provision need only be made for multiplication of rectangular matrices. Distillation at Total Reflux
Th e maxim um separation possible in a given num ber of stage is attained at total reflux when the overhead vapor is completely condensed and returned to the column. This situation is also attractive for ex perim ental purpose s, as no material is lost from the column. Using a m balan ce aroun d any top section of the column it is easy to show that uij+
(31)
This suggests that eithe r the set of com ponen t vapor flows, u,,, or
AlChE Journal (Vol. 31, No
the component liquid flows lij, can be removed from the set of variables representing t he jt h nonequilibrium stage. The component material balance equations or M f ; . )would then be removed from the set of equations for that stage. We will thenbe left with 4c variables and functions per stage. Single Phase Control
In many separation processes the resistance to mass transfer lies predominantly in one phase. Distillation s a process that sometimes is controlled by the vapor phase resistance (Vogelpohl, 1979; Kayihan et a)., 1977; Krishna et al ., 1977; Dribika and Sandall, 1979). we assume that the liquid phase resistance can neglected, then the interface liquid compositio n is that of the bulk liquid and it is no longer necessary to consider the as independent variables. We are then left with 4c variabfes per stage, which, of course, means that t he number f equations must also be reduced by 1. As there is no resistance to mass transfer in the liquid phase, it is no longer possible to calculate the mass transfer rates in the liquid from the rate equations which, therefore, ar e dropped from the set f equations for the jt stage. would be eliminated For liquid phase control the yi and from the sets of variables and equations.
ij
S” SL
WI [V
[WI
Vapor-Phase-Controlled Distillation at Total Reflux
a vapor-phase-controlled distillation at total reflux, we have to deal with only 3c equations and variables per nonequilibrium stage: (~1)’
(01j.
ucj,T)’,TF,Nlj.
(M;. M z , E Y , E f , R ; . (31)~
Ncj>Yij.
Subscripts
component number stagenumber
YS-l,j>TjI)
R,V-l,j,E;,Qlj..
Qcj)
Superscrlpts
We note that this situation is encountered in many of the experimental simulations more fully described in Part
interface liquid bulk liquid feed vapor bulk vapor feed partial molar property approximation
Other Applications
The model developed in this paper can be applied to many other situations:distillation with chemical reaction, dynamic simulations, and distillation systems in which t he presence of two liquid phases are of particular interest. These particularcases require som modification and/or extension of the model equationse presented here (rather than simple elimination of equations and variables a in the examples considered above). We will consider these variations on theme in fu tur e parts of this series. ACKNOWLEDGMENT This paper is based on work supp orted by the N ational Science Foun dation through Grant Number CPE 8105516.
NOTATION
[A1 [C
[J
interfacial area (m2) approximated part of Jacobian number of components computed part Jacobian point energy flux (J/m2s energy balance function (J/s total interphase energy transfer rate ( J/s component feed rate (kmol/s) total feed rate (kmol/s) discrepancy function vector heat transfer coefficient (J/m%K) enthalpy (Jik mol) Jacobian matrix of partial derivatives mass transfer coefficient (kmol/m2-s) equilibrium ratio
AlChE Journal (Vol. 31, No. 3
component liquid flow (kmol/s) total liquid flow (kmol/s) mass balance function number of stages point molar flux (kmol/ m%) total molar flux (kmol/m2-s) total interphase mass transfer rate (kmol/s) pressure (Pa) conductive heat flux (J/m2.s) heat removal from liquid phase (J/ s) heat removal from vapor phase (J/s interface equilibrium function (dimensionless) ratio side stream to interstage flow rate relation functions summation functions side stream vapor flow (kmol/s) side stream liquid flow (kmol/s) temperature (K subdiagonal block submatrix diagonal block submatrix component vapor flow (kmol/s) total vapor flow (kmol/s) superdiagonal block submatrix liquid phase compositions vector f variable vapor phase compositions
Matrices
column matrix row matrix square or rectangular matrices
LITERATURE CITED Billingsley, D. S., an Chirach avalla, “Num erical Solution of Non equilibrium Multicomponent Mass Transfer Operations,”AIChE J., 7, 968 (1981). Bird, R. B., E. Stew art, and E. N. Lightfoot Transport Phenomena, Wiley, New York (1960). Broyden, C. G. Class of M ethods for Solving Nonlinear Simultaneou Equations,” Math. Comp., 1 9 ,5 7 7 (1 9 6 5 ) J. MorB, “Quasi-Newton Methods, Motivation and Dennis, J. E. an Theory,” SIAM Reo., 1 9 ,4 6 (1 9 7 7 ) Deiner, an J. A. Gerster, “Point Efficiencies in Distillation of Acetone-Methanol-Water,”Ind. Eng. Chem . Process Design Deve lop., 339 (1968). Dribika, M. M ., and C. Sanda ll, “Simultaneou s Heat and Mass Transfer for Multicomponent Distillation in a Wetted-wall Column,”Chem. Eng. Sci., 3 4 ,7 3 3 (1 9 79 ) Hausen, H ., “The Definition the‘ Degr ee Exchange on Rectifying Plates for Binary and Ternary Mixtures,” Chem. Ing. Tech., 2 5 , 5 9 5 (1953). Henley, E. J., an J. D. Seader, Egul&rium-Stuge Separation OperatZons in Chemical Engineering, Wiley, Ne w York (1981). Holland, C. Fundamentals and Modeling of Separation Processes, Prentice-Hall, Englewood Cliffs, NJ (1 975)
March, 1985
Page 455
Holland, C. D., Fundamentals of Multicomponent Dis tillation, McGraw-Hill, New York (1981). S. McMahon, “Comparison of Vaporizatio n EfHolland, C. D., and ficiencies with Murphree-type Efficiencies in Distillation-I,” Chem. Eng. Sci., 25,431 (1970). Kayihan, F., 0. C. Sandall, and D. Mellichamp, “Simultaneous Heat and Mass Transfer in Binary Distillatiort-11 Experimenta l,” Che m. Eng Sci., 32,747 (1977). King, C. J. Separat ion Processes, 2nd ed., McCraw-Hill, New Y o r k (1980). Krishna, R., H. F. Martinez, R. Sreedhar, and G. L. Standart, “Murphree Point Efficiencies in Multicomponent Systems”, Trans. Inst . Chem. Engrs., 55 178 (1977). Krishna, R. and L. Standart, “Mass and Energy Transfer in Multicomponent Systems,”Chem. Eng. Commun., 3,20 1 (1979). Lucia, ,, and S. Macchietto,“A New Approac h to the Approxim ation o Quantities Involving Physical Properties Derivatives in Equation Oriented Process Design,” AZChE 29.705 (1983). Lucia, A,, and K. R. Westmau, “L ow Cost Solutions to Multistage , Multicomponent Separations Problems by a Hybrid Fixed Point Algorithm,” Foundations of Comp uter -ai ded Process Design, 2nd International Conferenc e, Snowmass, CO (1983). Medina, A. G. N. Ashton,and C . McDer mott, “Murphre e and Vaporization Efficienciesin Multicomponent Distillation,” Chem. En g. Sci.,33,331 (1978). Medina, A. ., N. Ashton, and C. McDermott, “Hausen and Murphree Efficiencies in Binary and Multicomponent Distillation,” Chem . Eng Sci., 34,1105 (1979). Murphree, E. V., “Rectifying Column Calculations,” Ind . Eng. Chem., 17,74 7 (1925). Naphtali, L. M., and D. P. Sandholm, “Multicomponent Separation Calculations by Linearization,” AlChE 17,1 48 (1971) Schlunder, E. V. ed. Heat-Exchanger Design Handbook, Hemisphere Publishing Corp., Washington, D.C. (1983).
Schubert, L. K., “Modificatio n of a Quasi-Ne wton Method for Nonlinear Equations with a Sparse Jacobian,” Math. Co mp., 24,27 (1970). Sherwood, T. K., R. L. Pigford, and C. R. Wilke, Mass Transfer, McGraw-Hill, New Y o r k (1975). Standart, G. L., “Studies on Distillation-V. Generalize d Definition of Theoretical Plate or Stage f Contacting Equipment, Chem. Eng. Sci., 20,611 (1965). Standart, G. L., “Comparison Murphree-type Efficiencies with Vaporization Efficiencies,” Chem. En g. Sci., 26,985 (1971). Toor, H. L., “Diffusion in Three-Component Gas Mixtures,” AlChE 3, 198 (1957). Toor, H. L., “Prediction of Efficiencie s and Mass Transfer on a Stage with Multicomponen t Systems,” AIChE J. 10,545 (1964). Toor, H. L., and J. K. Burchard, “Pla te Efficiencies in Multicomponent Distillation,” AZCh J., 6,202 (1960). Vogelpohl, A,, “Murphree Efficienci es in M ulticomponent Systems,” In Chem. Eng. Symp. Ser. No 56, 2,2 5 (1979) Waggoner, R. C., and L. E. Burkhart, “Non-equilibrium Computation for Multistage Extractors,” Compu t. Chem. Eng., 2, 169 (1978). Waggoner, R. C., and G. D. Loud, “A lgorithms or the Solution f Material Balance Equations for Non-Conventional Multistage Operations,” Comput. Chem. Eng., 1,4 9 (1977). Zuiderweg, F. J. “Sieve Trays-A View f the State of the Art,” Chem Eng. Sci., 37,1441 (1982). Supplementar ymaterial has &n deposited as Document No 04262 with the National Auxiliary Publications Service (NAPS), c/o MicrofichePublications, 4 North P earl St., Portchester, NY 10573, and can be obtained for $4.00 for microfiche or $7.75 for photocopies.
Manusmipt received July 29,1983; eoisim received Dec. 6, and accepted January
28,1984.
Nonequilibrium Stage Model of Multicomponent Separation Processes Part
Comparison with Experiment
nonequilib rium stage model of count ercurrent sep aration processes is use composition and temperature profiles during binary and multicomponent distillation in wetted-wall and bubble-cap-tray columns. The profiles predicte by the model are compared with experimental data for the binary systems benzene-toluene, ethanol-wa ter, and acetone-chloroform a nd for the t ernar y systems benzene-toluene-ethylbenzene, acetone-met hanol-water, methanol-isopropanolwater, acetone-met hanolethano l, and benzene-toluene-m-xylene. he model does very commendable job of predicting th e composition profiles measured for these systems; average absolute differences between predicted and measured mole fractions ar e seldom greater tha mole percent and are often very much less LO predict
R. KRISHNAMURTHY
an R. TAYLOR Department of Chemical Engineering Clarkson University Potsdam, NY 13676
SCOPE In Part this series, a nonequilibrium stage model of countercurr ent multicomponent separation processes was described. The motivation for developing the mode l was the ob__
R. Krishnamurthy is presently with The BO Group, Murray Hill, NJ 07974.
Page 456
March, 1985
servation that t he equilibriu m stage model modified by the inclusion of a s tage efficiency that has th e same value for all components (the model presently used to simulate separation processes) is unab le to pr edict accurately th e composition profiles in processes separating multicomponent mixtures in which
AlChE Journal (Vol. 31, No. 3)