- Email: [email protected]
Process dynamics of binary distillation
Chemical Engineering Science, 1969, Vol. 24, pp. 1687- 1698.
Pergamon Press.
Printed in Great Britain.
Process dynamics of binary distillation V. J. POHJOLA and H. V. NORDEN Technical University, Department of Chemistry, Otaniemi, Helsinki, Finland (Fint received 16 September 1968; in revisedform 30 April 1969) Abstract-The dynamic behaviour of binary distillation in the vicinity of an equilibrium state is investigated. The genera1 linear state equations and the principles of solution are given. These are applied to a more detailed study of a theoretical nonlinear process for which the response of compositions to smal1 perturbations in the feed flow rate is solved. The effects of certain approximations needed to render the differential equations linear with constant coefficients are studied qualitatively. A feature in the present model which may be of interest is the dividing of the process into two subprocesses, mass transfer and liquid flow, and the combining of the independently obtained responses for each subprocess by superposition. As the procedure suggests the use of the Laplace transform, a simple numerical Laplace inversion method is also given. INTRODUCTION
THE AUTOMATIC control and optimization of processes and the advanced design of plants and equipment are al1 at their best based on detailed knowledge of dynamic behaviour of the process in question. In distillation, particularly, many difficulties arise when it is desired to predict how the process wil1 act as a function of time. There are phenomena the exact physical mechanisms of which are not known or there are no ways to characterize them simply enough. Such features of practical distillation as the large number of plates and components, the large number of process and the nonlinear relationships variables between them make the analytical attempts hopeless to solve synthetic models, and render the solutions on digital computers time-consuming and expensive. Though computers are becoming faster and capable of more precise simulation, it is felt by several authors [ 1-31 that there wil1 be a growing need of relatively simple short-cut methods for the purposes mentioned above. Several reviews on the vast literature covering the control and process dynamics of distillation have been made [4-61. The secondary effects in distillation, by which are meant physical phenomena other than mass transfer which may have some
effect on the response of compositions, are discussed by Rosenbrock [7]. He discusses also certain distinctive problems in distillation from automatie control’s point of view [3]. Articles dealing primarily with the mathematics involved are reviewed by Mah et al. [SI. This paper is concerned with a theoretical study of composition response of a nonlinear binary distillation, where the hydraulic lag, the most common secondary effect, is present. Starting from genera1 dynamic state equations it is arrived finally at a linear readily solvable model which is expected to describe the response in the neighbourhood of the initial steady state with an accuracy sufficient for control purposes.
GENERAL
MODEL
A linear mathematica1 model of binary distillation wil1 be developed in genera1 terms. A matrix notation wil1 be used to expose the basic structure of the model and to make it easier to follow the-manipulation of equations. Particularly the response of compositions to perturbations in the feed flow rate wil1 be studied. State equations
The dynamics of binary distillation, a typical multivariable process, can be described by a
1687
V. J. POHJOLA
set of ordinary first-order generally nonlinear, y
= i(t)
and H. V. NORDÉN
differential equations,
= f [a(t),p(t>].
(1)
Here a is a state vector and it includes al1 the variables that are needed to specify the state of the process. p is a vector of input variables. Consider an uncontrolled distillation process the state of which in certain range of conditions can be described by its compositions and liquid flow rates. Let x = [x,] be a column vector of compositions and 1 = [IJ be a column vector of liquid flows, the elements of the vectors being deviations from certain initial steady state values of compositions and flow rates of the liquid streams leaving plates n and m, respectively. (n and m may get values from 0 to N + 1 where N is the total plate number of the distillation column; 0 refers to the condenser and N + 1 to the reboiler.) Let the elements of p be deviations from certain initial steady state values of compositions and flow rates of the streams entering the column. Suppose that at least in some vicinity of an initial steady state the functional dependence of Eq. (1) can be approximated linear. It can be written then more explicitly
[‘i]=.[;]+RP,
=A[;]+K
where 0 signifies a zero matrix. In the same way b=
;.
[l
(5)
Thus
or written in an expanded form X=Px+QI+sf
i=m+tf.
(7)
03)
(2) The solution of state equations
where A and B are matrices which are in genera1 dependent on time. To determine the response of the process to certain input p requires the solution of Eq. (2) with a set of initial conditions. Suppose next that there is only one input variable and that this is some perturbation in the feed flow rate, J? The state equation can be written then [:]
complexity only the time-dependent elements are thought to be included in 1. The major simplification that can be made at this stage is to assume that the changes in liquid flows do not depend on changes in compositions. This assumption is justified if there are no considerable changes in latent heats during the transient, and if the changes in compositions do not affect the rate of heat transfer in the condenser and reboiler. With this simplification A can be partitioned into submatrices as follows
(3)
where b is a column vector. To avoid any extra
The solution of Eqs. (7-8) in the case of timevarying coefficients is tedious. A step-by-step integration procedure is suggested by Mah et al.[8]. We concentrate here, however, on a situation where the assumption of constant coefficients is justified. Such a situation arises when the deviations in the state variables are so smal1 that the coefficients practically are constant, or when the nonlinearity of the original problem is so strong that Eq. (2) is for any practica1 purpose an equally poor approximation of Eq. (1) whether the dependence on time of the coefficients is taken into account or not. In the case of constant coefficients the solutions are
1688
Process dynamics of binary distillation
XQ) =
J;ew
I(t) = ILexp
P(~-T)I[Q~(T)
+s.(T)ld~
[R(t-T)]tf(T)
dr.
(9) (10)
The evaluation of the vector convolution integral of Eq. (9) wil1 be examined in the following. The vector Ql sfcan be written in the form
+
[lf1
where
Z is a constant
-f ,
matrix.
If in addition
is denoted for brevity by v, Eq. (9) becomes X(t) = 10 eXp [P(?-T)]
zV(T)&.
(12)
Consider temporarily a solution where the hydraulic lag is neglected and f is a unit step change. Then each element of 1 is also a unit step change. Let this particular solution beu: t
exp [P(t-T)]
u(t) =
1 Z 1 dT.
0
(13)
Ijl
111 Now the expression du(t-7) dt
x(t)
=
du(t--7)
vtT)dT
dt
Thus it has been shown that the original problem to determine the response of compositions to smal1 perturbations in the feed flow rate can be divided into two sub-problems: the determination of composition response to a unit step in the feed flow rate neglecting hydraulic lag, U, and the determination of liquid flow response v to the actual changes in J This procedure allows to treat the mass transfer and the liquid flow phenomena separately. Further the two responses can be superposed by (17) which is valid. Let f be any function of time.
il !
MODEL FOR PROCES
A PARTICULAR
In the following the composition response of a theoretical binary distillation to perturbations in the feed flow rate wil1 be studied in detail. The unsteady state material balances around the reboiler, condenser and each plate wil1 be linearized and solved along the main lines brought out in the preceding chapter. Dejinition of the process
1
= exp [P(t-i)]Z
DETAILED
(14)
1
The pr.incipal nomenclature of the process is given in Fig. 1. The column is assumed to operate
represents the response of compositions to a unit impulse in the feed flow rate when the hydraulic lag is neglected. It can be written also in the form dU(iyT)
= eXp
[P(t-T)]
z.
(15)
The relation between the elements of u and U is u,=
(17)
(11)
[l
Ql+sf=Z
Substitution of (15) into (12) gives the final composition response in terms of U and v:
E Unk, k=l
(16) Fig. 1. The nomenclature of the distillation process.
where K is the total number of elements in v. 1689
V. J. POHJOLA
and H. V. NORDÉN
with a partial reboiler and a total condenser. The following additional assumptions are made: (1) Constant molal overílow on the plates in the steady state. (2) Pressure changes in the column during the transient are negligible. (3) The vapor hold-up in the column and the liquid hold-up’s between the plates are negligible. (4) The liquid hold-up’s on the plates are perfectly mixed. (5) The feed and the reflux enter the column as a liquid at its boiling point; any additional condensation or vaporization due to potential temperature differente between the entering liquid and the liquid on the feed plate or the top plate is negligible. (6) The vapor leaving the reboiler is in equilibrium with the bottom product. Thus no enthalpy balances wil1 be considered.
Material
where the smal1 letters h, 1, p,f, x and y are used for the deviations from the initial steady state AH, AL, AP, AF, AX and AY, respectively. As no changes occur in the vapor flows, the deviations of the VS are zero. In addition the vapor flow is constant throughout the column, so the subscripts of the vapor flow terms wil1be dropped A total material balance is dH n = Lpl dt
%
(20)
+fn.
(21)
= In_1- (l,+p,)
Combining ( 19) and (2 1) we obtain (H,O+ h,) %
= (Li-1 + In-&-I - (L,o + p,o+ l,_, +fn)x,
A material balance of the more volatile component around a typical plate in the unsteady state is
d(Hdf”)= L._,X,_,
- (L,+
P,)X,-
v,r,
+J’n+~y,+~+Fttx~.
(18)
Here P, is a product stream from plate n and F, is a feed onto plate n. (Po = D, PN+l = B, and FM = F. P, = F, = 0 otherwise.) Multiplying both sides of (18) by a forward differente operator, A, an expression is obtained where the variables are deviations from the initial steady state and in addition the constants are defined in the initial steady state and can be assumed known. Eq. (18) becomes (X,O+x,)
%
= (LE-1 + ln_l)X,_l - (L,o+ 1, + P,o+p,)x, - v,“y?z + C+1Y,+, --KI”(h+P,)
P,)+ F,,
or
balances
(H,O+ h,) %+
- (L, +
+ x-ll,-1
+xF”fn7
(19)
+ I/oY,+, + cx;+
- vy,
- x,o)I,-1
+ (X,O - XnO)fn.
(22)
Equation (22) is a differential equation with time-varying coefficients. Its nonlinearity is determined by the nonlinearity of the dependence of y upon x and 1. A straightforward means to render the coefficients constant is to neglect the products of differences in (22). This is based on the assumption that the deviations from the initial steady state are small. With this simplification Eq. (22) becomes H,’ %
= L:_~x,_~ - (L,O + P,“)x, - voy,
+ WY,+, + (X;-I
-X,t”)f,-~
+
(XF” - X,O>f,. (23)
The dependence of y on x and possibly on 1 is determined by the phase equilibrium relationship and the Murphree plate efficiency. An explicit expression is obtained when the efficiency for the liquid phase and the equilibrium curve in the form Xn =
1690
Proces
E
xn-1-&_
=
dynamics of binary distillation
xn-1-x,
n x,-,-x:
(24)
- x,-1-qY(
Rearranging and multiplying by A we obtain
where again the products of differences have been neglected. Eq. (32) can be written more briefly Yn = %“&l-I
&(Y,z) = (l-B,z”-Sn)xa-l+
(E,O+&)xn
- (XL - XnO)Sn, 1 where Z = -
(25)
[=AEandx=AX.
E’
AS
Ady,) - AYn = cp’tyn”), where Y,,mis defined by Y”O=z Y,” S Y?l,
(27)
Acp(Y,) =
(28)
we may put
where y = A Y. If the equilibrium curve is available in the form YZ = $(X,), the following relationship may be used
where oo, po and yo are constants defined in the initial steady state. It is interesting to note that after neglecting the products of differences it is the behavior of the plate efficiency which determines the nonlinearity of (33) and consequently of (23). Unfortunately there are no data at hand on how the plate efficiencies do behave in transient conditions. Therefore no attempt is made here to replace ,$ by some linear function of x and/or 1. Rather we content ourselves to assume the plate efficiencies to remain constant. Thus (33) reduces to y, = (Y,OX,_l+&OX n.
H = F%+zv.
h?h l+X&,-E,
5oco ei51 Cl a2
(30)
6,
x0
x1
PX=
x2
c2
.
where
.
.
.
9(36) .
h
=E,
n
L,
=IY,*-Kt
n
L,X,-X,*’
$7ö)[( 1-
Yn =
n
-
aN
(31)
Here en is the vapor efficiency and m, is an approximation of the slope of the phase equilibrium curve on plate n. Combining (25) and (28) and solving for y, we obtain finally
(35)
The first term on the right-hand side is
There is also a simple relationship between the Murphree plate efficiencies for liquid and vapor phases: =
(34)
After combining Eqs. (23) and (34) the material balances of the more volatile component for n = 0 . . . N + 1 can be expressed as a set of linear differential equations with constant coefficients
(29)
E,
(33)
+ PRox, + Yn0L
. bN
CN
x,
aN+l
bNfl
xN+l I
where
I
a, =
-L (Lz-, - V%xno) Hn”
(37)
b,=-~[L.O+Pr+r(-a~+,+B.“)]
(38)
c, = J-
(39)
ZnO)x,-1+ ~,Ox,
Gc-l-&“)~nl~
(32) 1691
Hn”
vop:,, .
V. J. POHJOLA
and H. V. NORDÉN
It is to be noted that L-, = LN+l = 0, (Y,,= & = 0, CI,,,+~ = &+2 = 0, and as the plate efficiency in the reboiler is equal to unit y (assumption (6)), (llN+l= 0. The second term is 0 ____---------_
-_’
zv= ZM
>
Response of compositions to a step change in the feedflow rate neglecting hydraulic lag It is obvious from the expression of the final composition response, (17), that most of the advantages of the model presented herein are taken by the use of the Laplace transform. Giving the elements of v in Eq. (35) a constant valuef and taking the Laplace transform with respect to time, Eq. (35) becomes
(40)
1.M
Z h4+1
(SI-P)X Z N+l
= zs-’
(46)
f,
1N
where (41) Z,=$,(X;_l-X.o),
M < ut < N+l.
(42)
As the changes in the feed flow rate,f, affect the liquid flows below the feed plate only, the liquid flows in the enriching section remain unchanged. So the elements from 0 to M - 1 have been omitted in 1. In the total material balance, (2 l), the relation between h and 1 must be known. Suppose [9] that l,, = &h,,
(SI - P)ü = zs-1.
(47)
LiquidJEow response Equation (44) is Laplace-transformed respect to time to give
(43)
where K, is a constant. The total material balances for 12= M . . . N can be written now i = Kl+kf,
where the bar signifies the transformation and s is the transform variable. Equation (46) is a set of linear algebraic equations with a tridiagonal coefficient matrix and is easily solved by Gaussian reduction or by Crout’s method, a modification of the former. The solution contains s as a parameter. An expression analogous to ( 15) is
(SI - K)Ï = kf:
(48)
This can be solved analytically. we obtain
(44)
1, =
with
1
On plate n
(49)
f,
ifiM(&-‘s + 1)
where
(-KM)
KM
i 1M+ 1
KV+1(-K44+1) Kw+z(-&,+z ) KI + kf=
0 0
+* . I
KN(-KN>
1692
f.
(45)
Process dynamics of binary distillation
where f is the transform of an arbitrary bance in the feed flow rate.
distur-
Superposition
Equation transformed
(17) can be written notation simply x = SÜY.
in a Laplace-
(50)
a special emphasis is put on a universal compact formulation of the equations.
and
AcknowledgmentsThe authors wish to thank 0. A. Asbjernsen, Dr. Techn., for valuable criticism. Also, thanks are due to the Computing Center of Helsinki University for giving free machine time for calculations on an electronic computer IBM 1620.
For example on plate 12: NOTATION K x,
=
s
I: k=l
Unrl&.
Laplace inversion
The final composition response is obtained by the Laplace inversion of (50): x(t) = L-’ Z(S) . 1 1
coefficient matrices state vector coefficient column vector elements of P bottom product flow rate distillate flow rate Murphree liquid efficiency El F F.U feed flow rate onto plate n F, AF AF, f relationship in vector ; functional notation H5 weight in Eq. (53) liquid hold-up on plate n H, hz AH, 1 unit matrix J, K integers K average of time constants for flow response time constant for flow response, K, defined in Eq. (43) K coefficient matrix, defined in Eq. (45) k coefficient column vector, defined in Eq. (45) flow rate of liquid leaving plate 12 L 1, AL, 1 = [[“l column vector of liquid flows L Laplace transform operator M feed plate of the slope of the m, approximation equilibrium curve on plate n N total number of plates zero matrix 0 P,Q,R coefficient matrices A, B
(51)
a b a,, b,, c, B D
(52)
Nordén [ 101 has developed a simple numerical inversion formula
Gt(ti)
=
i
C
HjGt(sj-i)9
(53)
*3=1
where the arguments are chosen from geometric series ti = tod
(54)
and s._.=!&!d-i 3 I to
The H;s are weights independent of t. v is a constant, slightly larger than 1. qo and t. are suitably chosen constants. A detailed discussion of the method is given in [lol. DISCUSSION
The present paper reports some of the fundamental results of a larger research program of the authors and some co-workers directed towards the process dynamics and computer control of continuous distillation. In the text 1693
V. J. POHJOLA
and H. V. NORDÉN
flow rate of product stream from plate n PZ AP, P input vector s Laplace transform variable S coefficient column vector t time variable t coefficient column vector constant in Eq. (54) u = rz7.k: composition response matrix which neglects hydraulic lag u = [4ll U in column vector form V vapor flow rate V, V liquid flow response v = [Gtl augmented
-
p,
vector =
Superscripts 0
- -
fractional response of composition fractional response of composition Wl to a step change in the feed flow rate when the hydraulic lag is neglected, defined in the Appendix fractional response of liquid flow, w2 defined in the Appendix liquid composition on plate n X?l &I A-Kl [&l column vector of compositions Y?l vapor composition on plate n Yfl Ar,, Z coefficient matrix W
x =
f1 11
Greek symbols % Cl- %)/
m-l-x,“)/d(y,)
forward differente operator, deviation from the initial steady state value Murphree vapor efficiency constant in Eq. (55) dimensionless time variable = &2°1H,o) t constant defined in Eq. (3 1) constant in Eq. (54) LI& A? dummy time variable phase equilibrium relationships
f *
Subscripts B D F M n P
initial steady state value final steady state value phase equilibrium composition
bottom product distillate feed feed plate plate, numbering to reboiler product stream
condenser
Other notations A bar signifies a Laplace-transformed variable A prime is used for a derivative Matrices are denoted by bold-face capita1 letters whereas lower case bold-face letters are used for vectors
REFERENCES
[ll SARGENTR.
[2] 131 [4] [5] [61 [7] [8] [9] [lO] [l 11
from
W. H., Chem. Engng Prog. 1967 6371. LEFKOWITZ I., Chem. Engng Prog. Symp. Ser. 1963 59 178. ROSENBROCK H. H., Chem. Engng Prog. 1962 58 43. ARCHER D. H. and ROTHFUS R. R., Chem. Engng Prog. Symp. Ser. 196157 2. ROSEBROCK H. H., Trans. Instn chem. Engrs 1962 40 376. WILLIAMS T. J., Chem. Engng Prog. Symp. Ser. 1963 59 1. ROSENBROCK H. H., Trans. Insm chem. Engrs 1962 40 35. MAH R. S. H., MICHAELSON S. and SARGENT R. W. H., Chem. Engng Sci. 1962 17 619. SMITH B. D., Design ofEquilibrium Stage Processes. McGraw-Hill 1963. NORDÉN H. V.,ActaAcod.Abo. Mathematicaet Physica 196122 1. NORDÉN H. V. and SEPPÄ I.,Actupolytech. stand. Ch. No. 36, 1965.
1694
Process dynamics of binary distillation
APPENDIX It has been shown that the response of compositions to perturbations in the feed flow rate can be obtained by superposing the response of compositions to a unit step change in the feed flow rate when the hydraulic lag is neglected, and the response of liquid flows to the actual perturbations. To arrive at this some simplifications have been made which concern the balance equations of the more volatile component. Thus it might be of some interest to solve numerically this part of problem, that is Eq. (46), to get some idea about the effects of the approximations which led to Eq. (34). As the exact solutions of the original nonlinear equations are not known, we content ourselves to compare the response curves from (46) on large values of time with the final steady state values of changes of compositions which are computable. Giving the deviation terms appearing in the coefficients of (22) their final steady state values one ought to get solutions which approach asymptotically the final steady state values of the changes of compositions if (34) is valid. If, however, the plate efficiencies in the final steady state differ from those in the initial steady state, these asymptotic values deviate from the truc final values by amounts which are related to the errors caused by neglecting Sn in (33), (and to a less extent to the errors which are caused by neglecting the changing of the slope of the equilibrium curve, VA, dur+g the transient). The numerical values which wil1 be given next have been chosen somewhat arbitrarily. However, the initial and final values of the plate efficiencies have been chosen so that the error in yn caused by neglecting e. in (33) is about the same magnitude and into the same direction on every plate. This can be regarded as the most unfavorable situation as no cantelling of errors occurs. It is emphasized that no genera1 conclusions from this kind of study can be drawn and al1 what follows is meant purely as qualitative.
Table 1. The constants used in computing the initial and final steady state compositions F” = 100 kmol/hr f= 5 kmol/hr Do = 50 kmol/hr Loo/Do x,0 N M
= = = =
3.75 0~500 50
optimum in the initial steady state (= 24) H,O= 1 kmol, n= l...N Ho0 = 5 kmol HN,, = 5 kmol
h,/kmol = 0.004 I,,/kmol/hr,
n= M .. N
h ,V+1-0 Y* = 1.5x-0.5x
l=-X+0.8, 0~ XS 0.2 0.2 s X =S 1 1E = -0.25X+0.65,
Response
curves
The response curves for condenser (n = O), reboiler (n = 51), feed plate (n = 24), and plates 6, 15, 33, and 42 are given in Figs. 2-8. The notation is as follows: 1
Steady state
Taking into account the assumptions made earlier any steady state of the process can be defined by fixing six process degrees of freedom[9]. The six variables chosen in this study are: feed flow F, feed composition X,, distillate flow D, reflux ratio LdD, total number of plates N, and feed plate location M. Instead of giving M a numerical value it is required that the feed plate location be optimum in the initial conditions. In addition, the following quantities have to be known: the phase equilibrium curve, the plate efficiencies for each plate in the initial and final steady states, the liquid hold-ups in the reboiler, condenser and on each plate, their dependence on the liquid flow rate, and the magnitude of the step change in the feed flow rate. The numerical values are given in Table 1. The analytic expressions for the phase equilibrium curve and for the dependence of the Murphree vapor efficiency on composition are assumed to be valid both in the initial and in the final steady states. The steady state compositions are obtained by iterative plate-by-plate calculation. The optimal feed plate location is attained in the initial steady state with M = 24. This value of M maximizes the differente XDo- XBO.
1695
CESVd.?.4NaIl-D
log 0 Fig. 2. Composition response of the bottom product, w,(0), and flow response of the liquid stream entering the reboiler, w,(0).
Fig. 3. Composition response of the liquid stream leaving plate 42, w,(fI), and flow response of the liquid stream entering plate 42, w*(B).
V. J. POHJOLA
and H. V. NORDÉN
Curves a. Solutions of Eq. (46) with the elements of P evaluated giving the time-varying coefficients of Eq. (22) their final values instead of the initial ones appearing in Eq. (23). Curves b and c. Solutions which are obtained as above but instead of evaluatina the time-varving coefficients of Eqs. (25) and (28) in theinitial steady state the coefficients are given their final steady state values (curves c), and the averages of the initial and final values (curves b). (6” in the last term of Eq. (25) is assumed zero) Curves d. Solutions used to check the accuracy of the present numerical methods. These are as solutions c, but &, in the last term of Eq. (25) is also evaluated in the final steady state. These curves approach asymptotically the final steady state.
Fig. 6. Composition response of the liquid stream leaving plate 15. 3-
2 5-
Fig. 4. Composition response of the liquid stream leaving plate 33, ~~(0). and flow response of the liquid stream entering plate 33, W(8).
2-
5-
lufl
05
I
aw-
I
0
Of
4
-c )5Fig. 7. Composition response of the hquid stream leaving plate 6.
Fig. 5. Composition response of the liquid stream leaving the feed plate 24. 1696
Process dynamics of binary distillation Table 2. The constants used in the Laplace inversion
y= qo
=
100” 10-0”
oo= lOO‘5 J= Hl = H2 = H3 = H4 = H, = HB = H, = Hs =
10 - 10.9758588954
164.049656959 - 1070.20258017
3981.52528520 - 9247.18974157
13773XtO56045 - 13003.6575318
7432.85149041 H,=-2317.82126813 Hl0 = 298.414943523
number of plates is as large as 50, a strong accumulation of errors occur. The asymptotic values of the solutions d were found to be equal to the final steady state values to three significant figures. This was found to be also the accuracy of the Laplate inversion when the resuhs obtained by use of ten-point and nine-point inversion formulae were compared. Approximative
methods of solution
As an illustration some curves of liquid flow response to a step change in the feed flow rate have been sketched in Figs. 2-4. Here a symbol wy is used in the meaning I
1
w*=rf=T. To obtain the response directly in the time domain without inverting Eq. (49) some statistical distribution curves can be used. In distillation the assumption of equal time constants for the flow response is in genera1 valid and the Poisson distribution is useful. For platen we obtain
Fig. 8. Composition response of the distillate.
Kf
Al1 the curves have been constructed using the dimensionless time 0 = (L,“/H.o)t, where the numerical vahtes of LB0 and H.O are those in the enriching section. (Thus t/min = 0-320). Each response is defined as a fractional approach of composition to the correct final steady state value and is denoted by w,:
w,=E
.lJ’
The constants used in the Laplace inversion are given in Table 2. The weights have been computed by Nordén er al. [1 ll.
By comparing the solutions a. *. c it can be seen that the estimating of the time-varying coefficients in Eqs. (25) and (28) towards their final vahtes reduces in most cases the deviations of the curves from the final steady state. In addition the curves c show that the neglecting of the last term in Eq. (25) can cause very large deviations on some plates though the error in each individual yn is very smalt. This is a consequente of our choice of numerical values. As the error in each yn is into the same direction and the total
wz(t) =
(n_lM)!o TnmMexp (-
r)dr,
where K is the average of ah K,, (i = M . . . N). In constructing the curves of Figs. 2-4 the dimensionless time 0 is used and according to Table 1ah the time constants are equal. The response curves for composition and liquid flow in the reboiler have been redrawn in Fig. 9 to show how the initial part of the final composition response w forms in the case of a step change in the feed flow rate. It can be shown that on smal1 values of time w is wel1 approximated by superposing wr and w, which are determined for the plate in question. Thus on any platen * M Eq. (5 1) reduces to - w = SW,W*, where the symbols of fractional responses time domain this equation becomes
1697
w(e) =
10” wz(r)dw(B-T),
are used. In the
V. J. POHJOLA
and H. V. NORDÉN
Ftg. 9. Initial part of the final composition response of the bottom product, w(0), obtained from the responses of Fig. 2. which can be solved graphically as is done to obtain the curve in Fig. 9.
The numerical calculations in this work have been carried out on an electronic computer IBM 1620.
Résmné-On étudie le comportement dynamique d’une distillation binaire dans le voisinage dun état d’équilibre. On donne les principes de la solution et les équations générales de l’état linéaire. Ceux-ci sont appliqués à une étude plus détaillée d’un processus théorique non lineaire pour lequel on donne la réponse des compositions à de faibles perturbations du taux d’ahmentation. Les effets de certaines approximations nécessaires pour donner les équations différentielles linéaires avec des coefficients constants sont étudiés qualitativement. Une des caractéristiques intéressantes de ce modèle est la division du processus en deux processus intermédiaires, le transfert de la masse et l’écoulement liquide, et la combinaison, par superposition, des réponses obtenues indépendamment pour chacun des sous-processus. Ce procédé suggérant l’utilisation de la transformation de Laplace, une simple méthode numérique d’inversion de Laplace est aussi donnée. Zusannnenfassung-Das dynamische Verhalten binärer Destillation in der Nähe eines Gleichgewichtszustandes wird untersucht. Die allgemeinen linearen Zustandsgleichungen und die Grundlagen einer Lösung werden gegeben. Diese werden auf die detailliertere Untersuchung eines theoretischen, nichtlinearen Vorganges angewendet, fiir welchen die Reaktion der Zusammensetzungen auf geringe Störungen in der Geschwindigkeit der Zufuhr gelöst wird. Die Wirkungen gewisser Annäherungen, die benötigt werden, urn die Differentialgleichungen linear mit konstanten Koeffizienten zu machen, wird qualitativ studiert. Ein Merkmal des gegenwärtigen Modelles, das von Interesse sein dürfte ist die Aufspaltung des Vorganges in zwei Teilvorgänge, die Stoffibertragung und die Fhi~sigkeitsströmung, und die Kombination der fiir jeden Teilvorgang erhaltenen Reaktionen durch Uberlagerung. Da das Verfahren auf die Verwendung der LaplaceschenTransformation schliessen Iasst, wird eine einfache, numerische Laplacesche Inversionsmethode ebenfalls angegeben.
1698
Pergamon Press.
Printed in Great Britain.
Process dynamics of binary distillation V. J. POHJOLA and H. V. NORDEN Technical University, Department of Chemistry, Otaniemi, Helsinki, Finland (Fint received 16 September 1968; in revisedform 30 April 1969) Abstract-The dynamic behaviour of binary distillation in the vicinity of an equilibrium state is investigated. The genera1 linear state equations and the principles of solution are given. These are applied to a more detailed study of a theoretical nonlinear process for which the response of compositions to smal1 perturbations in the feed flow rate is solved. The effects of certain approximations needed to render the differential equations linear with constant coefficients are studied qualitatively. A feature in the present model which may be of interest is the dividing of the process into two subprocesses, mass transfer and liquid flow, and the combining of the independently obtained responses for each subprocess by superposition. As the procedure suggests the use of the Laplace transform, a simple numerical Laplace inversion method is also given. INTRODUCTION
THE AUTOMATIC control and optimization of processes and the advanced design of plants and equipment are al1 at their best based on detailed knowledge of dynamic behaviour of the process in question. In distillation, particularly, many difficulties arise when it is desired to predict how the process wil1 act as a function of time. There are phenomena the exact physical mechanisms of which are not known or there are no ways to characterize them simply enough. Such features of practical distillation as the large number of plates and components, the large number of process and the nonlinear relationships variables between them make the analytical attempts hopeless to solve synthetic models, and render the solutions on digital computers time-consuming and expensive. Though computers are becoming faster and capable of more precise simulation, it is felt by several authors [ 1-31 that there wil1 be a growing need of relatively simple short-cut methods for the purposes mentioned above. Several reviews on the vast literature covering the control and process dynamics of distillation have been made [4-61. The secondary effects in distillation, by which are meant physical phenomena other than mass transfer which may have some
effect on the response of compositions, are discussed by Rosenbrock [7]. He discusses also certain distinctive problems in distillation from automatie control’s point of view [3]. Articles dealing primarily with the mathematics involved are reviewed by Mah et al. [SI. This paper is concerned with a theoretical study of composition response of a nonlinear binary distillation, where the hydraulic lag, the most common secondary effect, is present. Starting from genera1 dynamic state equations it is arrived finally at a linear readily solvable model which is expected to describe the response in the neighbourhood of the initial steady state with an accuracy sufficient for control purposes.
GENERAL
MODEL
A linear mathematica1 model of binary distillation wil1 be developed in genera1 terms. A matrix notation wil1 be used to expose the basic structure of the model and to make it easier to follow the-manipulation of equations. Particularly the response of compositions to perturbations in the feed flow rate wil1 be studied. State equations
The dynamics of binary distillation, a typical multivariable process, can be described by a
1687
V. J. POHJOLA
set of ordinary first-order generally nonlinear, y
= i(t)
and H. V. NORDÉN
differential equations,
= f [a(t),p(t>].
(1)
Here a is a state vector and it includes al1 the variables that are needed to specify the state of the process. p is a vector of input variables. Consider an uncontrolled distillation process the state of which in certain range of conditions can be described by its compositions and liquid flow rates. Let x = [x,] be a column vector of compositions and 1 = [IJ be a column vector of liquid flows, the elements of the vectors being deviations from certain initial steady state values of compositions and flow rates of the liquid streams leaving plates n and m, respectively. (n and m may get values from 0 to N + 1 where N is the total plate number of the distillation column; 0 refers to the condenser and N + 1 to the reboiler.) Let the elements of p be deviations from certain initial steady state values of compositions and flow rates of the streams entering the column. Suppose that at least in some vicinity of an initial steady state the functional dependence of Eq. (1) can be approximated linear. It can be written then more explicitly
[‘i]=.[;]+RP,
=A[;]+K
where 0 signifies a zero matrix. In the same way b=
;.
[l
(5)
Thus
or written in an expanded form X=Px+QI+sf
i=m+tf.
(7)
03)
(2) The solution of state equations
where A and B are matrices which are in genera1 dependent on time. To determine the response of the process to certain input p requires the solution of Eq. (2) with a set of initial conditions. Suppose next that there is only one input variable and that this is some perturbation in the feed flow rate, J? The state equation can be written then [:]
complexity only the time-dependent elements are thought to be included in 1. The major simplification that can be made at this stage is to assume that the changes in liquid flows do not depend on changes in compositions. This assumption is justified if there are no considerable changes in latent heats during the transient, and if the changes in compositions do not affect the rate of heat transfer in the condenser and reboiler. With this simplification A can be partitioned into submatrices as follows
(3)
where b is a column vector. To avoid any extra
The solution of Eqs. (7-8) in the case of timevarying coefficients is tedious. A step-by-step integration procedure is suggested by Mah et al.[8]. We concentrate here, however, on a situation where the assumption of constant coefficients is justified. Such a situation arises when the deviations in the state variables are so smal1 that the coefficients practically are constant, or when the nonlinearity of the original problem is so strong that Eq. (2) is for any practica1 purpose an equally poor approximation of Eq. (1) whether the dependence on time of the coefficients is taken into account or not. In the case of constant coefficients the solutions are
1688
Process dynamics of binary distillation
XQ) =
J;ew
I(t) = ILexp
P(~-T)I[Q~(T)
+s.(T)ld~
[R(t-T)]tf(T)
dr.
(9) (10)
The evaluation of the vector convolution integral of Eq. (9) wil1 be examined in the following. The vector Ql sfcan be written in the form
+
[lf1
where
Z is a constant
-f ,
matrix.
If in addition
is denoted for brevity by v, Eq. (9) becomes X(t) = 10 eXp [P(?-T)]
zV(T)&.
(12)
Consider temporarily a solution where the hydraulic lag is neglected and f is a unit step change. Then each element of 1 is also a unit step change. Let this particular solution beu: t
exp [P(t-T)]
u(t) =
1 Z 1 dT.
0
(13)
Ijl
111 Now the expression du(t-7) dt
x(t)
=
du(t--7)
vtT)dT
dt
Thus it has been shown that the original problem to determine the response of compositions to smal1 perturbations in the feed flow rate can be divided into two sub-problems: the determination of composition response to a unit step in the feed flow rate neglecting hydraulic lag, U, and the determination of liquid flow response v to the actual changes in J This procedure allows to treat the mass transfer and the liquid flow phenomena separately. Further the two responses can be superposed by (17) which is valid. Let f be any function of time.
il !
MODEL FOR PROCES
A PARTICULAR
In the following the composition response of a theoretical binary distillation to perturbations in the feed flow rate wil1 be studied in detail. The unsteady state material balances around the reboiler, condenser and each plate wil1 be linearized and solved along the main lines brought out in the preceding chapter. Dejinition of the process
1
= exp [P(t-i)]Z
DETAILED
(14)
1
The pr.incipal nomenclature of the process is given in Fig. 1. The column is assumed to operate
represents the response of compositions to a unit impulse in the feed flow rate when the hydraulic lag is neglected. It can be written also in the form dU(iyT)
= eXp
[P(t-T)]
z.
(15)
The relation between the elements of u and U is u,=
(17)
(11)
[l
Ql+sf=Z
Substitution of (15) into (12) gives the final composition response in terms of U and v:
E Unk, k=l
(16) Fig. 1. The nomenclature of the distillation process.
where K is the total number of elements in v. 1689
V. J. POHJOLA
and H. V. NORDÉN
with a partial reboiler and a total condenser. The following additional assumptions are made: (1) Constant molal overílow on the plates in the steady state. (2) Pressure changes in the column during the transient are negligible. (3) The vapor hold-up in the column and the liquid hold-up’s between the plates are negligible. (4) The liquid hold-up’s on the plates are perfectly mixed. (5) The feed and the reflux enter the column as a liquid at its boiling point; any additional condensation or vaporization due to potential temperature differente between the entering liquid and the liquid on the feed plate or the top plate is negligible. (6) The vapor leaving the reboiler is in equilibrium with the bottom product. Thus no enthalpy balances wil1 be considered.
Material
where the smal1 letters h, 1, p,f, x and y are used for the deviations from the initial steady state AH, AL, AP, AF, AX and AY, respectively. As no changes occur in the vapor flows, the deviations of the VS are zero. In addition the vapor flow is constant throughout the column, so the subscripts of the vapor flow terms wil1be dropped A total material balance is dH n = Lpl dt
%
(20)
+fn.
(21)
= In_1- (l,+p,)
Combining ( 19) and (2 1) we obtain (H,O+ h,) %
= (Li-1 + In-&-I - (L,o + p,o+ l,_, +fn)x,
A material balance of the more volatile component around a typical plate in the unsteady state is
d(Hdf”)= L._,X,_,
- (L,+
P,)X,-
v,r,
+J’n+~y,+~+Fttx~.
(18)
Here P, is a product stream from plate n and F, is a feed onto plate n. (Po = D, PN+l = B, and FM = F. P, = F, = 0 otherwise.) Multiplying both sides of (18) by a forward differente operator, A, an expression is obtained where the variables are deviations from the initial steady state and in addition the constants are defined in the initial steady state and can be assumed known. Eq. (18) becomes (X,O+x,)
%
= (LE-1 + ln_l)X,_l - (L,o+ 1, + P,o+p,)x, - v,“y?z + C+1Y,+, --KI”(h+P,)
P,)+ F,,
or
balances
(H,O+ h,) %+
- (L, +
+ x-ll,-1
+xF”fn7
(19)
+ I/oY,+, + cx;+
- vy,
- x,o)I,-1
+ (X,O - XnO)fn.
(22)
Equation (22) is a differential equation with time-varying coefficients. Its nonlinearity is determined by the nonlinearity of the dependence of y upon x and 1. A straightforward means to render the coefficients constant is to neglect the products of differences in (22). This is based on the assumption that the deviations from the initial steady state are small. With this simplification Eq. (22) becomes H,’ %
= L:_~x,_~ - (L,O + P,“)x, - voy,
+ WY,+, + (X;-I
-X,t”)f,-~
+
(XF” - X,O>f,. (23)
The dependence of y on x and possibly on 1 is determined by the phase equilibrium relationship and the Murphree plate efficiency. An explicit expression is obtained when the efficiency for the liquid phase and the equilibrium curve in the form Xn =
1690
Proces
E
xn-1-&_
=
dynamics of binary distillation
xn-1-x,
n x,-,-x:
(24)
- x,-1-qY(
Rearranging and multiplying by A we obtain
where again the products of differences have been neglected. Eq. (32) can be written more briefly Yn = %“&l-I
&(Y,z) = (l-B,z”-Sn)xa-l+
(E,O+&)xn
- (XL - XnO)Sn, 1 where Z = -
(25)
[=AEandx=AX.
E’
AS
Ady,) - AYn = cp’tyn”), where Y,,mis defined by Y”O=z Y,” S Y?l,
(27)
Acp(Y,) =
(28)
we may put
where y = A Y. If the equilibrium curve is available in the form YZ = $(X,), the following relationship may be used
where oo, po and yo are constants defined in the initial steady state. It is interesting to note that after neglecting the products of differences it is the behavior of the plate efficiency which determines the nonlinearity of (33) and consequently of (23). Unfortunately there are no data at hand on how the plate efficiencies do behave in transient conditions. Therefore no attempt is made here to replace ,$ by some linear function of x and/or 1. Rather we content ourselves to assume the plate efficiencies to remain constant. Thus (33) reduces to y, = (Y,OX,_l+&OX n.
H = F%+zv.
h?h l+X&,-E,
5oco ei51 Cl a2
(30)
6,
x0
x1
PX=
x2
c2
.
where
.
.
.
9(36) .
h
=E,
n
L,
=IY,*-Kt
n
L,X,-X,*’
$7ö)[( 1-
Yn =
n
-
aN
(31)
Here en is the vapor efficiency and m, is an approximation of the slope of the phase equilibrium curve on plate n. Combining (25) and (28) and solving for y, we obtain finally
(35)
The first term on the right-hand side is
There is also a simple relationship between the Murphree plate efficiencies for liquid and vapor phases: =
(34)
After combining Eqs. (23) and (34) the material balances of the more volatile component for n = 0 . . . N + 1 can be expressed as a set of linear differential equations with constant coefficients
(29)
E,
(33)
+ PRox, + Yn0L
. bN
CN
x,
aN+l
bNfl
xN+l I
where
I
a, =
-L (Lz-, - V%xno) Hn”
(37)
b,=-~[L.O+Pr+r(-a~+,+B.“)]
(38)
c, = J-
(39)
ZnO)x,-1+ ~,Ox,
Gc-l-&“)~nl~
(32) 1691
Hn”
vop:,, .
V. J. POHJOLA
and H. V. NORDÉN
It is to be noted that L-, = LN+l = 0, (Y,,= & = 0, CI,,,+~ = &+2 = 0, and as the plate efficiency in the reboiler is equal to unit y (assumption (6)), (llN+l= 0. The second term is 0 ____---------_
-_’
zv= ZM
>
Response of compositions to a step change in the feedflow rate neglecting hydraulic lag It is obvious from the expression of the final composition response, (17), that most of the advantages of the model presented herein are taken by the use of the Laplace transform. Giving the elements of v in Eq. (35) a constant valuef and taking the Laplace transform with respect to time, Eq. (35) becomes
(40)
1.M
Z h4+1
(SI-P)X Z N+l
= zs-’
(46)
f,
1N
where (41) Z,=$,(X;_l-X.o),
M < ut < N+l.
(42)
As the changes in the feed flow rate,f, affect the liquid flows below the feed plate only, the liquid flows in the enriching section remain unchanged. So the elements from 0 to M - 1 have been omitted in 1. In the total material balance, (2 l), the relation between h and 1 must be known. Suppose [9] that l,, = &h,,
(SI - P)ü = zs-1.
(47)
LiquidJEow response Equation (44) is Laplace-transformed respect to time to give
(43)
where K, is a constant. The total material balances for 12= M . . . N can be written now i = Kl+kf,
where the bar signifies the transformation and s is the transform variable. Equation (46) is a set of linear algebraic equations with a tridiagonal coefficient matrix and is easily solved by Gaussian reduction or by Crout’s method, a modification of the former. The solution contains s as a parameter. An expression analogous to ( 15) is
(SI - K)Ï = kf:
(48)
This can be solved analytically. we obtain
(44)
1, =
with
1
On plate n
(49)
f,
ifiM(&-‘s + 1)
where
(-KM)
KM
i 1M+ 1
KV+1(-K44+1) Kw+z(-&,+z ) KI + kf=
0 0
+* . I
KN(-KN>
1692
f.
(45)
Process dynamics of binary distillation
where f is the transform of an arbitrary bance in the feed flow rate.
distur-
Superposition
Equation transformed
(17) can be written notation simply x = SÜY.
in a Laplace-
(50)
a special emphasis is put on a universal compact formulation of the equations.
and
AcknowledgmentsThe authors wish to thank 0. A. Asbjernsen, Dr. Techn., for valuable criticism. Also, thanks are due to the Computing Center of Helsinki University for giving free machine time for calculations on an electronic computer IBM 1620.
For example on plate 12: NOTATION K x,
=
s
I: k=l
Unrl&.
Laplace inversion
The final composition response is obtained by the Laplace inversion of (50): x(t) = L-’ Z(S) . 1 1
coefficient matrices state vector coefficient column vector elements of P bottom product flow rate distillate flow rate Murphree liquid efficiency El F F.U feed flow rate onto plate n F, AF AF, f relationship in vector ; functional notation H5 weight in Eq. (53) liquid hold-up on plate n H, hz AH, 1 unit matrix J, K integers K average of time constants for flow response time constant for flow response, K, defined in Eq. (43) K coefficient matrix, defined in Eq. (45) k coefficient column vector, defined in Eq. (45) flow rate of liquid leaving plate 12 L 1, AL, 1 = [[“l column vector of liquid flows L Laplace transform operator M feed plate of the slope of the m, approximation equilibrium curve on plate n N total number of plates zero matrix 0 P,Q,R coefficient matrices A, B
(51)
a b a,, b,, c, B D
(52)
Nordén [ 101 has developed a simple numerical inversion formula
Gt(ti)
=
i
C
HjGt(sj-i)9
(53)
*3=1
where the arguments are chosen from geometric series ti = tod
(54)
and s._.=!&!d-i 3 I to
The H;s are weights independent of t. v is a constant, slightly larger than 1. qo and t. are suitably chosen constants. A detailed discussion of the method is given in [lol. DISCUSSION
The present paper reports some of the fundamental results of a larger research program of the authors and some co-workers directed towards the process dynamics and computer control of continuous distillation. In the text 1693
V. J. POHJOLA
and H. V. NORDÉN
flow rate of product stream from plate n PZ AP, P input vector s Laplace transform variable S coefficient column vector t time variable t coefficient column vector constant in Eq. (54) u = rz7.k: composition response matrix which neglects hydraulic lag u = [4ll U in column vector form V vapor flow rate V, V liquid flow response v = [Gtl augmented
-
p,
vector =
Superscripts 0
- -
fractional response of composition fractional response of composition Wl to a step change in the feed flow rate when the hydraulic lag is neglected, defined in the Appendix fractional response of liquid flow, w2 defined in the Appendix liquid composition on plate n X?l &I A-Kl [&l column vector of compositions Y?l vapor composition on plate n Yfl Ar,, Z coefficient matrix W
x =
f1 11
Greek symbols % Cl- %)/
m-l-x,“)/d(y,)
forward differente operator, deviation from the initial steady state value Murphree vapor efficiency constant in Eq. (55) dimensionless time variable = &2°1H,o) t constant defined in Eq. (3 1) constant in Eq. (54) LI& A? dummy time variable phase equilibrium relationships
f *
Subscripts B D F M n P
initial steady state value final steady state value phase equilibrium composition
bottom product distillate feed feed plate plate, numbering to reboiler product stream
condenser
Other notations A bar signifies a Laplace-transformed variable A prime is used for a derivative Matrices are denoted by bold-face capita1 letters whereas lower case bold-face letters are used for vectors
REFERENCES
[ll SARGENTR.
[2] 131 [4] [5] [61 [7] [8] [9] [lO] [l 11
from
W. H., Chem. Engng Prog. 1967 6371. LEFKOWITZ I., Chem. Engng Prog. Symp. Ser. 1963 59 178. ROSENBROCK H. H., Chem. Engng Prog. 1962 58 43. ARCHER D. H. and ROTHFUS R. R., Chem. Engng Prog. Symp. Ser. 196157 2. ROSEBROCK H. H., Trans. Instn chem. Engrs 1962 40 376. WILLIAMS T. J., Chem. Engng Prog. Symp. Ser. 1963 59 1. ROSENBROCK H. H., Trans. Insm chem. Engrs 1962 40 35. MAH R. S. H., MICHAELSON S. and SARGENT R. W. H., Chem. Engng Sci. 1962 17 619. SMITH B. D., Design ofEquilibrium Stage Processes. McGraw-Hill 1963. NORDÉN H. V.,ActaAcod.Abo. Mathematicaet Physica 196122 1. NORDÉN H. V. and SEPPÄ I.,Actupolytech. stand. Ch. No. 36, 1965.
1694
Process dynamics of binary distillation
APPENDIX It has been shown that the response of compositions to perturbations in the feed flow rate can be obtained by superposing the response of compositions to a unit step change in the feed flow rate when the hydraulic lag is neglected, and the response of liquid flows to the actual perturbations. To arrive at this some simplifications have been made which concern the balance equations of the more volatile component. Thus it might be of some interest to solve numerically this part of problem, that is Eq. (46), to get some idea about the effects of the approximations which led to Eq. (34). As the exact solutions of the original nonlinear equations are not known, we content ourselves to compare the response curves from (46) on large values of time with the final steady state values of changes of compositions which are computable. Giving the deviation terms appearing in the coefficients of (22) their final steady state values one ought to get solutions which approach asymptotically the final steady state values of the changes of compositions if (34) is valid. If, however, the plate efficiencies in the final steady state differ from those in the initial steady state, these asymptotic values deviate from the truc final values by amounts which are related to the errors caused by neglecting Sn in (33), (and to a less extent to the errors which are caused by neglecting the changing of the slope of the equilibrium curve, VA, dur+g the transient). The numerical values which wil1 be given next have been chosen somewhat arbitrarily. However, the initial and final values of the plate efficiencies have been chosen so that the error in yn caused by neglecting e. in (33) is about the same magnitude and into the same direction on every plate. This can be regarded as the most unfavorable situation as no cantelling of errors occurs. It is emphasized that no genera1 conclusions from this kind of study can be drawn and al1 what follows is meant purely as qualitative.
Table 1. The constants used in computing the initial and final steady state compositions F” = 100 kmol/hr f= 5 kmol/hr Do = 50 kmol/hr Loo/Do x,0 N M
= = = =
3.75 0~500 50
optimum in the initial steady state (= 24) H,O= 1 kmol, n= l...N Ho0 = 5 kmol HN,, = 5 kmol
h,/kmol = 0.004 I,,/kmol/hr,
n= M .. N
h ,V+1-0 Y* = 1.5x-0.5x
l=-X+0.8, 0~ XS 0.2 0.2 s X =S 1 1E = -0.25X+0.65,
Response
curves
The response curves for condenser (n = O), reboiler (n = 51), feed plate (n = 24), and plates 6, 15, 33, and 42 are given in Figs. 2-8. The notation is as follows: 1
Steady state
Taking into account the assumptions made earlier any steady state of the process can be defined by fixing six process degrees of freedom[9]. The six variables chosen in this study are: feed flow F, feed composition X,, distillate flow D, reflux ratio LdD, total number of plates N, and feed plate location M. Instead of giving M a numerical value it is required that the feed plate location be optimum in the initial conditions. In addition, the following quantities have to be known: the phase equilibrium curve, the plate efficiencies for each plate in the initial and final steady states, the liquid hold-ups in the reboiler, condenser and on each plate, their dependence on the liquid flow rate, and the magnitude of the step change in the feed flow rate. The numerical values are given in Table 1. The analytic expressions for the phase equilibrium curve and for the dependence of the Murphree vapor efficiency on composition are assumed to be valid both in the initial and in the final steady states. The steady state compositions are obtained by iterative plate-by-plate calculation. The optimal feed plate location is attained in the initial steady state with M = 24. This value of M maximizes the differente XDo- XBO.
1695
CESVd.?.4NaIl-D
log 0 Fig. 2. Composition response of the bottom product, w,(0), and flow response of the liquid stream entering the reboiler, w,(0).
Fig. 3. Composition response of the liquid stream leaving plate 42, w,(fI), and flow response of the liquid stream entering plate 42, w*(B).
V. J. POHJOLA
and H. V. NORDÉN
Curves a. Solutions of Eq. (46) with the elements of P evaluated giving the time-varying coefficients of Eq. (22) their final values instead of the initial ones appearing in Eq. (23). Curves b and c. Solutions which are obtained as above but instead of evaluatina the time-varving coefficients of Eqs. (25) and (28) in theinitial steady state the coefficients are given their final steady state values (curves c), and the averages of the initial and final values (curves b). (6” in the last term of Eq. (25) is assumed zero) Curves d. Solutions used to check the accuracy of the present numerical methods. These are as solutions c, but &, in the last term of Eq. (25) is also evaluated in the final steady state. These curves approach asymptotically the final steady state.
Fig. 6. Composition response of the liquid stream leaving plate 15. 3-
2 5-
Fig. 4. Composition response of the liquid stream leaving plate 33, ~~(0). and flow response of the liquid stream entering plate 33, W(8).
2-
5-
lufl
05
I
aw-
I
0
Of
4
-c )5Fig. 7. Composition response of the hquid stream leaving plate 6.
Fig. 5. Composition response of the liquid stream leaving the feed plate 24. 1696
Process dynamics of binary distillation Table 2. The constants used in the Laplace inversion
y= qo
=
100” 10-0”
oo= lOO‘5 J= Hl = H2 = H3 = H4 = H, = HB = H, = Hs =
10 - 10.9758588954
164.049656959 - 1070.20258017
3981.52528520 - 9247.18974157
13773XtO56045 - 13003.6575318
7432.85149041 H,=-2317.82126813 Hl0 = 298.414943523
number of plates is as large as 50, a strong accumulation of errors occur. The asymptotic values of the solutions d were found to be equal to the final steady state values to three significant figures. This was found to be also the accuracy of the Laplate inversion when the resuhs obtained by use of ten-point and nine-point inversion formulae were compared. Approximative
methods of solution
As an illustration some curves of liquid flow response to a step change in the feed flow rate have been sketched in Figs. 2-4. Here a symbol wy is used in the meaning I
1
w*=rf=T. To obtain the response directly in the time domain without inverting Eq. (49) some statistical distribution curves can be used. In distillation the assumption of equal time constants for the flow response is in genera1 valid and the Poisson distribution is useful. For platen we obtain
Fig. 8. Composition response of the distillate.
Kf
Al1 the curves have been constructed using the dimensionless time 0 = (L,“/H.o)t, where the numerical vahtes of LB0 and H.O are those in the enriching section. (Thus t/min = 0-320). Each response is defined as a fractional approach of composition to the correct final steady state value and is denoted by w,:
w,=E
.lJ’
The constants used in the Laplace inversion are given in Table 2. The weights have been computed by Nordén er al. [1 ll.
By comparing the solutions a. *. c it can be seen that the estimating of the time-varying coefficients in Eqs. (25) and (28) towards their final vahtes reduces in most cases the deviations of the curves from the final steady state. In addition the curves c show that the neglecting of the last term in Eq. (25) can cause very large deviations on some plates though the error in each individual yn is very smalt. This is a consequente of our choice of numerical values. As the error in each yn is into the same direction and the total
wz(t) =
(n_lM)!o TnmMexp (-
r)dr,
where K is the average of ah K,, (i = M . . . N). In constructing the curves of Figs. 2-4 the dimensionless time 0 is used and according to Table 1ah the time constants are equal. The response curves for composition and liquid flow in the reboiler have been redrawn in Fig. 9 to show how the initial part of the final composition response w forms in the case of a step change in the feed flow rate. It can be shown that on smal1 values of time w is wel1 approximated by superposing wr and w, which are determined for the plate in question. Thus on any platen * M Eq. (5 1) reduces to - w = SW,W*, where the symbols of fractional responses time domain this equation becomes
1697
w(e) =
10” wz(r)dw(B-T),
are used. In the
V. J. POHJOLA
and H. V. NORDÉN
Ftg. 9. Initial part of the final composition response of the bottom product, w(0), obtained from the responses of Fig. 2. which can be solved graphically as is done to obtain the curve in Fig. 9.
The numerical calculations in this work have been carried out on an electronic computer IBM 1620.
Résmné-On étudie le comportement dynamique d’une distillation binaire dans le voisinage dun état d’équilibre. On donne les principes de la solution et les équations générales de l’état linéaire. Ceux-ci sont appliqués à une étude plus détaillée d’un processus théorique non lineaire pour lequel on donne la réponse des compositions à de faibles perturbations du taux d’ahmentation. Les effets de certaines approximations nécessaires pour donner les équations différentielles linéaires avec des coefficients constants sont étudiés qualitativement. Une des caractéristiques intéressantes de ce modèle est la division du processus en deux processus intermédiaires, le transfert de la masse et l’écoulement liquide, et la combinaison, par superposition, des réponses obtenues indépendamment pour chacun des sous-processus. Ce procédé suggérant l’utilisation de la transformation de Laplace, une simple méthode numérique d’inversion de Laplace est aussi donnée. Zusannnenfassung-Das dynamische Verhalten binärer Destillation in der Nähe eines Gleichgewichtszustandes wird untersucht. Die allgemeinen linearen Zustandsgleichungen und die Grundlagen einer Lösung werden gegeben. Diese werden auf die detailliertere Untersuchung eines theoretischen, nichtlinearen Vorganges angewendet, fiir welchen die Reaktion der Zusammensetzungen auf geringe Störungen in der Geschwindigkeit der Zufuhr gelöst wird. Die Wirkungen gewisser Annäherungen, die benötigt werden, urn die Differentialgleichungen linear mit konstanten Koeffizienten zu machen, wird qualitativ studiert. Ein Merkmal des gegenwärtigen Modelles, das von Interesse sein dürfte ist die Aufspaltung des Vorganges in zwei Teilvorgänge, die Stoffibertragung und die Fhi~sigkeitsströmung, und die Kombination der fiir jeden Teilvorgang erhaltenen Reaktionen durch Uberlagerung. Da das Verfahren auf die Verwendung der LaplaceschenTransformation schliessen Iasst, wird eine einfache, numerische Laplacesche Inversionsmethode ebenfalls angegeben.
1698