Academia.eduAcademia.edu

Soft-sensors for process estimation and inferential control

1991, Journal of Process Control

This paper presents two adaptive estimators (software based sensors or 'soft-sensors') for inferring process outputs that are subject to large measurement delays, from other (secondary) outputs which may be sampled more rapidly. In other words, these estimators utilize plant data sampled at different rates. The parameters of the estimators can be continuously estimated and updated on-line, thus enabling the tracking of slow variations in process characteristics. The estimators employ either an input-output or a state space process description as the starting points for algorithm synthesis. In contrast to mechanistic model-based estimators such as Kalman filters, the proposed adaptive techniques require minimal design effort. The key contribution of the paper is thus the formulation of applications independent, adaptive multi-rate algorithms, which can provide accurate estimates of infrequently measured process outputs, from other more rapidly sampled secondary outputs. Theoretical developments are supported by results of recent applications to a variety of industrial scale processes: estimation of biomass concentration in an industrial mycelial fermentation; top product composition of a large industrial distillation tower and melt flow index on an industrial polymerization reactor. Measurements from established instruments such as off-gas carbon dioxide in the fermenter; overheads temperature in the distillation column and hydrogen concentration in the reactor were used as the secondary variables for the respective processes. The range of the applications is an indication of the utility of the techniques. A significant improvement in overall process control performance is also possible when estimated plant outputs, rather than the infrequently obtained measurements, are used for feedback control. This is demonstrated by non-linear simulation studies.

