498
Ind. Eng. Chem. Process Des. Dev. 1986, 25, 498-503
windup algorithm presented here. The simulation used in Figure 1 was rerun with a Proportional Integral controller using four different reset windup protection algorithms. The controller used a gain of 3.0 and a reset rate of 0.01 min-'. Figure 2 shows the results when using these different PID controllers. Line A shows the response when using the Velocity form algorithm with the reset windup protection described here. Line C shows the response when using the Position form algorithm presented here. These two algorithms show the best response of the four curves. A Position form algorithm where the integration is stopped as soon as the input to the process saturates is shown as line D. The response is seen to lag behind the first two curves since the control action is significantly lower when it comes off of the saturation limit a t time 71. This is because the integration is stopped as soon as the saturation limit is reached rather than when the Integral action begins to cause saturation. Line B shows the response of a Velocity form controller where the incremental contributions to the output are discarded if the controller is saturated. This algorithm again shows very poor response as the control action a t time 51 pulls the control input prematurely off of the limit. This is the same damaging problem that caused the poor control action in Figure 1. This last algorithm throws away the increments of the output that are calculated beyond saturation, hence, interfering with both the Proportional action and the Integral action. This algorithm still shows an offset of 4.1% at 200 min, whereas the normal position
form algorithm shows an offset of 1.1%,the improved position form algorithm presented here shows an offset of 0.3%,and the improved velocity algorithm presented here shows an offset of 0.2%. Conclusions The antireset windup protection algorithms presented in this paper allow saturation of PID and Minimum Variance controllers to be handled correctly. The algorithm described for the Velocity form of these controllers is significantly different from that generally used and allows performance equivalent to Position form algorithms even for a purely Proportional controller. Similar conclusions apply to the Smith Predictor and Dahlin Controllers since they are special cases of the Minimum Variance controller. Literature Cited Astrom, K. J.; Wittenmark, B. "Computer Controlled Systems Theory and Design"; Kailath, T., Ed.; Prentice-Hall: Englewood Cliffs, NJ, 1984. Dahlin, E. B. Instrum. Control Syst. 1988, 4 , 77. Harris, T. J.; MacGregor, J. F.; Wright, J. D. Can. J . Cbem. fng. 1982, 6 0 , 425. Khandheria, J.; Luyben, W. L. Ind. Eng. Cbem. Process Des. Dev. 1978, 15, 278. MacGregor, J. F. Department of Chemical Engineering, McMaster University, Personal Communication, 1981. Segall, N. L. M.Eng. Thesis, McMaster University, Hamilton, Canada, 1983. Smith, C. L. "Digital Computer Process Control"; Intext: Toronto, Canada, 1972. Smith, 0. J. M. Cbem. Eng. Prog. 1957, 5 3 , 217.
Received f o r review December 17, 1984 Accepted September 5, 1985
Design of MultHoop SISO Controllers in Illhmivariable Processes Cheng-Chlng Yu and Wllllam L. Luyben" Process Modeling and Control Center, Department of Chemical Engineering, Lehlgh University, Bethlehem, Pennsylvania 180 15
A practical design procedure is proposed for determining the
structure, variable pairing, and tuning of multiloop
SISO controllers in a multivariable-processenvironment. The selection of controlled variables relies primarily on engineering judgment. The choice of manipulated variables is based on Morari's Resiliency Index (minimum singular value of the process transfer function matrix). Variable pairing is decided by first eliminating those pairings that give negative RGA's or negative Niederfinksi or Morari Indexes of integral controllability. The final variable-pairing selection is based on the Tyreus load rejection criterion. The simple BLT controller tuning method, recently proposed by Luyben, is used to give consistent, logical comparisons of alternatives. The method is illustrated by applying it to a multivariable ternary distillation column example.
A recent paper (Luyben, 1985) proposed a practical, easy-to-use method (BLT tuning) for tuning SISO (single-input-single-output) controllers in a multivariable environment. The technique is based on detuning all controllers equally from Ziegler-Nichols settings until the desired degree of closed-loop stability of the system is achieved as indicated by a multivariable Nyquist plot. The structure of the system must be specified to use this tuning method: what variables are controlled, what variables are manipulated, and how they are paired in the SISO structure. This paper addresses these structural questions. A logical, systematic design procedure is proposed which produces a stable, workable multiloop SISO system. 0196-4305/06/1125-0498$01.50f0
It should be emphasized that multivariable controllers are not considered in this paper. There are no claims that the multiloop SISO structure gives the best control system. It is possible that a multivariable controller could improve control performance. Multivariable controllers (such as Internal Model Control [IMC], Dynamic Matrix Control [DMC], Linear Quadratic Control [LQ], etc.) will be compared to multiloop SISO systems in a future paper. What is claimed about the proposed procedure is that it will produce a workable, stable, simple SISO system with only a modest amount of engineering effort. This conventional SISO system can then serve as a realistic benchmark, against which more complex multivariable controller structures can be compared. This procedure has 0 1986 American
Chemical Society
Ind. Eng. Chem. Process Des. Dev., Vol. 25, No. 2, 1986
been tested on open-loop stable processes only. The application to open-loop unstable processes will be the subject of future research. Outline of Procedure Before the proposed design procedure is discussed in detail, a brief outline of the steps is given below. 1. Select Controlled Variables: Use primarily engineering judgment based on process understanding. , 2. Select Manipulated Variables: Find the set of manipulated variables that gives the largest Morari Resiliency Index (MRI). 3. Eliminate Unworkable Variable Pairings: Eliminate pairings with negative RGA’s or negative measures of integral controllability (Niederlinski and Morari Indexes). 4. Find the Best Pairing from the Remaining Sets: (a) Tune all combinations using BLT tuning and check for stability (characteristic loci) and robustness (Doyle-Stein Index). (b) Select the pairing that gives the lowest magnitude closed-loop load transfer functions: Tyreus LoadRejection Criterion (TLC). Discussion of Procedure The starting point is plant transfer functions relating all controlled variables to all manipulated variables. It is assumed that these W transfer functions are available ( N is the number of SISO loops in the multivariable process). The N transfer functions relating the controlled variables to the load disturbance must also be known. This aspect of the problem is by no means a trivial job, particularly for the many chemical engineering systems that exhibit strong nonlinearities (reactors, high-purity distillation columns, etc.). Plant test data can sometimes be used to develop these transfer functions, but noise and changing conditions make these transfer functions difficult to obtain reliably. Pulse testing of a nonlinear mathematical model of the process is reasonably successful in most systems, but analytical or numerical linearization methods must be used in some nonlinear processes. A. Choice of Controlled Variables. Usually the controlled variables are fairly easy to select. Engineering judgment, based on a good understanding of the process, leads in most cases to a logical choice of what needs to be controlled. Considerations of economics, safety, constraints, availability and reliability of sensors, etc., must be factored into this decision. We offer no new guidance on this question in this paper. It should be remembered that these controlled variables can be directly measured (temperature, pressure, composition, etc.) or computed from other directly measured variables (ratios, heat removal rates, inferred product compositions, mass flow rates, extensive variables, etc.). Occasionally some guidance can be obtained from singular-value decomposition (Downs, 1984). B. Choice of Manipulated Variables, Once the controlled variables have been specified, the “control structure” depends only on the choice of manipulated variables. For a given process, selecting different manipulated variables will produce different “control structure” alternatives. These “control structures” are independent of the “controller structure”, i.e., multiloop SISO controllers or one multivariable controller. For example, in a distillation column, the controlled variables might be chosen as distillate and bottoms compositions ( X , and XB).Now we could choose the manipulated variables reflux and vapor boilup (R and V). This choice will give one possible control structure. Had we chosen to manipulate distillate flow and vapor boilup (D and V), we would have a different control
499
structure for the same basic distillation process. Morari has studied methods for assessing the inherent resilience of alternative control structures (Morari, ,1983). Based on this work, we propose to use the “Morari Resiliency Index” (MRI) to guide the selection of manipulative variables. MRI = a[G(iw)] (1) The MRI is the minimum singular value (a)of the plant transfer function matrix G(io). The set of manipulated variables that gives the largest minimum singular value over the frequency range of interest is the best. This selection of control structure (manipulated variables) is independent of variable pairing (controller structure) and controller tuning. The MRI is a measure of the inherent ability of the process (control structure) to handle disturbances, model plant mismatches, changes in operating conditions, etc. The larger the value of MRI, the more resilient the control structure. The problem of the effect of scaling on singular values is handled by expressing the gains of all the plant transfer functions in dimensionless form. The gains with engineering units are divided by transmitter spans and multiplied by valve gains. This yields the dimensionless gain that the control system has to deal with. C. Eliminate Unworkable Variable Pairings. 1. Eliminate Pairings with Negative RGA’s. Pairing on negative RGA’s is undesirable because it produces either an unstable system or a system that lacks integrity (will go unstable with some loops on manual: Grosdidier et al., 1985). 2. Eliminate Pairings with Negative Niederlinski Indexes. The Niederlinski Index (NI) is defined NI =
det [G(O)I
(2)
bikf0) k=1
where G(o) = steady-state gain matrix of plant transfer functions and gii(o)= diagonal elements of G(o). If all the SISO controllers con”tainintegral action and have positive loop gains, a negative value of the Niederlinski Index shows that the system will be closed-loop unstable with this variable pairing for any controller tuning (Niederlinski, 1971). 3. Eliminate Pairings with Negative Morari Indexes of Integral Controllability (MIC) MIC = A[.:,]
(3)
where X[G] = eigenvalues of the matrix G and G& = plant steady-state gain matrix with the signs adjusted so that all diagonal elements have positive signs. If all the SISO controllers contain integral action and have positive loop gains, a negative value of any of the eigenvalues of Gf,) shows that the variable pairing will produce an unstable closed-loop system. Note that the application of MIC is limited to multiloop SISO controllers and strictly proper process transfer function matrices. MIC is a special case of Theorem 7 presented by Grosdidier et al. (1985). D. Final Variable-Pairing Selection. All the remaining pairing possibilities must be examined. For each case, controllers are tuned by using the BLT procedure. Then the closed-loop load-rejection performances are compared in order to find the pairing which does the best job of keeping controlled variables constant in the face of load disturbances. 1. BLT Tuning (Luyben, 1985). This simple procedure uses a multivariable Nyquist plot and equal detuning
500
Ind. Eng. Chem. Process
Des. Dev., Vol.
25, No. 2, 1986
of all loops from the single-loop Ziegler-Nichols settings. (a) Calculate the Ziegler-Nichols settings for each P I controller by using the diagonal elements of G. (b) Assume a detuning factor “F”,and calculate controller settings for all loops. Kci = K z N ~ / F TIi
= FTZN,
(4) (5)
n h M 4 weir
(c) Define the function W(iu,
W(iu)= -1
+ det [I + G(iu)B(iuJ
weir
(6)
where B = diagonal matrix of SISO feedback controllers. (d) Calculate the closed-loop function L,(,)
n C,
.52
t
(e) Adjust the detuning factor F until the peak in the L, log modulus curve (Lcm)is equal to 2N (8)
where L,, is the biggest log modulus in the L, curve and N is the number of SISO loops in the multivariable system. Luyben (1985) demonstrated that this method gives reasonable, conservative estimates for the controller settings on 10 case studies. It provides a method of obtaining settings for all possible pairings on a consistent basis so that realistic comparisons can be made. 2. Tyreus Load-Rejection Criterion (TLC). Tyreus (1984) proposed the use of the load-rejection capability of the closed-loop system as a rational method for selecting the best variable pairing. The closed-loop relationship between the load variable L and the controlled variables X is
X = [I + GB]-l&L
(9)
where X = vector of N controlled variables, G L= vector of open-loop plant transfer functions relating the controlled variables to the load disturbance, aQd L = a single-load disturbance. The attenuation of the load disturbance by the process and the control system can be shown in a plot of the magnitude of each element of Xas a function of frequency. A perfect control system would keep all the XL’sa t zero. The best variable pairing is the one that gives the smallest magnitudes of the XL’s.This is the TLC criterion. We have formulated the procedure above for a singleload disturbance. Obviously the method can be applied to systems with several load disturbances provided appropriate weighting factors can be defined for the effects of each load input. E. Robustness of Closed-Loop System. The ability of the closed-loop system to remain stable despite changes in process parameters is called robustness. Doyle and Stein (1981) proposed two useful measures of robustness. DSO = a [ I + (GB)-’]
(10)
DSI = (T[I + (BG)-l]
(11)
The final variable pairing and controller tuning should be checked for robustness by plotting DSO and DSI as functions of frequency. Singular values below 0.3-0.2 indicate a lack of stability robustness. Further controller detuning may be required.
Distillation Example The ternary distillation column, studied by Yu and Luyben (1984) for single-loop control, is presented as an
w 127 pria
(7)
L,, = 2N
6 35 cm length I 6 1 m
height
B
432
c3
1Ox1O6
iC4
,163
nC4
,837
Figure 1. Steady-state design of ternary deisobutanizer column. R -V- I
D-V- I
RR-V
Figure 2. Alternative control schemes.
example. Several other multivariable distillation systems have been studied by using this procedure and will be reported in later papers dealing with these process applications: heat-integrated columns (Chiang, 1985) and complex distillation configurations (Alatiqi, 1985). The steady-state design of the system is given in Figure 1. The ternary propane/isobutane/normal butane system is studied at 120 psia. Propane is the lighter-than-light key component, and its concentration in the feed is 6 mol %. The column has 32 trays, a total condenser, a partial reboiler, saturated liquid reflux, and theoretical trays. A. Controlled Variables. Since we wish to control the distillate and bottoms product compositions, X D N and XBN will be chosen as controlled variables where X D N = normalized mole fraction of iC4 in the distillate = mol of iC, in distillate/(mol of iC4 + mol of nC, in distillate), XBN = normalized mole fraction of nC4 in bottoms = mol of nC4 in bottoms/(mol of iC4 + mol of nC4 in bottoms). These normalized compositions are used so as to avoid the problems that occur when trying to control the absolute composition in the presence of non-key components. Composition analyzers are used with two different sampling times ( T , = 1 and 4 min). B. Manipulated Variables. Many possible choices of manipulative variables exist. The conventional alternatives are reflux flow and vapor boilup ( R and V), distillate flow and vapor boilup (Dand V), and reflux ratio and vapor boilup (RR and V). Figure 2 shows some possible control schemes. Table I gives the open-loop process transfer functions G and G Lfor this process with several choices of manipulated variables. The load disturbance
Ind. Eng. Chem. Process Des. Dev., Vol. 25, No. 2, 1986 501 Table I. Ouen-Loop Transfer Function Matrix and Load Transfer Functions g11 g1z gz1
gzz
R-V" 16.3/(18.18 1) -18.0/(21.49 1) -26.2~'.~~/[(47.6S + 1)(1.2S + l ) ] 28.63/[(47.69 + 1)(4.4S + l ) ]
+ +
gL1
gLZd
RR- V' 1.5/[(34.28 + 1)(4.1S l ) ] -1.12/(1.4S 1) -2.06e-2.3S/[(23.6S 1)2] 4.73/[(19.79 + 1)(0.7S + 111 0.60e-1.9S/[(17.7S+ 1)(0.9S + l ) ] -0.53e-5.2S/ [ (16.4s + 1)2]
D- Vb -2.3/[(29.49 1)(0.2S l)] 1.04e-19S/(1.35S 1) 3.1e-1.5s/[(19S + 112] 1.8(77S + 1)/[(21.2s + 1)(0.9S + l ) ] 0.55e-2.3S/[(15.85S+ 1)(0.9S + l ) ] -0.48e-5.7S/[(17.7S + 1)2]
+
+
+
+
+
+
DTransmitterspan of 0.2 mol fraction for XDNand XBN, valve gain of 1457 g-mol/s for V, and valve gain of 1420 g-mol/s for R. "Valve gain of 160 g-mol/s for D. "Valve gain of 17.75 for RR. With units of mole fraction/mole fraction. R-V D-V ............RR-V
Table 11. Steady-State Analysis of Alternative Control Structures control variable MRI Q O ) RGA NI structure pairing 0.114 -86.1 -1/86.1 R-V R-V-1 (XDN-R,XBN-V) R-V-2 0.114 87.1 1/87.1
~
-.. ..
MIC hG&$
45.1, -0.16 44.1, 0.24
1.86
0.56
1/0.56
2.04 f 1.89i
1.86
0.44
1/0.44
2.06 iz 1.75i
0.89
1.48
1/1.48
5.33, 0.9
is considered to be the propane concentration in the feed
ZF (LLK)
X = GM + &ZF(LLK)
(12) Table I1 gives the MRI values, the RGA values, the Niederlinski Indexes (NI), and Morari Indexes of Integral Controllability (MIC). The D-V structure has a bigger MRI value, and therefore this choice of manipulated variables is recommended (Table I1 and Figure 3). The R-V structure has a very small MRI (0.114 at w = 0). Therefore, this structure is inherently sensitive. For example, a 1% change in steady-state gains could make G(o)singular. Also notice that the R-V scheme with the conventional pairing XDN-R and XBN-V (R-V-1) shows a negative RGA. This unusual result is due to the multicomponent nature of the system. As can be seen in the gains of the open-loop transfer functions (Table I), V has more of an effect on XDNthan does R. At this point, the manipulated variables would be chosen as D and V. However, for purposes of illustration, the R-V and RR-V structures will be studied further in order to demonstrate that the MRI is an effective criterion. C. Eliminate Unworkable Variables Pairings. (1)Values of Values of RGA indicate both pairings (DV-1 and D-V-2) for the D-V structure give stable systems (Table 11). As mentioned earlier, the conventional R-V-1 pairing gives a negative RGA. Therefore, only the R-V-2
..", I
,001
.a
I
I
I
*I
I
IO
I
100
w
Figure 3. Morari Resiliency Indexes MRI for R-V, D-V, and RR-V control structures.
pairing can be used. The RR-V structure also has one stable pairing as indicated by RGA. (2) For a 2 X 2 system, NI contains essentially the same information as RGA (Table 11). However, this is not true for higher-order systems. (3) Values of MIC also show that both D-V-1 and D-V-2 are integral-controllable, and there is only one possible pairing for R-V and RR-V structures. For a 2 X 2 system, RGA, NI, and MIC give essentially the same results. But for higher-order systems, there are cases where MIC gives more information than RGA and NI. Before getting into the final variable pairing selection, the R-V structure will be investigated briefly. Table I11
Table 111. Dynamic Analyses of Alternative Control Structures control scheme
TS, min
R- V- 1
0
R-V-2
0
tuning method Buckleyd (SISO) Buckley
F factor
K,,/K,,
T ~ I / T ~ ~
DSO
DSI
stabilityb (linear)
stability' (nonlinear)
Oe
no
Yes
gain" margin
1
1.84/0.39
4.2/30.0
2
-1.95/-0.46
4.3/19.8
0.50
0.70
5.5
Yes
no
3.35 2.02 3.88 3.21 2.34 3.98 1.73 2.48
-2.45/0.54 -1.17/0.31 0.25/0.60 0.30/0.35 5.32/1.01 3.12/0.59 2.0/0.42 1.39/0.29
13.1/18.4 26.9/28.7 138.2/100.2 130.2/124.4 24.3/11.9 41.4/20.3 39.2/24.6 56.0/35.2
0.63 0.67 0.49 0.42 0.45 0.87 0.68 0.89
0.26 0.36 0.30 0.38 0.10 0.22 0.20 0.30
4.1 2.2 1.7 1.6 3.0 5.4 2.4 3.5
yes Yes yes yes yes yes yes Yes
yes Yes no no no yes no Yes
(SISO) D- V -1 D- V-2 RR- V
1
4 1 4 1 1
4 4
BLT-4 BLT-4 BLT-4 BLT-4 BLT-4 BLT-2 BLT-4 BLT- 1
'Smallest gain margin of characteristic loci. *Stability on linear transfer function model. Buckley's single-loop tuning (Luyben, 1973). e System with positive feedback.
Stability on nonlinear rigorous column model.
502
Ind. Eng. Chem. Process Des. Dev., Vol. 25, No. 2, 1986 R - V - I (T,*O)
-4J-
7k
-
R V- 2 ( T8*O)
I
Figure 4. Characteristic loci plots for R-V-1 and R-V-2 pairings.
Figure 6. Characteristic loci plots for D-V-1 and D-V-2 schemes with BLT-4 tuning. 01
R-V-2
b
__
R-V-I
20
40
I
SO
80
'\
too
TIME ,mln
Figure 5. Nonlinear step responses for L-V-1 and L-V-2 schemes for a step change in ZF(LLK) (ZF(LLK) = 0.10).
gives the controller parameters for both R-V-1 and R-V-2 pairings with no analyzer dead time ( T , = 0). Characteristic loci plots (Figure 4) show that the R-V-1 scheme is closed-loopunstable as indicated by MIC. Characteristic loci plots also indicate that R-V-2 is stable and has a gain margin of 5.5. However, the nonlinear simulations, Figure 5 , show that the R-V-1 is stable and the R-V-2 is unstable. This structure is very sensitive to nonlinearities and was eliminated from further study. D. Final Variable-PairingSelection. (1) There are three possible control schemes left: D-V-1, D-V-2, and RR-V. All the controllers were tuned with BLT tuning for L , = 4 (BLT-4). Table I11 summarizes the controller parameters. Reasonable stability margins were obtained as shown in characteristic loci plots (Figure 6) and gain margins in Table 111. (2) The Tyreus Load-Rejection Criterion (TLC) was used to determine the final variable pairing. Figure 7
t
m
.'a
.'I
w
I
1
I
IO
1
KK)
Figure 7. Load-rejection plot (TLC) for D-V-1, D-V-2, and RR-V schemes.
shows that the D-V-1 scheme attenuates the load disturbance (ZF(LLK)) much better than the D-V-2 scheme. Therefore, D-V-2 is eliminated according to TLC. There is no clear preference between the D-V-1 and RR-V schemes (Figure 7) in terms of TLC. E. Robustness of Closed-Loop System. Doyle and Stein Indexes (DSO and DSI) give good measures of robustness. Figure 8 shows the lack of robustness for the RR-V scheme (DSI = 0.2 at w = 0.2). To get the same degree of robustness, the RR-V scheme had to be detuned to a +1 Lcm. Therefore, D-V-1 is the control scheme recommended. Figure 9 shows the responses of the column for step changes in feed composition. The D-V-1 scheme gives good response. The RR-V scheme is unstable for the +4 L,, tuning. It should be emphasized that we can determine the best control structure at the very beginning (from MRI).
Ind. Eng. Chem. Process Des. Dev., Vol. 25, No. 2, 1986 503
This procedure gives a workable, stable, simple control system with a modest amount of engineering manpower. It has been successfully tested on a number of multivariable distillation examples.
D -V- I (Ts* 4) RR-V (TS-4) - - - -
-8
1'
54 ,001
.h
I
I
I
.I
I
lo
W
Figure 8.
Doyle-Stein Indexes ( D S O a n d D S I ) RR-V schemes.
i
loo for D-V-1 a n d
,
Greek Letters X = eigenvalue CT = minimum singular value B = maximum singular value rIL= reset time of ith controller, min
-.-I-
D-V-I(BLT-4)
Nomenclature B = bottoms flow rate B = controller transfer function matrix D = distillate flow rate F = feed flow rate F = detuning factor in BLT G = open-loop process transfer function matrix g = ijth element of G dL= load transfer function vector I = identity matrix K,, = controller gain of ith controller KZN,= Ziegler-Nichols gain of ith controller L = load variable L, = closed-loop log modulus M = vector of manipulated variables L,, = maximum closed-loop log modulus N = order of system (number of SISO controllers) R = reflux flow rate RR = reflux ratio = R J D T , = analyzer dead time V = vapor boilup rate W = -1 + det (I + GB) X = vector of controlled variables XB = bottoms composition XBN= normalized bottoms composition (mole fraction nC,) XD = distillate composition XDN= normalized distillate composition (mole fraction iC4) ZF(LLK) = feed composition of lighter-than-light key
R R - V (BLT-4) ----
= Ziegler-Nichols reset time of ith controller, min = frequency, rad/min
7ZN,
w
Literature Cited
t
0
I
90
40 40 TIME, min
do
do
Figure 9. N o n l i n e a r step responses of D-V-1 a n d RR-V for a step change in ZF(LLK) (ZF(LLK) = 0.10).
schemes
However, for purposes of illustration, several alternatives were included throughout this example. Conclusion A practical design procedure is presented for determining the control structure, variable pairing, and tuning of multiloop SISO controllers in a multivariable process.
Alatiqi, I. Ph.D. Thesis, Lehigh University, 1985. Chiang, T. P. PhD. Thesis, Lehigh University, 1985. Doyle, J. C.: Stein, G. IEEE Trans. 1981, AC-26, 4. Downs, J., paper presented at the Lehigh University Distillation Control Short Course, Bethlehem, PA, May 1984. Grosdidier, P.: Holt, 8. R.: Morari. M. Ind. Eng. Chem. Fundam. 1985, 24, 221. Luyben, W. L. "Process Modeling, Simulation, and Control for Chemical Engineers": McGraw-Hill: New York, 1973; p 420. Luyben, W. L. Ind. Eng. Chem. Process Des. Dev. 1986, 25,326. Morari, M. Chem. Eng. Sci. 1983, 38, 1881. Niederiinski, A. Aufomatica 1971, 7, 691. Tyreus, B. D., paper presented at the Lehigh University Distillation Control Short Course, Bethlehem, PA, May 1984. Yu, C. C.: Luyben, W. L. Ind. Eng. Chem. Process Des. Dev. 1984, 23, 590.
Received for review April 8, 1985 Revised manuscript received August 30, 1985 Accepted September 18, 1985