FU N Pa pe rs Soft-sensors for process estimation and inferential control Ming T. Tham, Gary A. Montague, A. Julian Morris, and Paul A. Lant Department of Chemical and Process Engineering, University of Newcastle, Newcastle upon Tyne, NE1 7RU UK (Received 13 June 1990; revised 20 August 1990) This paper presents two adaptive estimators (software based sensors or ‘soft-sensors’) for inferring process outputs that are subject to large measurement delays, from other (secondary) outputs which may be sampled more rapidly. In other words, these estimators utilize plant data sampled at different rates. The parameters of the estimators can be continuously estimated and updated on-line, thus enabling the tracking of slow variations in process characteristics. The estimators employ either an input-output or a state space process description as the starting points for algorithm synthesis. In contrast to mechanistic model-based estimators such as Kalman filters, the proposed adaptive techniques require minimal design effort. The key contribution of the paper is thus the formulation of applications independent, adaptive multi-rate algorithms, which can provide accurate estimates of infrequently measured process outputs, from other more rapidly sampled secondary outputs. Theoretical developments are supported by results of recent applications to a variety of industrial scale processes: estimation of biomass concentration in an industrial mycelial fermentation; top product composition of a large industrial distillation tower and melt flow index on an industrial polymerization reactor. Measurements from established instruments such as off-gas carbon dioxide in the fermenter; overheads temperature in the distillation column and hydrogen concentration in the reactor were used as the secondary variables for the respective processes. The range of the applications is an indication of the utility of the techniques. A significant improvement in overall process control performance is also possible when estimated plant outputs, rather than the infrequently obtained measurements, are used for feedback control. This is demonstrated by non-linear simulation studies. zyxwvutsrqponmlkjihgfedcbaZYXWVU (Keywords: estimation; feedback control; sensors) In many process control situations, due to sampling limitations, the infrequent measurement of some process outputs prevents the early detection of load disturbances. This can result in large deviations from setpoints and long disturbance recovery times. Often, these adverse effects cannot be acceptably overcome even by the use of existing advanced control algorithms and can lead to unsatisfactory system performance. The problem of controlling infrequently measured process outputs has long been studied and publications in this area date back to the early 1970s. There are essentially two approaches to the problem. One is to formulate special algorithms to control the infrequently sampled outputs. For example, Soderstrom’ formulated a number of minimum variance controllers enabling the manipulated input to be changed between the sampling intervals of the primary process output. However, these control algorithms were only developed for first-order plant models and no comparative results were presented. Parrish and Brosilow2 also proposed a controller design method based upon the reconstruction of disturbance effects. The controller parameters are determined on-line by heuristic tuning rules. Simulation results were presented to demonstrate the superior performance of their control strategy over conventional PID control. A second approach is to use the information provided by other easily measurable variables to provide an estimate of the controlled output. The estimated outputs can then be used for overall plant control. Control schemes based on the feedback of output estimates are often termed ‘inferential control’ schemes. An ideal situation arises when the plant states are completely observable from the secondary outputs. Kalman filtering techniques can then be employed to estimate plant states using secondary output measurements, and hence estimates of the controlled output computed from its relationship with the states. Control of the plant is achieved by feedback of either the state estimates or the output estimates to appropriate controllers. The literature on the above methods is extensive, and previous publication+ are examples of the more recent contributions describing the application issues of these techniques. However, the application of Kalman filtering techniques is limited to those processes that are completely observable from the secondary outputs. For most plants, such a set of J. Proc. Cont. 1991, Vol zyxwvutsrqponmlkji 1,January 3 O95%1524/9l/OlOOO~l2 h ,041 Rllt,rrurnrth_"PinPmlnn 1 ,A Soft-sensors for process estimation and inferential control: M. T Tham secondary outputs can be difficult to determine and, in some cases, may not even exist. Brosilow et ~1.~3~ have suggested an estimator design technique based upon an input-output representation of the plant. The design is approached by obtaining a least squares static estimator, which can be used to infer the controlled output from secondary measurements at the steady state. The estimator is then additionally made applicable to transient periods by incorporating heuristically derived lead-lag elements into its structure. Methods for minimizing the steady state estimation error by appropriate choice of the secondary measurements were also proposed. However, a set of secondary measurements for which this error is satisfactorily small may not exist. In fact, the evaluation work carried out by Patke et al.8 indicated that an inferential control scheme based on the output estimation technique proposed by Brosilow et al. can result in significant offsets. Moreover, to apply the technique, the gains and approximate time constants of the controlled output and all secondary outputs, to all plant disturbances and manipulated inputs, must be known. The control strategies reviewed above rely either on infrequent measurements of the controlled output or on the use of secondary measurements. It is, however, possible to use both types of measurements. One such technique is to set up a ‘parallel cascade’ control strategy9, in which measurements of the controlled output are fed to a controller whose output acts as the setpoint to a secondary output controller. Although no offset problems arise with this control configuration, transient performances can at times be pooI.8. Another approach is to use both measurements of the controlled output and secondary output within an adaptive inferential control framework. In such a scheme, the infrequent measurements of the controlled output are used only for parameter estimation while output estimation, and hence plant control, is achieved using a secondary output at its faster sampling rate. To our knowledge, very little investigative work has been carried out in this area, although a similar algorithm has been proposed previously’o. The algorithm’O is, however, restricted to the assumption of a first-order plant model, an equal number of disturbance inputs and secondary measurements, and the use of a dead beat control law. There are a number of industrial situations where infrequent sampling of the ‘controlled’ process output can present potential operability problems. Common examples are in product composition control of distillation columns@ and chemical reactors”. In these cases, the sampling delay is a direct result of the long cycle time of on-line composition analysers. Due to the potential problems associated with infrequent sampling, the control of product quality of many industrial multicomponent or high purity columns, for example, is commonly achieved by regulating an a priori chosen tray temperature at its setpoint. However, single-temperature feedback control is not always effective as maintaining a constant tray temperature does not necessarily result in constant product compositions. A similar problem arises in polymerization reactors, in which measurements of 4 J. Proc. Cont. 1991, et al. reactor feed rate and one or more reactants can be used to infer polymer properties that would otherwise only be available by relatively slow on-line or off-line analysis. In fermentation processes, important internal variables such as biomass, substrate and secondary prod&t concentrations characterize the state and progress of the fermentation. However, they present measurement difficulties. In spite of recent encouraging developments in ion selective and enzyme sensors, optical and high frequency based methods, most of the concentration variables in a fermentation liquid phase cannot be measured accurately and reliably on-line. Therefore, laboratory analyses/assays are usually required to support fermentation supervision and control’2. Nevertheless a number of successful strategies have been reported’3-23. It is interesting to observe that mechanistic model-based approaches and extended Kalman filtering are mainly employed, even though sufficiently accurate biochemical models are required. This paper presents two algorithms that utilize plant data sampled at different rates to estimate the controlled output. For example, overheads temperature is used to estimate top product composition of a distillation column; and biomass concentration of a fermentation is estimated using off-gas carbon dioxide (CO,) evolution rate. The first estimation algorithm is based upon a general input-output representationz4, while the second, summarized for comparison purposes, is derived from a state space representation of the plant*s. The original ideas behind the two methods are also briefly assessed in Guilandoust and Morri@. The extensions to Kalman filtering suggested by Kreisselmeier*‘, and the self-adaptive Kalman filtering approach suggested by Young28, use a similar model form to that proposed in this work. This major difference, however, lies in the multi-rate aspects addressed here. zyxwvutsrqponmlkjihgfedcbaZYXWVUTS Adaptive estimation using an input-output process representation Process description The process dynamics are assumed to be represented by the following input-output (I/O) relationships: r&) = W- ‘)u(t - m’) + I&- ‘)W) (1) v,(t) = W- ‘)u(t - (2) Y(O =.Yo(t - v(t) = vo(0 + 4 + m2) + I,G- ‘MO v,(t> v*(t) (3) (4) where: t is a time index; vO(t) is the controlled output (e.g. distillation column product composition, fermenter biomass concentration); and v,(t) is the secondary output (e.g. column temperature, fermenter off-gas CO,).Their corresponding observations, y(t) and v(t), are contaminated by measurement noises, v,(t) and v,(t) respectively, which are assumed to be zero mean white sequences. The Vol zyxwvutsrqponmlkjihgfedcbaZYXWVUTSRQPONMLKJIHGFEDCBA 1,January Soft-sensors for process estimation and inferential control: M. T. Tham et al. and output measurement delay is denoted by zyxwvutsrqponmlkjihgfedcbaZYXWVUTSRQPONMLKJIHGFEDCBA ‘d’, and is expressed as an integer multiple of the secondary output sampling interval. Note that the discretization time step T(t) = I&- ‘)w(t) + v*(t) (9) is based upon the smaller sampling interval of the As q% ,(t) and v,(t) are zero mean white sequences, it secondary variable. Thus, values of y(t) only become follows that c(t) is also a zero mean and white sequence. available every d time steps whereas v(t) is measureable Additionally, since w(t) is a row vector of stationary at each time step. The time delays in the responses of the signals, according to the spectral factorization theorem30, controlled and secondary outputs, to changes in the Equations (8) and (9) may be expressed as: manipulated input u(t) (e.g. column reflux flow, bioreactor feed rate) are ‘m,’ and ‘m2’, respectively. The term zyxwvutsrqponmlkjihgfedcbaZYXWVUTSRQPONMLKJIHGFED S(t) = H,(q-‘)W (10) war’ is a vector of unmeasureable disturbances. G,(q-I) and G,(q-I) and all elements of the row vectors and I,(q-I) and I,(q-I) are discrete time transfer functions with the order of numerator less than or equal to that of i(t) = KG- W(t) (11) the denominator. The term q-l is the backward shift operator. It is assumed that changes in any one of the where H,(q- 1) and &(q- 1) are stable and proper polynocomponents of w(t) affect both yO(t) and v,(t), otherwise mial ratios and the sequence w(t) is zero mean and white. v,(t) will not be suitable for estimating y,(t). Therefore an It is noted that both k(t) and c(t) are generated by the element of I,(q-I) cannot be zero, unless the correspondsame disturbance sequence o(t). If this was not so, then ing element in I,(q- I) is also zero, and vice versa. the two signals would clearly be uncorrelated. This cannot be the case, because variations in load disturbances Process disturbances and their representation are assumed to affect both controlled and secondary In deriving the estimators, it is important to consider outputs-see Equations (Z+-(9). Similarly, neither how disturbances are to be modelled. It has been H,(q-I) nor H,(q-1) can be zero without the other being observed4 that process disturbances are characterized by zero also. infrequently occurring random changes, lasting for perSubstitution of Equations (10) and (11) into Equations iods of the order of at least the process response time. (5) and (6), and elimination of o(t), results in: This behaviour has been modelled as Brownian motion processes3~4~7J9, and was chosen because it can be pictured y(t + d) = [G,(q-‘)q-“lfm2 - G&‘)H,(q-‘)/ as a sequence of step changes with independent random H2(q-‘)lu(t--2)+[H,(q-1)lH2(q-‘)lv(t)+~(t+d) 02) amplitudes, occurring at times described by a Poisson distribution. Even if this disturbance model may well for m, 2 m2. For the case in which m, 5 m2: represent the physical situation, its use poses an observability problem. The Brownian motion process descripy (t + d)= [G,(q-‘) - q-“2+“lG2(q-‘)H,(q-‘)/ tion contains a pole on the unit circle and hence is nonH2(q-‘)lu(t-m,)+[H,(q-‘)lH2(q-‘)lv(t)+E(t+d) (13) stationary. For the process to be observable from any set of measurements, the number of measurements should Both Equations (12) and (13) relate measurements of the therefore exceed the number of disturbance inputs4-a primary output, y(t), to the secondary output, v(t), and situation rarely met in process control. The random stepthe values of the manipulated input, u(l). If G,(q-I), like behaviour of process disturbances can alternatively G,(q-I), H,(q-I) and H,(q-I) were all known, u(t) and v(t) be modelled by the stationary exponentially correlated could be used in Equation (12) or (13) to compute: noise (ECN) process, which does not introduce the F(t+d)=y (t+d)- E(t+d) observability problem associated with the Brownian (14) motion process. Therefore, in the sequel, all process diswhich is an estimate of the controlled output, y ,,(t). From turbances will be assumed to be stationary random Equation (7), it is observed that the deviations of P(t + d) sequences so that each may be represented as the resfrom the controlled output is determined essentially by ponse of a stable filter to a white input sequencejo. the variations of v,(t) and v,(t). Nevertheless, noise effects can be reduced by judicious filtering of the Derivation of the I/O estimator respective variables. Thus, given a sufficiently accurate Equations (1) to (4) can be rewritten as: plant model, Equation (12) or (13) may be used to provide delay free estimates of the controlled output at y (t+d)=G,(q- ‘)u(t- m,)+~(t)+E(t+d) (5) the sampling rate of v(t). In situations where there is a large number of disturbance inputs or there are complex v(t) = GAq- ‘)u(t - m2) + i(t) (6) or even unknown plant dynamics, the polynomial ratios G,(q-I), G,(q-I), H,(q-I) and H,(q-I) can only be where approximately obtained using off-line identification. In this work, however, these transfer functions are not E(t + d) = q% ,(t) - v*(t) (7) assumed to be known apriori. Instead, it is proposed that their parameters be identified on-line. S(t) = I,(q- W (t) + v*(t) (8) J. Proc. Cont. 1991, Vol 1, January 5 Soft-sensors for process estimation and inferential control: M. T. Tham et al. To facilitate on-line parameter estimation, and based on the definitions of G’(q-‘), Gz(q-‘), H,(q-‘) and H,(q- ‘), Equations (12) or (13) can be rewritten as: using a suitable algorithm. The value of y(t) = OT4(t - d) is then calculated and inserted in +(I). The following relationship is then used to provide estimates of the primary output at those time instants when its measurements are not available: j(r + d) = OT$(f) (20) = -c@(t)-..-a&Q-(n1)d) where zyxwvutsrqponmlkjihgfedcbaZYXWVUTSRQPONMLKJIHGFEDCBA m = Min(m,,m3 and, + j$u(t-m1)+ ..+ l&&~-m-@ A(q-‘)= 1 +a,q-,+azq-*+...+a,q-” +r,v(t)+ .. +r,dv(t-nd) The parameters of Equation (19) are again updated at time t + d when a new value of the primary output again (16) becomes available. The expression for estimating y(t) at C(q- ‘)=c,+c,q- ‘t- c*q- 2+...+c,q- ” the faster sampling rate of the secondary output is thus given by Equation (20). A convergence proof for this algorithm, in the form of a multi-rate predictor for the The value of the integer ‘n’ is chosen to achieve the ‘slow’ process output given ‘fast’ manipulated input varirequired accuracy of representation. It should be noted able data, has been published3,. that Equation (15) is not of the normal ARMAX type, Equations (19) has n(2d+ 1) + 1 unknown parameters although its form may, at first sight, suggest otherwise. and ‘d’ may typically lie in the range 2 to 6. A first order This is because v(t) is not an unknown signal as in the estimator (n= 1) may therefore have between 6 and 14 case of an ARMAX relationship. Indeed, both the primunknown parameters and a second order estimator, ary output, y(t), and secondary output, v(t), respond to between 11 and 27. These figures may, at first sight, give changes in control input and disturbances. In its present the impression that the tuning-in period for larger values form, the parameters of Equation (15), cannot be estiof ‘d’ will be inconveniently long. This is, however, not mated in a straightforward manner beause y(t) is only the case. Although the number of unknown parameters available at every ‘d’ time steps. This problem can, increases with ‘d’, each time the controlled output is nevertheless, be resolved by making use of Equation (14) sampled, the data being supplied to the estimation algorand rewriting Equation (15) as: ithm is also increased [see Equation (19)]. Simulation studies, investigating the effects of analyser delays up to 6 A(q- ‘)j(t+d)=B(q- ‘)u(t- m)+C(q- ‘)v(t) (17) sample intervals, have in fact confirmed that the initial tuning-in period increases only marginally for increasing For a given value of ‘n’, it can be shown that the followvalues of ‘d.’ ing relationship: It has also been observed that estimates of the controlled output can often be obtained with sufficient (l+a~q- ~+...+u”(‘q- “d)j(t+d)= accuracy using a first order estimator. Additionally, the (p’q-’ + &q- *+ . . . + f&q- nd)u(t- m) parameter ‘ad of the first-order estimator may be shown +(r,+r,q- ‘+...+&dq- “d)v(f) (18) to be related to a, by a,, = (-a#. As may be deduced from Equation (12), Iall is normally less than unity, and Equation (17) are equivalent. However, Equation except for the special case in which H,(q-‘) has a domi(18) has the desired features because all of its parameters nant unstable zero. For larger values of ‘d’, it therefore are associated with measured data values, and may therefollows that cl,,-C < 1, especially when la, 1is significantly fore be estimated. Rearrangement of Equation (18) less than unity. As a result, ‘ud’ becomes negligible, and results in: Equation (20) simplifies to: y(t)= -Qj@-d)-...-u,&t-nd)+ fI,u(t- m- d- I)+...+p”&(t- m- (n+l)d)+ ~,v(t- d)+...+~ndv(t- (n+l)d)+~(t) (19) which can then be expressed more concisely as: y(t) = WQ(t - d) + e(t) ,..., u(t- m- d- (21) zyxwvuts Adaptive estimation based on a state space process representation x(t+ l)=Ax(t)+Bu(t-m)+Lw(t) (22) y(t) = Dx(t - d) + v,(t) (23) v(t) = Hx(t) + v,(t) (24) 1),..., v(t-d) ,... ] and e(t) is the equation error. When a measurement of y(t) is available, the parameter vector 0 is estimated 6 l)+...+~,$(f-m-6) Instead of adopting an input-output design approach, an estimator can also be derived based on the following state space model: where: $(I- d)=[- j$t- d) j(t+d)=fi,u(t-m+r,v(t)+...+rdV(t-6) J. Proc. Cont. 1991, Vol 1, January Soft-sensors for process estimation where r is again a time index, and y(t) and v(t) are sampled controlled primary and secondary outputs respectively. The corresponding measurement noises, v,(t) and vi(f), are assumed to be independent random sequences with zero mean and finite variances. The state vector is x(r)eRn. It has an initial value x(O), which has the following properties: E@(O)) = x, (25) and E{[x(O)- xJ[x(O) - xJT} = P, w(r)eR” is a vector of random and unmeasurable load disturbances. The primary output measurement delay is ‘d‘, while ‘m’ is the smallest of the time delays in the responses of the y(t) and v(r) to changes in the manipulated input, u(t). Any difference between the two delays can be included in the model by extending the state vector. Here, the aim is to estimate the system states using observations of v(r) and thence obtain an estimate of y(t) from Equation (23). If the matrices A, B, L, H and D, as well as the noise covariances were known, then the state estimates, %(r), could be computed from v(t) by Kalman filtering and used to provide estimates of y(t) using the expression: j(r + 6) = D%(t) (26) The Kalman filter for Equations (22) and (24) is give@ R(t+ l)=Ari(t)+Bu(t-m)+K(t)c(t) (27) where n(0) = iz, E(2)= v(r) - qr> and inferential ;(t)=HjZ(t) X(*+1)= [I’_, ; ‘I #)+ [j_)+ [lj&;) v(t)=[l,O,...,O]P(t)+c(t) From the definition, becomes: (33) e(r) = v(r)-D(r), Equation i(r) = i.,(r) (33) (34) Elimination of the state estimates i,(r) to Z,(t) from the right-hand side of Equation (32) by successive substitution yields: 1) .?-,(r+ l)= - azti(r)- .. - a,i(t- n+ 2)+ b&t-m) +..+b,u(r-m-n+2)+k&r)+..+k,e(t-n+2) .?“(r+ l)= -a,?(t) C(t) is the filtered value of v(t); K(r) is the Kalman gain determined via solution of the Riccati equations and c(t) is the innovations sequence. Although this method is one of the more popular approaches adopted for output estimation, it requires a rather detailed a priori knowledge about the process. This requirement can, however, be relaxed if the estimation is performed within an adaptive framework. In this case, an innovations model is used in place of Equations (22) to (24): M. T Tham et al. more practical to use the former set of expressions. Since the state estimates can be computed directly, the filtering process need not be performed. The matrices of the model can be parameterized in various ways. Since no prior knowledge about the model structure is assumed, it is reasonable that the parameterization be performed such that the identified parameters would directly correspond to those of an input-output mode133. This is best achieved by adopting an observer canonical model structure. Having estimated the states, Equation (30) could then be used to obtain j(r + d) if the row vector D was known. In principle, D could be estimated using values of y(r) and Equation (26), if the state estimates were not correlated. Unfortunately, this is not the case as may be noted from Equation (29). The problem can, however, be resolved by noting that in the observer canonical structure, Equations (29) and (31) are of the form: &(t+ l)=t(t+ (28) control: (35) + b,u(r - m) + k,c(t) Equation (35) is then substituted into Equation (30) leading to an expression with the following form: y(r)=P,u(r-m-d- l)+..+~,_,u(r-m-n-d+ 1) +r,i(t-d)+..+t,_,it(r-n-d+l)+Q(t-d-1) +..+&,_,~(r-n-d+l)+v,(r) (36) This can be expressed more compactly as: S(t+ l)=A%(t)+Bu(t-m)+Ke(t) y(r)=oTQ(r-d)+e(r) (29) zyxwvutsrqponmlkjihgfedcbaZYXWVUTSRQPONMLKJIHGFEDCBA (37) y(t) = D%(t - 6) + v,(t) (30) v(t) = H%(r) i- c(t) (31) where: K is the time invariant (limiting) Kalman gain, and % is the optimal estimate of the states, x. Although Equations (29)-(31) and (22)--(24) are equivalents2, it is with QT= [b,,..,b,_ I,r O,..,r”-l, 6 I,..>6n-1I +(r-d)T=[~(r-m-d-l) ,.., ^v(r-d) ,.., ~(r-d-I),..] and e(r) is the equation error. Each time a primary measurement becomes available, the parameter vector 0 can be estimated. Since Equation (37) is independent of state J. Proc. Cont. 1991, Vol 1, January 7 Soft-sensors for process estimation and inferential control: M. T Tham estimates, this should not present any difficulties. After updating 0 at time t, it is used in: =l3,u(t-ml)+..+p,_,u(t-m-n+ +..+r,_,i(t-n+l)+&c(t-1)+..+6,_&-n+1)(38) l)-tz,i@) to provide estimates of the controlled output. The values of c(t) in Equations (36) and (38) can be computed from Equations (32) and (33) once the parameters of Equation (32) have been determined34. 0 is again updated at time t + d when a new measurement of the controlled output becomes available. Equation (38) thus describes the relationship which provides estimates of the controlled output at the faster secondary output sample rate. The number of parameters to be estimated in Equation (36) is 3n - 2. If secondary output measurements can be obtained with low noise levels, and load disturbances only occur infrequently, then the ‘8’ terms tend to be very small and may thus be ignored. The calculation of jj(t + 6) becomes correspondingly simpler as identification of a secondary output model is no longer necessary. The result is a significant reduction in the number of unknown parameters and Equation (38) simplifies to: jj(t+d)=P,u(t--mt,v(t)+..+T,_,v(t-n+ l)+..+p,_,u(t-m-n+ 1) l)+ (39) At this point, it is interesting to note the similarity between Equations (39) and (21). This agreement between the final form of two essentially different approaches provides additional confidence in the arguments presented. The above discussion has implicitly assumed that the process is completely observable with respect to v(t).The question arises, however, as to what happens if this is not the case. From Equations (22~(24), it is evident that r(t) shares some of its poles with v(t). On the other hand, the zeros of v(t) and v(t) are, in general, different from each other. Loosely speaking, substituting the state estimates into Equation (23) may be regarded as assigning the poles of v(t) to r(t). The zeros assigned to y(t) are then determined by the parameter estimation procedure. As may be expected, this operation places the zeros such that the effect of the difference between the poles of v(t) and y(t) is minimized. Thus, it is also possible that some of the poles of v(t) are effectively cancelled. Simulation studies carried out on systems with widely different primary and secondary output dynamics indicate satisfactory estimation despite the differences in dynamics. Adaptive inferential control With the availability of ‘fast’ primary output estimates, (e.g. fermenter biomass concentration, distillation column product composition), closed loop adaptive ‘inferential’ control of that variable becomes a feasible control option. The values of P(t + d) can be used as a feedback 8 J. Proc. Cont. 1991, Vol 1, January et al. signal for any constant parameter or adaptive control algorithm. Examination of Equations (19) and (36) indicates that unless a proportional controller is used, closed loop identifiability problems due to correlated data should not arise [e.g. Ref. 351. Thus, for a time invariant process, the estimated parameters should converge to final values after an acceptable number of sample intervals. At this stage, parameter estimation could be halted if desired, with accurate estimates of the controlled output still being available. In applications to a slowly time varying process, it may be preferable to continue with the parameter estimation at all times. In this case, it may also be desirable to use an adaptive controller since the use of the adaptive inferential estimators does not preclude the application of adaptive control algorithms. Given a suitable choice of secondary variable, the inferential control strategy will automatically provide for ‘anticipatory’ control action (in a feedforward sense). This is because disturbance effects will manifest first in measured secondary variable response, and thus the strategy should yield significant improvements in disturbance rejection performance. Industrial evaluation Commercial scale plant validation of both the statespace and input-output based adaptive estimators (soft sensors) was undertaken in association with our industrial collaborators. For reasons of commercial confidentiality, in all the figures to be presented, the ordinate scales have been removed. In the following applications, a UD-factorized recursive least squares algorithm (UDRLS)33 was used to estimate the parameters of the estimators. In all the results to be presented, the estimators were initialized with zero parameters at the start of each application. Application to a continuous fermenter In this study, the process is a continuous stirred tank fermenter producing a fungal mycelium. Following downstream processing, the mycelia are the final process product. However, specifications placed upon its quality dictate certain fermenter operating conditions. In particular, biomass concentration has to be tightly regulated. The current strategy is based upon 4 hourly laboratory analysis of biomass concentration and the dilution rate adjusted upon the return of assay information. Although process analysis reveals that a sampling interval of 1 h would be more suitable for maintaining a constant biomass concentration, it is however, not possible to sample biomass at this frequency. In fact, present policy is based upon increasing the sampling interval, so as to reduce the demand for 24 h laboratory support and operator supervision. A move to 8 or 12 h analysis would therefore be desirable. This has important consequences on the control strategy in that an alternative technique to offline analysis must be found to regulate biomass effectively. The adaptive estimators provide a means by which Soft- sensors for process estimation and inferential control: M . T. Tham Adaptive Dilution inferential et al. estimation rate 6 c 2 2 s Ei Carbon dioxide in off-gas (%I zyxwvutsrqponmlkjihgfedcbaZYXWVUTSRQPONMLKJIHGFEDCBA 2 2 .2 m 0 400 200 0 100 Time Figure ments 1 Fermenter CO, evolution Adaptive 100 Figure 2 Fermenter pared with estimated 400 (h) (h) rate and dilution inferential rate measure- estimation 200 Time 300 200 Time (h) dry weight assayed approximately biomass concentration 4 hourly com- frequent estimates of biomass can be obtained using CO2 evolution rate and dilution rate. Measurements of the latter variables can be obtained frequently and with little delays. Hence, an inferential estimator supplied with hourly CO* and dilution rate measurements, supported by infrequent and irregular biomass concentration assays, will enable hourly predictions of biomass concentration. In this application, the state-space based estimator was implemented. Figure I shows the COZ and dilution rate measurements over a 470 h period of continuous operation. A major plant disturbance occurred at around time 285 h, causing a rapid fall in CO2 evolution rate. Although this was unexpected, it provided the estimator with a test of its ability to track major process changes. Note that at approximately 360 h, CO:! measurements became unavailable and was assumed constant for a period of about 20 h. The laboratory biomass assays (step-like response) and inferred biomass concentration over the fermentation period are compared in Figure 2. It demon- Figure 3 Fermenter pared with estimated dry weight assayed approximately biomass concentration 8 hourly com- strates the ability of the ‘state-space based’ estimator to provide very acceptable predictions of biomass transient behaviour. The biomass assay results can be seen to occur at intervals that are not completely regular. With the state-space based estimator, this is not important since the algorithm can accommodate such variations in assay times (P.A. Lant, PhD Thesis in preparation, University of Newcastle-upon-Tyne, UK). The performance of the estimator during the plant upset is also encouraging. It predicted the down-turn well, with only a slight undershoot at around 280 h. This is due to continuing estimator parameter adaptation to the new operating conditions. At time 248 h, the laboratory assay indicated a large increase in biomass concentration. It is interesting to observe, however, that the estimator continued to predict a fall based upon measurements of dilution rate and CO* evolution rate. Other possibly ‘suspect’ biomass assays can also be observed at times 58, 100 and 183 h respectively. Nevertheless, the ability of the estimator to predict the general trend of biomass concentration response was not compromised. The effect of increasing the interval of sampling for biomass concentration analysis from 4 to 8 h is shown in Figure 3. It can be observed that there was only a marginal deterioration in the quality of biomass estimates. However, and as expected, the initial estimator ‘tuning-in’ period is slightly longer (see Figure 2). Application to a high purity distillation column Distillation columns can also exhibit control problems, due to the delays introduced by on-line product analysers. The problem can be exacerbated in the case of columns with high purity products, where composition dynamics can be extremely non-linear. In one such column, an industrial demethanizer, the aim of the control strategy is to regulate the top product composition by varying reflux flow. Here, all process variables are measured every 5 min except for top product composition. Being measured by an analyser, composition values are available only at 20 min intervals. The aim, therefore, is J. Proc. Cont. 1991, Vol 1, January 9 Soft-sensors for process estimation and inferential control: M. T. Tham et al. Feed I I I 10 8 # flowrate I 12 Figure 4 1 I 14 Time Overheads I I 16 (min) 18 I 10 I 1 12 I I 14 Time x lo3 Figure 6 Demethanizer column feed flow rate measurement temperature (min) I I Time Figure 5 16 14 (min) I I 18 x lo3 and estimate I T 12 16 Demethanizer column temperature measurement Measurement I 10 I 1 18 Time (min) x lo3 x lo3 Demethanizer column reflux flow rate measurement Figure 7 Demethanizer column top product composition compared with estimated composition over 8250 min analysis degraded performances occurred only during large proto use ‘fast’ measurements of column overheads vapour cess transients. Recall that high purity columns possess temperature, together with reflux flow rate, to provide strongly non-linear composition dynamics, and thus estimates of product concentration at 5 min intervals. As represent a severe test of the resilience of the estimator. a secondary variable, a tray liquid temperature would be The lapse in performance quality is attributed to the preferable. The problems with the use of a vapour temestimator having to ‘tune’ to relatively frequent changes perature is that it is overly sensitive to vapour disturin plant operating conditions. Nevertheless, when new bances; changes due to column boil-up control and steadier operating conditions are approached, the estinumerous other factors which need not necessarily affect mation error is significantly reduced. This is in marked product composition. However, the choice of a vapour contrast to the steady state errors reported in the evalutemperature was dictated by its availability at the time of ation study of other inferential control schemes by Patke the tests. et a1.8Figure 8 shows the performance of the estimator Results are presented for a 137 h segment of a 20 day between 15000 and 16000 min, and provides a more application using the input-output form of the estimadetailed picture of estimator responses during a period tor. The data for feed flow rate, reflux flow rate and when process operating conditions are changing. It can column temperature are shown in zyxwvutsrqponmlkjihgfedcbaZYXWVUTSRQPONMLKJIHGFEDCBA Figures 4-6 respectibe observed that the estimates of top product compovely. A comparative plot of actual product composition sition are very close to the composition analyser results. with its estimates is shown in Figure 7. In spite of the Indeed, the errors are within the tolerances demanded by ‘bad’ choice of secondary variable, the estimator perthe plant operating personnel. formed remarkably well. It can be observed that 10 J. Proc. Cont. 1991, Vol 1, January Soft-sensors for process estimation and inferential control: M. T Tham et al. Measurement and estimates (15000 min < t < 16000 min) I I I I 2 0 f I 15 I I 15.2 Time Figure 8 compared I I 15.4 Demethanizer with estimated (min) I I 15.6 I I Figure 9 ment 2I I 41 Polymerization Application I 6, reactor ’ 8’ ’ ’ 10 fresh propylene to a polymerization Polymerization Figure 10 ’ ’ 12 I I (min) II I 8 I 10 12 I ’ zyxwvutsr 14 x lo3 reactor coolant flow rate measurement analysis 0 0 1 6 16 x lo3 column top product composition composition over 1000 min I1 4 Time I 15.8 I ’ I I 2 I I 4 I 6 Time ’ 14 feed rate measure- I Figure ment Polymerization 11 ’ (min) reactor 8I I I 10 I 1 12 ’ ’ 14 x lo3 hydrogen concentration measure- reactor As a further test of the utility of the estimator, the input-output algorithm was applied to an industrial polymerization process. The melt flow index (MFI) of the product was estimated using measurements of reactor feed rate, reactor coolant flow rate and hydrogen concentration above the reaction mass (shown in Figures 9-11 respectively). The secondary variables were measured at 10 min intervals, while measures of MFI were obtained from laboratory analyses every 2 h. Figure 12 shows the results of the very first trial application to this process and compares the laboratory assayed melt flow index (step-like response) and its estimates. The tuningin period, from zero initial parameter values took approximately 300 h (equivalent to 150 laboratory samples). Although the algorithm was implemented with minimal process knowledge used to ‘condition’ the estimator, its performance is good. In particular, the ability to track the change in polymer grade at around 1085 h is particularly encouraging. I 0 I I1 I,,,,,,,, 2 4 ,I 6 Time 8 (min) 10 12 14 x lo3 reactor melt flow index analysis Figure 12 Polymerization with estimated melt flow index over 14500 min J. Proc. Cont. 1991, compared Vol 1, January zyxwvutsrqpo 11 Soft-sensors for process estimation and inferential control: M. T Tham et al. Cooling water 1 I_ Implementation aspects In the application of the adaptive estimators, several important implementation issues arose. The data used for parameter estimation have to be suitably ‘conditioned’. Data scaling to ensure that input-output data ranges were of the same order, is particularly important. Care must also be taken to ensure that the measured signals are filtered for noise rejection, anti-aliasing, sudden changes or spikes and ‘rogue’ data. Failure to achieve satisfactory data conditioning will lead to poor estimator performance. In the above case studies, all process measurements were filtered using low phase-shift algorithms and transformed to represent deviations about a mean operating level. Data spikes were attenutated via the use of logic filters designed based on on-line calculated statistics of each data variable. Additionally, to promote robust estimator performance, whenever a primary output analysis (e.g. fermenter biomass, column product composition) becomes available, it is used in the algorithm to replace the corresponding estimated value. It is in ensuring the integrity of the algorithms to specific process problems that proved to be most demanding. For instance, periodic re-calibration of the infrared carbon dioxide analyser in the bioreactor application (or the on-line chromatograph in the case of the distillation column), unless accounted for, can cause significant problems. Due to the adaptive nature of the algorithm, however,it is capable of dealing with the slow drifts in calibration usually experienced. The techniques adopted here reflect the experiences with our work on CC adaptive and self-tuning contro118J9J6. zyxwvutsrqponmlkjihgfedcbaZYXWVUTSRQPONMLKJIHGFEDCBA l-0 + Bottom product Adaptive inferential control of a distillation column Schematic diagram of pilot plant distillation column: CR, analyser recorder; FRC, flow recorder/controller; GC, gas chromato- Figure 13 The adaptive estimators described above can also be graph; LC, level controller included as part of a feedback control system to provide an adaptive inferential control strategy. The practicality variable and is sampled every 0.5 min. The UDRLS of such a scheme to effect distillation product compoalgorithm was again used to estimate the parameters of a sition control has been evaluated by non-linear simula1st order input-output estimator. Although the perfortion. The process under consideration is a 10 stage pilot mances of estimators of various orders have been investiplant column, installed at the University of Alberta, Canada. It separates a 5&50wt% methanol-water feed gated, in this application, it was generally found that the advantages to be gained from the use of high order mixture, which is introduced at a rate of 18.23gs-’ into estimators were marginal. the column on the 4th tray. Bottom and top product compositions are to be maintained at 5wt% and 95wt% To provide control, fixed parameter PI controllers were used and their settings determined using Zieglermethanol respectively. A schematic diagram of the pilot plant column is shown in Figure 13. The column is Nichols ultimate sensitivity tuning rules, followed by fine-tuning to obtain responses with minimum integral of modelled by a comprehensive set of dynamic heat and mass balance relationships37. Both the pilot plant and the absolute error. Comparative performances of PI feedback control and adaptive inferential control, under non-linear model have been widely used by many investigators to study different advanced control schemes36*3a. major disturbance conditions, were assessed by subjecting the column to f 25% feed flow changes from steady The control objective is to use steam flow rate to the state operating conditions. This disturbance sequence reboiler to regulate bottom product composition, when the column is disturbed by step changes in feed flow rate. consisted of a 25% increase from steady state; a 25% Bottom product composition is measured by a gas chrodecrease back to steady state; a 25% decrease from matograph, which is subject to an analyser delay of 3 steady state; and finally a 25% increase back to steady state. min. The liquid temperature of tray 3, i.e. the tray immeThe open loop responses of bottom product diately below the feed tray, was chosen as the secondary zyxwvutsrqponmlkjihgfedcbaZYXWVUTSRQPONMLKJIHGFEDCBA 12 .I.Proc. Cont. 1991, Vol I, January Soft-sensors for process estimation and inferential control: M. T. Tham et al. 3 _ r n_, z I$ 5_ JT 300 _ zyxwvutsrqponmlkjihgfedcbaZYXWVUTSRQPONMLKJIHGFEDCBA _ (min) TJ I I 03 ” zyxwvutsrqponmlkjihgfedcbaZYXWVUTSRQPONMLKJIHGFEDCBA zyxwvutsrqponmlkjihgfedcbaZYXWVUTSRQPONMLKJIHGFEDCBA 300 (min) Methanol distillation column bottom product composition (PI control of analyser composition) f 25% feed flow rate disturbances Figure 14 I- Figure 16 Methanol distillation column bottom product composition (Adaptive inferential control using reboiler temperature, f 25% feed flow rate disturbances) I I 300 (min) Figure 15 Methanol distillation column bottom product composition. (Adaptive inferential control using tray 3 temperature, *25% feed flow rate disturbances) sition and tray 3 liquid temperature to f 10% step changes in steam and feed flows (not shown), indicate that both bottom product composition and tray 3 temperature exhibit quite different gains and time constants for both positive and negative steps in steam and feed flow rates. Indeed, the dynamics of the column become severely non-linear especially when it is subject to large (> f 20%) changes in reflux, steam or feed flow rates. The linear estimator therefore has to cope with the resultant non-linearities of all these responses. This is an important factor when assessing the adaptation properties of the estimator to changing process characteristics. Note that in figures to be presented, the composition plots correspond to ‘actual’ bottoms composition, y,(t) rather than its analyser delayed value, v(t). The performance of the PI controller, using the measured product composition as the feedback signal is shown in Figure 14, while Figure 1.5 shows the performance of the adaptive inferential control scheme. Deviations of actual composition from its steady state value (full line) are plotted along with the estimated composition (broken line). Particular note should be made of the difference in ordinate scales in Figures 14 and 15. Clearly, an overall improved transient behaviour is obtained with the adaptive inferential control strategy. The PI settings used in the conventional control strategy were K= - 60 and T,= 20 min. It may be argued that the proportional gain is too high, hence the results in Figure 14. However, reducing the gain of the PI controller, and increasing the integral time to account for the measurement delay, would result in a more sluggish response to that shown in Figure 14. Furthermore, there was an increase in peak overshoot of about 0.6 wt% (response not shown). The use of the estimates for feed- back control can be regarded as implicitly providing for dead-time compensation in the closed loop, and thus allows the use of a controller with higher gains. This was indeed the case, and the PI settings used in the inferential scheme were K= - 175 and T,= 2.5 min. As a result, disturbance rejection responses are superior to those obtained using conventional feedback control. It is also suggested that fine tuning of the PI settings in the adaptive inferential control scheme could well provide for further improvements. It is noted that although the estimator tracks the actual composition response closely during the first part of the transient, following the + 25% increase in feed rate, the actual composition response exhibits slow offset removal as the estimator parameters tune in under low excitation conditions. This is particularly apparent during the two major disturbances occurring at times 20 and 160 min, respectively. In practice, a continuously repeating sequence of disturbances, such as used in this study, would not be expected to occur. Therefore, a longer and richer adaptation period would be available. The use of various tray temperatures as the secondary output have also been investigated. Tray 3 temperature was selected as the secondary variable to demonstrate the performance of the estimator (soft-sensor) based inferential control strategy using a ‘badly’ chosen secondary measurement. If the secondary measurement was made closer to the sampling point of the primary output, then the performance of the estimator would be significantly enhanced. This is illustrated in Figure 16, which shows the performance of the adaptive inferential control scheme when the temperature of the liquid in the reboiler was used as the secondary variable. As expected, the estimated composition values (broken line) accurately reflects the state of the actual composition (full line). Other comparative studies have been performed using the multicomponent distillation column model described in Patke et ai.* A much improved, offset free control performance, has been reported using the estimators described in this paper2”26. The ability to anticipate disturbance effects results in a significantly improved control behaviour. This has been achieved here despite a poor choice of secondary measurement location. The studies also hopefully indicate the potential of these estimators in providing for improved product quality regulation of multicomponent columns. It is also worth noting that in J. Proc. Cont. 1991, Vol I, January 13 Soft-sensors for process estimation and inferential control: M. T. Tham et al. Chemicals and Polymers, the UK SERC and the Univerpractice, it is unlikely that all available secondary outsity of Newcastle upon Tyne. puts will exhibit highly non-linear dynamics. Neither would the disturbances affecting the plant change as drastically as the steps used in this study. Hence the References ability of the adaptive algorithm to track slow changes in plant dynamics should be adequate. As has been shown, 1 Soderstrom, T, in ‘Methods and Applications in Adaptive Control’, (Ed. H. Unbehaven), Springer-Verlag, 1980 the adaptive inferential control technique described 2 Parrish, J. R. and Brosilow, C. B. Automatica 1985,21 (S), 527 above appears to possess such capabilities. zyxwvutsrqponmlkjihgfedcbaZYXWVUTSRQPONMLKJIHGFEDCBA 3 Morari, M. and Stephanopoulos, G. AIChE Journal 1980, 26 (2), 247 Conclusions In this paper, two adaptive estimators for inferring infrequently measured plant outputs, using more frequently measured secondary variables, have been derived. In one approach, a general input-output representation of the plant formed the starting point, while the other utilizes a state space model. These algorithms have been referred to as ‘soft-sensors’ since they are software based rather than being hardware based. Both estimators approximate the complex relationship between the two sets of output variables. The parameters of the estimators are identified on-line using values of the infrequently sampled plant output and the relatively faster sampled secondary outputs and inputs. The technique enables the controlled output to be inferred at the same rate at which the secondary output is being sampled. The algorithms have been successfully applied to several industrial processes. Fermenter biomass concentration was estimated using off-gas CO2 and fermenter dilution rate measurements; top product composition of a demethanizer column was estimated using column temperature and reflux flow rate measurements; and polymer melt flow index was estimated using measurements of reactor feed rate, reactor coolant flow rate and hydrogen concentration above the reaction mass. The results highlight the potential of the proposed adaptive estimation schemes for solving ‘difficult’ industrial control problems caused by instrumentation limitations. With regard to closed loop ‘adaptive inferential’ control, early detection of load disturbances effects (in a feedforward sense) permits the achievement of significantly faster disturbance rejection. Consequently, it is expected that control schemes using the estimated output in conjunction with fixed parameter PI or PID controllers will out-perform other conventional feedback PI/ PID controllers, or indeed self-tuning control schemes, which rely solely on infrequent output measurements. This potential for adaptive inferential control was also demonstrated. 4 Morari, M. and Stephanopoulos, G. Inf. J. Control 1980,31(2), 367 5 Morari, M. and Fung, K. W. Computers and Chemical Engineering 1982, 6 (4), 271 6 Cwiklinski. R. R. and Brosilow. C. B. Proc. JACC 1977. 1530 7 Joseph, B.,‘and Brosilow, C. B. AIChE Journal 1978,&l (3), 485 and 500 8 Patke, N. G., Deshpande, P. B. and Chou, A. C. Ind. Eng. Chem. Process Des. Dev., 1982, 21, 266 9 Luyden, W. L. Ind. Eng. Chem. Fundam 1973,12 (4), 463 10 D’Hulster. F. M. and van Cauwenberahe. A. R.. Proc. 2nd World Congress on Chemical Engineering, M&meal, Canada, 1981,378 11 Wright, J. D., MacGregor, J. F., Jutan, A., Tremblay, J. P. and Wong, A. Proc. JACC 1977, 1516 12 Montague, G. A. Morris, A. J. and Bush, J. R. IEEE Confrol Systems Mug. 1988, 44 13 Bastin, G. and Gevers, M. R. IEEE Trans. Auto. Conlrol 1988, 33 (7), 650 14 Dochain, D. and Bastin, G. Automatica 15 Dochain, D. and Bastin, G. Automatica 16 Halme, A. and Seklainaho, A. Proc. 6th. IFAC Symp. on Identification and System Parameter Estimation, Washington, USA, 1982 17 Halme, A. Proc. 1st. IFAC Symp. on Modelling and Control of Biotechnological Processes, Noordwijkerhout, Holland, 1985, 179 18 Montague, G. A., Morris, A. J., Wright, A. R., Aynsley, M. and Ward, A. C. Proc. IEE, Part D 1986,133 (5). 240 19 Montague, G. A., Morris, A. J., Wright, A. R., Aynsley, M. and Ward. A. C. Can. J. Chem. Enc. 1986.64.567 20 Nihtila, M., Harmo, P. and P&tula,’ M: Proc. 9th. IFAC World Congress, Budapest, Hungary, 1984 21 San, K. Y. and Stephanopoulos, G. Biorech. & Bioeng. 1984, 26, 1189 22 Stephanopoulos, G. and San, K. Y. Biotechnol. & Bioeng, 1984,26, 1176 23 Swiniarski, R., Ixsniewski, A., Dewski, M. A. M., Ng, M. H. and 24 25 26 2-l 28 29 30 31 32 33 34 35 36 Leigh, J. R. Proc. 1st. IFAC Workshop on Modelling and Control of Biotechnical Processes, Espoo, Finland, 1982, 23 1 Guilandoust, M. T., Morris, A. J. and Tham, M. T. Ind. Eng. Chem. Res. 1988, 27, 1658 Guilandoust, M. T., Morris, A. J. and Tham, M. T. Proc. IEE Part D 1987, 134 (3) 171 Guilandoust. M. T. and Morris. A. J. Proc. 35th Canadian Chemical Eng. Conf., Calgary, Canada, 1985, p. 213 Kreisselmeier, G. IEEE Trans. Auto. Control 1977, 22, 2 Youna. P. C. Electronic Letters. 1979. 15 (12). 358 Tuffs p. S. and Clarke, D. W. &oc. IBE, Par: D, 1985,128 (3), 93 Astrom, K. J. in ‘Introduction to Stochastic Control Theory’, Academic Press, London, UK, 1970 Lu, W., and Fisher, D. G. Int. J. of Control, 1988,48 (I), 149 Anderson, B. D. O.and Moore, J. L. in ‘Optimal Filtering’, Prentice Hall Publishers, Hemel Hempstead, Herts, UK, 1979 Ljung, L. and Soderstrom, T. ‘Theory and Practice of Recursive Identification’, MIT Press, London, UK, 1983 Tsay, Y. T. and Shieh, L.S. Proc. IEE, Part D, 1981, 128 (3), 93 Isermann, R. in ‘Digital Control Systems’, Springer-Verlag, FRG Morris, A. J., Nazer, Y. and Wood, R. K. Optimal Control and Applications Acknowledgements The authors would like to acknowledge support from ICI Biological Products Business, ICI Engineering, ICI 14 J. Proc. Cont. 1991, Vol 1, January 1984,20,621 1986,22, 705 1982, 3, 363 37 Simonsmier, U. ‘Non-linear Binary Distillation Column Models’, MSc Thesis, University of Alberta, Canada, 1977 38 Morris, A. J., Nazer, Y., Wood, R. K. and Lieuson, H. Proc. IFAC Conf. on Digital Computer App. to Process Control, (Ed. Isermann and Haltenecker), Dusseldorf, FRG, 1981, p. 345