CONTROL LOOP P E R F O R M A N C E ASSESSMENT A N D OSCILLATION DETECTION By Lahoucine Ettaleb B. Sc. Ecole Hassania des T. P. Casablanca, 1987 M . Sc. A. Ecole Polytechnique Montreal, 1994 A THESIS SUBMITTED IN PARTIAL F U L F I L L M E N T OF THE REQUIREMENTS DOCTOR FOR T H E DEGREE OF OF PHILOSOPHY in THE FACULTY OF GRADUATE STUDIES DEPARTMENT OF ELECTRICAL A N DCOMPUTER ENGINEERING We accept this thesis as conforming to the required standard T H E UNIVERSITY OF BRITISH C O L U M B I A July 1999 © Lahoucine Ettaleb, 1999 In presenting this thesis in partial fulfilment of the requirements for an advanced degree at the University of British Columbia, I agree that the Library shall make it freely available for reference and study. I further agree that permission for extensive copying of this thesis for scholarly purposes may be granted by the head of my department or by his or her representatives. It is understood that copying or publication of this thesis for financial gain shall not be allowed without my written permission. Department of ^lg-^VW <~»X gu^«k The University of British Columbia Vancouver, Canada DE-6 (2788) Cs^^ovW— wxWi^ Abstract In this thesis an analysis of process control loop performance based on the Harris's performance index is given. Good and poor performances will be defined and it will be shown that poor performance may be obtained when some loops are cycling. Then a procedure for localizing the oscillating loops in a cascade multiloop system will be introduced. The procedure uses data collected under normal operation and assumes knowledge of the different frequencies of the periodic component of the output signal. An oscillation index is introduced to characterize the oscillating loops. The estimation of time delay is a major step in assessing control loop performance. To overcome some difficulties of the available methods for its estimation, a new off line extended Least-Squares technique for parameter identification and time delay estimation will be introduced for processes where acting disturbances are white noise (i.e. C(q^ ) = 1). 1 The Harris performance index for single-input single-output (SISO) systems is a very useful tool, therefore its extension to multivariable systems becomes necessary. The multi-input multi-output (MIMO) performance indices available in the literature involve use of the MIMO time delay matrix, also called the interactor matrix, to perform the assessment. The method proposed in this thesis does not require knowledge of the interactor matrix. The method uses knowledge of the delays between different input/output pairs of the process and defines an absolute lower bound on the achievable output variance for each output. This bound is then used to define a performance index associated with each output. The sum of these bounds is used to characterize the overall control loop performance. The results are independent of the matrix time delay representation, or interactor. This method is evaluated by application to the loops controlling a lime kiln in a Kraft pulp mill. ii Table of Contents Abstract ii List of Tables vi List of Figures vii Acknowledgment 1 2 3 ix Introduction 1 1.1 Control Loop Performance Monitoring 1 1.2 Outline of this thesis 8 Control Loop Performance Monitoring: SISO case 11 2.1 Minimum Variance Control (MVC) 11 2.2 Evaluation of performance 15 2.2.1 Evaluation of I using time series 15 2.2.2 Evaluation of I using Laguerre network 18 2.2.3 Time delay estimation 19 v p A n Extended off-line Least-Squares Method and Time Delay Estimation 3.1 3.2 3.3 Introduction 22 22 Least squares parameter estimation 24 3.2.1 Least-squares method 24 3.2.2 Generalized least-squares 26 3.2.3 Extended least-squares 29 An ARMA model for a closed-loop system iii 30 3.4 A n off-line extended least-squares method 34 3.5 Off-line time delay estimation 39 3.6 3.7 4 5 Case of colored noise 45 Conclusion 46 Realistic Performance Assessment in the SISO Case 47 4.1 Introduction 47 4.2 Problem formulation 47 4.3 Satisfactory and unsatisfactory performance 50 4.4 Control chart 54 4.5 Tuning procedure 55 4.6 Simulation Example 4.7 Conclusions . 57 60 Monitoring Oscillations in a Multiloop System 61 i 6 5.1 Introduction 61 5.2 Background 63 5.3 Monitoring the oscillations 65 5.4 Accuracy of the Describing Function Approach 72 5.5 Some Limitations 73 5.6 Simulation example 73 5.7 Conclusion 83 Control Loop Performance Monitoring in the M I M O case 6.1 6.2 Assessment strategies based on MIMO minimum variance control (MVC) Ultimate output variance bound 84 . . . . 85 88 6.3 MIMO Performance Indices 91 6.4 Estimation of the lower bounds 95 6.5 Simulation examples 98 iv 6.6 7 8 6.5.1 First example 6.5.2 Second example 98 101 Conclusions 104 Performance Assessment of Control Loops on a Lime K i l n 106 7.1 Process description 106 7.2 Lime kiln process 108 7.3 Loop performance assessment for the lime kiln 112 7.4 Conclusion 115 Summary 116 8.1 Contributions of this thesis 116 8.2 Recommendation for future work 117 Bibliography 118 v List of Tables 3.1 Results 43 6.2 C i = 0.9 <TI = 1 CT = 1 6.3 Reference values for example 2 95 2 101 vi L i s t of Figures 2.1 Process control 12 2.2 Process control 13 3.3 Covariance of the residual 33 3.4 White noise and its estimate 34 3.5 Estimation value of the coefficient a 41 3.6 Estimation value of the coefficient b 42 3.7 Estimation value of the time delay d 42 3.8 Loss function V 43 3.9 Estimation value of the coefficient a by ordinary least-squares 44 3.10 Estimation value of the coefficient b by ordinary least-squares 44 3.11 Estimation value of the delay d by ordinary least-squares 45 4.12 Process control 49 4.13 I and steady state error versus P-gain 49 4.14 A single nonlinear loop 51 4.15 Saturation characteristic 51 4.16 Nyquist Diagram 52 4.17 Partition of Ip curve 52 4.18 Performance analysis 54 4.19 Process approximation 55 4.20 Process control 56 4.21 Process control 57 4.22 Comparison of true and approximate model 58 p vii 4.23 Estimated curves versus Controller parameters 59 4.24 Estimated curves versus Controller parameters 60 5.1 Process control 63 5.2 Process control 66 5.3 Process control 73 5.4 The output y(t) of the outer loop, example(5.6.1)-test 1 75 5.5 The output yi(t) of the inner loop, example(5.6.1)-test 1 76 5.6 The output out\(t), example(5.6.1)-test 1 76 5.7 The output outfit), example(5.6.1)-test 1 77 5.8 The output y{t) of the outer loop, example(5.6.1)-test 2 78 5.9 The output y\(t) of the inner loop, example(5.6.l)-test 2 78 5.10 The output out\(t), example(5.6.1)-test 2 . 79 5.11 The output out2(t), example(5.6.1)-test 2 79 5.12 The output y(t) of the outer loop, example(5.6.2) 80 5.13 The output yi(t) of the inner loop, example(5.6.2) 81 5.14 The output outi(t), example(5.6.2) 81 5.15 The output out (t), example(5.6.2) 82 5.16 The output outfit) versus iri2(t), example(5.6.2) 82 2 6.17 On line performance index estimation kl=0.08 and k2=0.05 104 7.18 Lime cycle 107 7.19 Lime Kiln 108 7.20 First bump test 110 7.21 Second bump test Ill 7.22 Validation of the model 112 7.23 On-line performance indices 114 viii Acknowledgment I am grateful to a number of individuals who greatly assisted me throughout the course of this thesis. Without their support this work would not have been be achieved. In particular I wish to thank: Both Dr. Michael S. Davies and Dr. Guy A. Dumont my thesis supervisors for giving me the opportunity to work under their guidance. The research reported in this thesis reflects true interactions, inspiration and encouragement from both of them. I am very grateful for their direction and their financial support. Dr. Ezra Kwok for his valuable suggestions and help with my research work. Rita Penco, the magic Librarian for providing any required documentation right on time. Georgina White, Ken Wong, Brenda Dutka , Lisa Brandly, Peter Taylor and Tim Paterson for their outstanding assistance. Thanks are also due to Brian D. McMillan for computer support and to my student colleagues and friends at the Pulp and Paper Centre. Ministere des T.P. of Morocco, Natural Sciences and Engineering Research of Canada and Networks of Centres of Excellence of Canada for providing the financial support. Dr. Abdelhaq Elhakimi, Dr. Mohamed Y. Tachafine, my colleagues at Ecole Hassania des T.P. for their moral support. John Ball from Canfor pulp and paper mills and Bruce Allison from Paprican for providing lime kiln data used in chapter 6. My wife Fatima Elmourabit, my children: Sarah, Fatima Azzahra and Hamza for their unconditional love and unselfish support. My large family, brothers, sisters and friends back in Morocco. I pay special tribute to Sandra Davies for her kindness and care. She was always there when myself or my family needed some one to help. I am indebted to remain grateful to her. ix Finally, it is to my mother Hajja Zahra Irramed that I dedicated this thesis for her eternal love and encouragement. x Chapter 1 Introduction 1.1 Control Loop Performance Monitoring The purpose of on-line surveillance of control loops is to obtain information on the performance of plant components while the plant is operating, without disturbing its normal operation. This information is used for a number of purposes such as 1) to provide information on the performance deterioration with time and its correlation with operating conditions so as to help the diagnosis of a problem; 2) to allow improved control of normal operating conditions so that the plant can run closer to its limits. The analysis of this information guides decisions on what repair work or other modifications will be needed at the next shutdown and whether continued operation is safe under the normal conditions. On-line surveillance is a general concept, and research in the area is widely treated in the literature because the problems addressed are as old as the process industry itself Mah [53]. When process measurements and data acquisition were regularly carried out by plant operators and engineers, the accuracy and consistency of the data were checked on the basis of their knowledge and experience of the process. Now these functions are performed automatically by computers; hence, a far greater volume of process data may be checked and processed without requiring of human vigilance and intervention. The surveillance in a wide sense includes many activities such as performance monitoring, fault detection and state or parameter estimation. A common characteristic among these techniques is the use of collected data. But the objective of each method is different. Performance monitoring attempts to state whether satisfactory system performance standards are met. If not, then fault detection procedure is to be initiated to locate the faulty component. If there 1 Chapter 1. Introduction 2 is none, then an engineering modification must be undertaken. Degradations in control performance can be caused by an inaccurate sensor reading due to bias, or an actuator which is not responding properly to control signals due to a non-linear phenomenon such as backlash or stiction, or any leaks in pipelines. These malfunctions are refered to as failures or faults. In some cases these faults can be modeled as an additive disturbance, and therefore they can be detected using Kalman filter-based state estimators. Observing changes in the mean and variance of the main process variable can also indicate that a failure has occured because faults introduce changes in the properties of system components. Taking this fact into account, some recent methods of fault detection are based on that principle and many detection algorithms using recursive parameter estimation are now available in the literature as detailed in [8]. A review of the most significant methods for fault detection is given in [46, 7]. Fault detection and parameter estimation are two wide areas of research well treated in the literature. The main focus in this thesis will be on control loop performance monitoring based on the use of the minimum variance controller as benchmark. Control loop performance has recently received great attention from researchers and manufacturers who realize that it is a key to manufacturing high quality product. In the past, statistical such as Shewart control charts were used. Recently, focus has been placed on the development of multivariate statistical control charts[55]. Principal Component Analysis (PCA) and Partial Least Squares are among the best known methods for processes with a large number of measurements. The idea behind these methods is to compress the highly correlated information contained in the data into a lower dimensional space which allows use of univariate control charts. These charts can also be generated using methods from chemometrics[69]. As mentioned before, the ultimate goal of performance monitoring is to improve product quality by reducing process variation, Smith [60]. The most commonly used measure of performance is the variance or standard deviation of key process variables. Standard deviation is used for monitoring because reduction of this Chapter 1. Introduction 3 quantity implies control and manufacture of a uniform product in the face of randomly varying raw materials and chemicals[9]. It was shown in Astrom et al. [3] that by using the concept of minimum variance control (MVC), it is possible to control the process variability, so that the remaining variations are due only to a random noise which, when dead time is present, may be correlated. In some cases it is not possible to implement MVC, and time series analysis tests exist which identify whether minimum variance has been approached. Pryor [57] discusses some aspects of the use of autocorrelation and power spectra to analyze paper machine dynamic data. The analysis of control performance monitoring from dynamic data was also addressed by Staton [63]. The approach of Harris [35] for control loop performance assessment is based on the use of routine operating data and little prior knowledge of the process. Although the control objective may not be to minimize process output variance cr , Harris proposes the use of the 2 minimum achievable output variance as a benchmark to assess the performance of a regulator. To determine the minimum variance <T^„ , the process time delay of Td sampling intervals must be known and a model of the closed loop output data must be estimated. An A R M A model of the form y(t) = H(q~ )e(t) is used. The first Td terms of the impulse response of the filter H 1 (hi,..,hT -i) d are feedback-invariant and are used to estimate c r ^ as 2 a 2 = (l + h + ... + h _ )al 2 mv 2 1 Td 1 (1.1) This last equation requires knowledge of the process time delay Td- The estimation of this entity is a well known problem in the literature and, many methods have been developed, among them the techniques described in Zheng et al.[71] or in Elnaggar and Dumont [22]. The method of the last reference is used in many applications and proved to give good results[52, 30, 56]. But, it has also some limitations[52]. In the third chapter, for processes where the acting disturbances are white noise i.e. the polynomial C(q~ ) = 1 a new method will be proposed. l Other performance measures have also been used in the literature. Notable is that used by Astrom [2] and Astrom et al.[6] for assessing achievable performance using PID control. Harris's technique has been used successfully to improve quality and uniformity of products in many applications, for example [56, 52]. Unfortunately, this method has some limitations. Chapter 1. Introduction 4 One of them is that the Harris's performance index is only valid for linear minimum phase systems. Tyler and Morari[66] proposed an extension to a nonminimum phase systems by taken into account unstable system poles and zeros. Another limitation is that when comparing the output variance to the minimum output variance cr^v Y the extreme situations are well defined: good variance performance when o m the performance index (ratio of the actual output variance to the minimum achievable output variance noted I ) is equal to one, and poor performance when the index is large. In practice, p the index takes on values between these extremes and it is necessary to determine whether the actual controller performance is acceptable, acceptable with a small tuning correction, or so unacceptable that redesign of the controller is needed. Another issue related to the minimum variance benchmark limitations is whether the performance index itself is appropriate. In practice, the ability of the controller to reduce the effect of noise on the output variance is only one of several desirable attributes. Therefore, having I ~ 1 for a certain controller may be too narrow as assessment of performance. In similar fashp ion, I = 1 is obtained when implementing a minimum variance control strategy (MVC) for a p minimum phase system. As discussed by Isermann [45], minimum variance strategy may result in slow output tracking of setpoint changes if the process zeros are too slow. This strategy may therefore, not always be useful. To overcome some limitations of the minimum variance benchmark, Kozub and Garcia[49] propose a performance index which allows an exponential decay of the closed-loop response. Horch and Isaksson[39] introduce a modified performance index which does not require more knowledge than the original one. Huang and Shah[41] propose a "user-defined benchmark" which requires knowledge of the disturbance model. They assume that the disturbance model is either an integrator or has to be estimated from data. However, a good control strategy is often obtained when compromising between output variance as a performance criterion and other measures such as steady-state rejection of disturbances. These points will be discussed in more detail in chapter 4 where it is also shown that poor performance may result when loops are cycling. Oscillations have an impact on product Chapter 1. Introduction 5 uniformity, and can cause increase energy consumption and waste of raw material. In a recent paper Clarke [16] showed that by eliminating oscillations in control loops, 6% saving in the steam consumption of a batch temperature controlled process was possible. Bialkowski [9], in a survey, indicated that 30% of all control loops in Canadian paper mills were oscillating because of sensing and valve nonlinearity problems. Most of the valve nonlinearities encountered in the pulp and paper industry are caused by backlash or stiction. Backlash is due to the difference in motion response between an increasing and a decreasing output usually caused by mechanical gearing. Backlash does not typically result in cycling, except in integrating processes, but it does degrade control performance by limiting the loop's ability to attenuate load disturbances and respond to set point changes. Stiction, on the other hand, affects the valve control element and is more severe than backlash, almost always resulting in loop cycling[68]. Stiction causes the actuator force to continue to increase without a change in the valve position until the force overcomes the friction force. Then, the valve will move suddenly to a new position. Stiction causes reduced valve resolution and results in a limit cycle due to the controller's inability to position the valve to maintain a setpoint[68]. It is important to mention that although these valve problems cannot be eliminated by controller tuning, they can however be minimized. But only repairing or replacing the valve can eliminate them. Oscillations may also be due to bad control tuning. Dealing with detection of oscillations due to non-linear problems is a subject of many articles in the literature. In Horch and Isaksson[38], extended Kalman filter was used to generate two sets of residuals: One describing a valve with an acceptable level of friction and another describing a valve with unacceptable amount of stiction. A generalized likelihood ratio is then used to decide which of the sets of residuals is more likely to correspond to the measured data. In Taha et al.[65] a friction index is proposed which measures the deviation of the valve's inputoutput characteristic in normal operating conditions from its static characteristic established by the manufacturer. The index measures the degradation of the valve due to its use. In chapter 5, the oscillating loops are first located and then the cause of the oscillations is found. The Chapter 1. Introduction 6 procedure uses prior knowledge of different frequencies present in the periodic component of the output signal. A n oscillation index will be introduced to characterize the oscillating loops. The prior knowledge of oscillation in the loops can be obtained by techniques such as the Discrete Fast Fourier Transform (DFFT) or Singular Value Decomposition (SVD)[47]. The Harris method has been extended to Multi-Input Single-Output (MISO) in Desborough et al.[18] and very recently to Multi-Input Multi-Output system (MIMO) by Huang et al.[42] and Harris et al.[36]. The last two references use minimum variance control as a benchmark. In fact with no explicit time delay structure, the key steps in the derivation of SISO and MIMO minimum variance control strategies are: • A predictor for the time shifted variable £(q)y(t), where £(q) stands for the matrix time delay representation, is developed. • The noise term ^(q)N(q~ ) is factored into future and past noise terms as 1 Z^Niq- ) 1 where F(q) = F q d d = F(q) + Riq' ) 1 (1.2) + ... + F q and ^ ( g ) is a proper transfer function. -1 t • Under M V control, and in steady state, the tracking error is shown to be a finite moving average polynomial of the form y(t)-y*(t) = cHq)F(qHt), (1.3) where {y*{t)} is a given reference trajectory and £(q) is the time delay representation considered. Note that the expression for this tracking error (equation 6.183) does not depend on controller parameters. For SISO systems £~ {q) l = q~ , therefore equation (1.3) reduces to equation (1.1) which Td expresses the minimum achievable output variance for SISO systems. Following from SISO assessment, one might be tempted to use the predicted performance of MIMO minimum variance controller as an absolute lower bound against which current performance can be assessed. A bound is indeed obtained when taking the variance of both sides Chapter 1. Introduction 7 of equation (1.3). However, this lower bound will depend on the type of matrix time delay representation used and so is not unique. The MIMO time delay representation was first proposed by Elliott et al. [21], in 1982, who showed that the notion of delay in the MIMO case is related to the system interactor matrix £(g) introduced by Wolovich et al.[70], in 1976, which satisfies the following two properties for a given system with the transfer function T: 1. l i m _ £ ( q ) T ( q ' ~ ) = K for any non singular matrix 1 g 00 2- £(<?) = H{q)dia,g{q , ...,q ) where H(q) is a unimodular matrix. dl dn The interactor matrix was introduced initially by Wolovich with no apparent justification other than the resulting computational procedure for its calculation. It was Elliot and Wolovich [21] who showed the usefulness of this entity and related it to time delay in the MIMO case. The interactor proposed by Bittanti et al.[ll], in 1994, satisfies (1.4) S(q)Z(q- ) = i, 1 T and also satisfies property (1) above. This interactor is called the spectral unitary interactor matrix. As can be deduced for the same system one can calculate at least two inter actors. The procedure for control loop performance assessment proposed by Huang[40] uses the unitary interactor. In Harris et al. [36], the procedure is based on the interactor - 7(9)^(9) where £(q) is the lower triangular form as defined by Wolovich, and 7(g) is obtained by solving the following spectral factorization equation b W f Q i W r 1 ) ] = m-TQiM*- )- ] 1 1 (1.5) where Q\ is a diagonal matrix. In fact, it is easy to show that (1.5) is a more general form of equation (1.4). Writing (1.5) as £(9)Qi£VT = Qi (1.6) Chapter 1. Introduction we let Q i = I then f ( g ) f ( g ) _1 8 T = I. If Qi # J and |(g) satisfies (1.4) i.e. | ( g ) | ( g ) _1 T = J, then by writing Q i = QI QI the interactor matrix £'(q) = Ql£{q)Qi h (1.7) is also a solution to (1.6). On the other hand if £'(g) satisfies (1.6) then the interactor £(g) defined in (1.7) satisfies (1.4). Therefore (1.6) serves as a link between (1.4) and (1.6), showing the equivalence of the inter actors. Knowledge of the interactor matrix is a drawback in the approaches of Huang and Harris. Though Huang et al.[43] showed that only a few Markov parameters are needed for its determination. If the assessment is meant to be on-line the estimation of the first few Markov parameters is not straightforward for the following reasons: • Most processes that require monitoring are under feedback control. • For most processes their linear model depends on the operating point which keeps changing therefore, the interactor matrix also does. For these and other reasons many authors avoid use of the interactor matrix for adaptive control (see for example Dugard et al.[19], Shah et al.[58], ...etc). However, the approach which is proposed in chapter 6 is based on defining an absolute lower bound on the achievable value for each output variance, thus avoiding the need to identify the interactor matrix. It also removes the ambiguity in the lowest achievable performance. 1.2 Outline of this thesis As mentioned before, this thesis is concerned with control loop performance monitoring based on minimum variance control as a benchmark. The research accomplished will be presented in six chapters. In the second chapter, the Harris approach for control loop performance monitoring is described along with different techniques for its implementation for single-input single-output systems. The Harris method requires knowledge of process time delay, therefore, in the third Chapter 1. Introduction 9 chapter, an off-line extended least-squares method to identify the parameters of a linear system under closed-loop conditions for which the acting disturbances are white noise (i.e. C{q~ ) = 1) l is proposed. An off-line time delay estimation method, as an application of the proposed method, will also be presented . These results have been reported in [24]. For some applications the Harris performance index does not lead to a clear cut assessment. In the fourth chapter, an analysis of process control loop performance based on Harris's performance index is given. Good and poor performance will be defined and it will be shown that poor performance may be obtained when some loops are cycling. A controller tuning procedure is given that is based only partly on output variance. These results have been reported in [27]. The important result from this third chapter is that poor performance results mainly from the fact that some loops are cycling. The fifth chapter, describes a procedure for localizing the oscillating loops in a cascade multiloop system. The procedure uses data collected under normal operation and assumes knowledge of the different frequencies of the periodic component of the output signal. An oscillation index is introduced to characterize the oscillating loops. These results have been reported in [26]. The detection and removal of oscillations caused by non linear effects is necessary if the linear methods of chapter 3 are to be used. For saturating gains the system can be brought into the linear region and then identified by controller adjustment. For more complex uncertainties it may not be possible to revert to linear operation. In the sixth chapter, a unified method for assessing single and multivariable control loop performance which does not require knowledge of the interactor matrix is presented. The method uses knowledge of the delays between different input-output pairs of the process and defines an absolute lower bound on the achievable output variance for each output. This bound is then used to define a performance index associated with each output. The sum of these bounds is used to characterize the overall control loop performance. The results are independent of the matrix time delay representation, or interactor. The application of the proposed method to an industrial process is presented in the seventh chapter, the proposed method is evaluated by application to the loops controlling a lime kiln in a Kraft pulp mill. In the last chapter Chapter 1. 10 Introduction conclusions and recommendations for future work will be given. The chronological use of the results of this thesis in a procedure to to improve control loop performance for a given system is as follows: 1. Oscillation detection: Investigate oscillating loops as a first step because oscillations cause poor performance and may lead to an incorrect calculation of the performance index. In the case of cascade loops, the results from chapter 5 are applied. 2. Performance indices calculation: Two cases are to be distinguished. • SISO systems for which methods from chapter 2 and chapter 3 are applied • MIMO systems for which results from chapter 6 are applied. 3. Interpretation of performance indices: For SISO systems, chapter 4 describes a realistic interpretation of the performance index and a tuning method that are based on the use of the performance index I po which can be estimated by the relay method. Chapter 2 Control Loop Performance Monitoring: SISO case The objective of this chapter is to describe the Harris single-input/single-output (SISO) approach for control loop performance monitoring and discuss the available methods for its on-line implementation. As mentioned before this technique relays on the use of minimum variance control as benchmark therefore the organization of this chapter is as follows. In the next section a brief discussion of minimum variance control will be given followed in section 2 by the definition of the performance index and how to determine it from operating data. The last section deals with time delay estimation. 2.1 Minimum Variance Control ( M V C ) The concept of minimum variance control was initially developed by Astrom et al.[3] to reduce process variability due to random disturbances. The derivation of M V C is outlined below. Let us consider the basic process control shown in figure (2.1). The linear system with delay is represented by the discrete time transfer function ^ G = ^&)^ < 2 8 > where A; is a positive integer, B and A\ are polynomials in q~ , the backward shift operator l (i.e q y ( i ) = y{t — 1)). The order UA of A\ is considered to be the same as or higher than n _1 B that of B. The disturbance model represented by ^ > = -£fk G <-> _1 2 The operator A is defined as A = 1 - q~~ . The polynomials A^q ) l -1 p and q respectively, and are stable. 11 9 and C(q~ ) are of order 1 Chapter 2. 12 Control Loop Performance Monitoring: SISO case u(t) + Gp(q') -f- y(t) o — Figure 2.1: Process control. The expression for the output is - m By long division ^^-IJA c a l ^ ^ w < + 2 - 1 0 ) be expressed in the following form (leading to a Diophantine n equation) . ntn-i\ . n(n~ \ l (2.11) where a unique solution is obtained if the polynomial F(q~ ) is monic and of order k — 1 1 (F(q~ ) = 1 + /i<Z +... + fk-iq~ ), 1 -1 k+1 and G(q~ ) is polynomial of order max[q — 1, p +rf— 1]. 1 Now, by using equations (2.10) and (2.11), the following expression gives the output at time t+k y ( t + k ) (^^< ) = t = + %V®) ( +Fe t + ) k y(t + k) + Fe(t + k). (- ) 2 12 (2.13) At time t + k, the best estimate of y(t + k) is y(t + k) because the term Fe{t + k) represents the information yet to come, as e is zero mean and e and y are uncorrelated. If the control law BA AFu(t) 2 + GAiy(t) = 0 (2.14) Chapter 2. Control Loop Performance Monitoring: SISO case 13 is used, then, from (2.12), y(t + k) = Fe{t + k) i.e. y(t + k) is a linear combination of the future values of the white disturbance sequence. It is not difficult to show that E[y (t + k)] is minimized in that case. The output variance is given 2 by a' mv = E[y\t + k)} = (1 + ft + ... + fU)<- (2.15) The control law defined by (2.14) represents the minimum variance control law. In closed loop form (Figure 2.2), the transfer function between the white noise e and the output y is given by Figure 2.2: Process control. (2.16) y(t) = i | r t(Bg«)(r') H{q~Mt). (2.17) By using the same manipulation as before one can easily derive the following expression of the Chapter 2. Control Loop Performance Monitoring: SISO case 14 output at time t rw x = / MG - BA FG A 2 H(q-i)e(t) \ C , , - (219) which shows that the first k Markov coefficients of H(q~ ) are feedback invariant. Then, it 1 follows that ol > o . (2.20) 2 mv Under (MVC) the variance of the output is minimized and it is given by 2 _ 2 which expresses the minimum output variance that can be achieved (limit of performance). The limit is independent of controller parameters. It depends only on the process noise model and the present system time delay. Therefore, the evaluation of a^ av gives insight into what performance may be expected from the actual system structure. If the performance is not satisfactory, then by referring to equation (2.15), change to the system structure can be made by either reducing the time delay by re-locating the sensors or reducing the disturbance variance. Equation (2.20) shows that minimum variance control is the best possible control in the sense that no regulator can have a lower variance. The connection between M V C and other strategies such as the PID controller and Dahlin's algorithm is discussed in Harris [37]. From equation (2.14) we can deduce that MVC control includes the inverse process model in the controller. It follows that for a non minimum phase process such a controller will result in an unstable control system. The problem can, however, be overcome by the technique discussed in detail by Astrom et al.[4]. It should also be noted that it is often undesirable to use M V C due to the excessive control action required. This matter will be discussed extensively in chapter 4. Chapter 2. Control Loop Performance Monitoring: SISO case 15 2.2 Evaluation of performance Even though the MVC cannot always be implemented, it still remains a convenient bound on achievable performance against which other controllers can be assessed. An index can be obtained by defining a performance index I based on the theoretical limit as p 1 < Ip < oo. Then the following comments can be made: • I ~ 1 : the actual performance is close to the best achievable performance. p • Ip « oo : the actual performance is significantly worse than optimum. To evaluate I one needs to identify the closed loop transfer function (equation 2.16) and p express it as a Diophantine equation to obtain the F(q~ ) required to determine the minimum 1 achievable output variance. I can then be computed [52]. An alternative to closed loop transfer p identification is to determine the transfer function between white noise and the output using Laguerre network[52] or using time series analysis [18]. These two methods are summarized in the following subsections. Evaluation of I using time series 2.2.1 p The procedure for evaluating the performance index I using time series analysis is detailed p in Desborough et al. [18]. This procesure is described below. It consists of fitting y(t) to an autoregressive moving average (ARMA) time series of the form y(t)-py = H{q- )e(t)l where p is the average of y(t). From the first section it was established (equation (2.18) that y y(t)-py = H{q- )e(t) l (2.21) (2.22) Chapter 2. Control Loop Performance Monitoring: SISO case = (2.23) iPi{q- )e{t) + i}>2{q- )e(t - k) l 16 l (2.24) Equation (2.24) is also equivalent to oo y(t) - ^ = e'(t) + oii(y(t - k - i + 1) - p ). y i=l In practice, the infinite series is truncated after m terms m y(t) - p y = e'(t) + oti(y(t - k - i + 1) - p) y i=l = m and, using matrix notation Y = Xa + e where y{n) y(n - 1) , a = Y = y(k + m) C*m y(n — k) y(n — k — 1) y(n — k — m + y(n — k — 1) y(n — k — 2) y(n — k — m) X = y{m) An estimate for y(m - 1) y(i) is given by the residual mean square error 1) Chapter 2. Control Loop Performance Monitoring: SISO case 17 By using the mean square error of y(t), an estimate of I is given by p s2 -(n — k — m + 1). (2.25) + f4 YY T If fiy is not known, it has to be estimated. For that, the closed loop model behaviour is modified to include a constant term an. The parameter vector a is augmented by cio a n ( i the matrix X is also augmented by a leading column of ones. These transformations give the following equations: Y = Xa + e (2.26) where y(n) OtQ y(n - 1) Y = , a= a.2 y(k + m) 1 y(n — k) 1 y{n — k —1) y(n — k- 1 ) ... y(n - k - 2 ) ... y(n- -k-m + l) y{n — k — m) X= 1 y(m) y(m - 1) y(i) The method described above gives an estimate off-line of Ip. The on line estimate can be obtained by using a recursive least squares algorithm to estimate the vector a in equation (2.26) and the following recursive equations to estimate the minimum variance performance as well as the moving mean square error of the output Smv(t) = S {t) = y \S v(t-l) m + (l-\)e\t) XS {t - 1) + (1 - A)y(t) y 2 Chapter 2. Control Loop Performance Monitoring: SISO 18 case where A is a forgetting factor (0.98 < A < 1) and e(t) is the prediction error at time t. An alternative method for evaluating I , to the one above is to use a Laguerre Network p Lynch and Dumont [52] as is summarized below. 2.2.2 Evaluation of I using Laguerre network p When using a state space representation of the relation between the white noise and the output x(t +1) = Ax(t) + be(t) (2.27) y(t) = Cx(t) + e(t), (2.28) the minimum achievable output variance is expressed in terms of the Markov parameters of the process Cb, CAb, CA b, 2 ... as -2 _ 2 l + k-1 (2.29) ^ipA^b) 2 i=l The same formula is still valid when using a Laguerre network instead of state space representation [52]. In this case, the system of equation (2.27)-(2.28) denotes the discrete Laguerre state space representation of Gunarsson [33], where the elements of matrix A and those of b and C are given by a,ij = a A = i = 3 = 0 i < j «u = (-«r (l-« ) i>3 j_1 bi = (-of-Vl-a C = 2 i = 2 ltoN [gi ... g \N The parameter a is the time scale, N is the number of Laguerre filters used (The optimal value of N is derived by Fu et al.[29]) and gw represents the N th filter xj^{q) is given by q-a \ q—aJ Laguerre gain. The N th Laguerre Chapter 2. Control Loop Performance Monitoring: SISO case 19 In the above analysis it was assumed that the time delay was known, but in most cases it is either unknown or time varying. The next section deals with the time delay's definition and estimation. 2.2.3 Time delay estimation The estimation of time delay is a familiar problem in control system design, in acoustics, and in radar remote control. Several good techniques have been developed for extracting this quantity from measurements. Most of those methods are based on analyzing the cross correlation function of u and y as the two signals of interest Ru (k) = E[u(t - k)y(t)} (2.30) V where E[.] denotes the expectation over time t to determine the time delay between u and y [48]. In some cases, when the signals are highly coloured, the cross-correlation function does not provide the correct time delay. This problem can be overcome by whitening the signals before performing the cross correlation Carter [14]. In Elnaggar et al. [22] the concept of variable regression estimation (VRE) to estimate the time delay is used. It is shown that this estimation can be added to any standard parameter estimation with minimal additional computation, which means faster convergence of the identification procedure. The model used to describe the process time delay is a first order given by _T a = —e -T s bq- ~ l+TS 1 + aq- e k d T 1 1 < b= 1 - e " then y(t) = ay(t - 1) + bu(t - k - 1) the prediction error e is defined as e(t,k) = y{t)-y{t, k). Chapter 2. Control Loop Performance Monitoring: SISO case where y(t, k) is the prediction of y given the estimate k. 20 The argument k that minimizes E[e(k) ] provides an estimate of time delay. 2 For b > 0 it is shown in [22] that this is equivalent to maximizing E\ given by Ex = Ru {k) - aRuy(k - 1). V The choice of a = 1 gives unbiased estimation of delay. E\ becomes E = E{{y{t + l)-y(t))u{t-k)}. 1 The on-line delay estimation procedure uses the following algorithm: • Ei(t, ki) = \Ei(t - 1, hi) + [y{t) - y(t - l)]u(t - ki) Vki e [k , k ] min max • choose the delay that maximizes E\ Ei(t, ki) = max Ei(t, ki) Vfcj € [k mini kmax\• In practice this procedure gives good results in open loop [52, 30, 56]. Unfortunately this technique has also some limitations [52]. In the next chapter a new method for time delay estimation will be proposed. Some methods for time delay estimation are based on Pade approximations which is often convenient as the approximation of the time delay is given in term of a rational transfer function instead of e~ . The approximation can be obtained from the following Pade approximation TdS of order n - s Td 6 „ 2 - T s + (-T s) /2 + ( - r a ) / 3 ! + .. + (-T s) /n\ ~ 2 + T s 4- {T sf/2 + (r a)3/3! + .. + {T s) /n\ ' 2 d n 3 d d d n d d d d The reverse problem of the Pade approximation is to obtain the equivalent time delay of a complex system described by \ TA + «i« + a s + .. 2 G(s) = K — 2 —5 . The equivalent time delay was derived by Matsubara [54] and is given by D = b\ — a\ n m i=l t=l Chapter 2. Control Loop Performance Monitoring: SISO case 21 where pi and Zi are respectively poles and zeros of the transfer function G(s). This powerful result was used by Isermann [44] to show that the following standard models k G l ( s ) = ke~ TdS (TT7^ G 2 ( s )= Gs(s) lTT7^ or k ES=i(i + ^ ) n can be used for describing a wide range of different low-pass models. These approximations or model reduction methods are useful for: • improving delay estimation. • Model reduction (reduction of parameter to be estimated). In the following chapter, an alternative method for time delay estimation will be presented, as an application of an off-line extended least-squares method. Chapter 3 A n E x t e n d e d off-line Least-Squares M e t h o d a n d T i m e Delay E s t i m a t i o n Most control loop performance assessment methods rely on some prior knowledge of the process under investigation. This knowledge may be the process time delay or the complete process model. Therefore, a major difficulty with these methods of assessment resides in the determination of the process model or the time delay. Because most performance assessment methods are expected to work with normal operating data, it is not clear that it is possible to identify the time delay in the absence of an external perturbation, e.g. setpoint change to ensure closed-loop identifiability In this chapter, it will be developed a technique that guarantees closed-loop identifiability of a process without using an external signal, as long as the process disturbance is assumed to be white. The case of colored disturbances although briefly discussed will not be studied here, this is reserved for future work. The organization of this chapter is as follows. In the following section, an introduction to closed loop parameter identification is given. In section 2, some existing variants of the least squares method are presented and an A R M A model is discussed for a regulatory closed loop system in section 3. In section 4, the proposed extended off-line least-squares method is presented. In section 5, the proposed method is used to determine the parameters and time delay of a system. Section 5 draws some conclusions. 3.1 Introduction In order to predict the future response of a dynamic process or to manipulate its outputs to follow desired trajectories, a mathematical model which describes the dynamics of the underlying 22 Chapter 3. An Extended off-line Least-Squares Method and Time Delay Estimation 23 process is necessary. In some cases this model can be obtained by analyzing the internal mechanism or physics of the system. In other cases the internal mechanism is not well understood or is highly complex, and an identification experiment has to be carried out. The model may be linear or non-linear though, in practice, linear models are more common since they are simple to analyze for identification, prediction or control purposes. Much work has been done in linear identification and control, and the theory is well known and documented. The estimation of the system parameters or equivalent of the process model can be performed in open-loop by many methods. A detailed description of these techniques is given, for example, in Ljung [51], Soderstrom and Stoica [62]. Under certain assumptions these methods can also be used to estimate the parameters of a process operating in closed-loop. If the process is open loop unstable, or the controller is adaptive, then it is preferable that the process model be found by closed-loop identification. Furthermore, it is now known that if the purpose of modeling is to design a controller, then it is advantageous to perform the identification in closed-loop[50]. For these reasons, closed-loop identification has become an important research area, and many methods dedicated to closed-loop identification are now available [28]. Most methods of closed-loop identification use a sufficiently exciting external signal added to the control signal or to the reference signal. In Van den Hof et al.[67], the sensitivity function is estimated first, then the process transfer function is determined. A recursive method is proposed by Landau and Karimi[50] using the sensitivity function and an external signal, but requires perfect knowledge of the controller. The use of an external signal to achieve identifiability under feedback conditions was proposed by Gustavsson et al.[34]. An alternative to using an external signal is to switch between different regulators as proposed by Soderstrom et al. [62]. The use of an excitation signal is the best known way to guarantee the estimation of the process parameters. However, for performance assessment, the use of such signals is rarely allowed. In that case, it will be shown in this chapter that in the case of a regulatory process Chapter 3. An Extended off-line Least-Squares Method and Time Delay Estimation 24 where the acting disturbance can be described by white noise, it is possible to estimate the process parameter without an external excitation signal. In the following an overview of least squares methods is given. 3.2 Least squares parameter estimation In the field of parameter estimation the least squares technique has reached a significant level of popularity and perfection. This technique was first introduced by Karl Friedrich Gauss in 1795 who used it for astronomical computations. Some variants of the least-squares method will be reviewed based on a survey article by Strejc[64]. More details and references are to be found in that article. The objective of this section is to show the difference between the proposed extended least-squares, which will be presented in the following section, and the existing leastsquares approaches. 3.2.1 Least-squares method Consider the problem of estimating the polynomials A(q~ ) and B(q~ ) of the system described 1 l by the following equation A{q- )y(t) = l B{q- )u{t)+e{t), l (3.31) where A(q~ ) = l Biq' ) = hq' l 1 q- l + ai 1 + ...+a q- , nA nA + ... + b q- nB nB (3.32) (3.33) y(t), u(t) and e(t) are respectively the process output, the process input and the process disturbance which is assumed to be a sequence of independent and identically distributed random variables with zero mean and variance a . Equation (3.31) is also equivalent to 2 y(t) = <p{t)0 + e(t) (3.34) Chapter 3. An Extended off-line Least-Squares Method and Time Delay Estimation 25 where ip{t) = [-y(t - 1)... - y(t - n ) u(t - l)...u(t - n )} A 6 = [ai ... a bi ... b ) . nB y(t) and the estimates £(1), (3.35) (3.36) T nA Given the measurements y(l), B e(t) of e(l), e(t) the least-squares estimate 9 of the vector 9 is determined by minimizing the following loss function N V(9) = We\t) =ee T (3.37) where e(t) is the equation error defined by e(t) = y(t) - v(t)9. (3.38) 0 = (3.39) The estimate 9 is given by {^NY^IYN if 4>%<J>N is invertible, which expresses the excitation condition (Astrom and Wittenmark[3]), where ' V(l) ' ' y(i) Y N and (J>N = = (3.40) <p(N) y(N) Some statistical properties of 9 are (Astrom and Wittenmark[3]): E{9) = f5and (3.41) cov(9) = <%{f <t> )-\ N N (3.42) Chapter 3. An Extended off-line Least-Squares Method and Time Delay Estimation 26 Recursive equations can be derived for the case when the observations are obtained sequentially as follows. Let 9(N) denote the least-squares estimate based on 7Y measurements, and assume that (<frJf(f>N) is a non-singular matrix. The least-squares estimate 9 then satisfies the recursive equations [64]: §(N + 1) = 9(N) + K [y(N + l)-<p {N + l)9(N)} T N K = P <p{N + + <p (N + l)P <p(N + 1)]" N+1 = + l)]P . N T N P (3.43) 1 N [I-K <p {N T N (3.44) (3.45) N This method is based on the assumption that e(t) is white noise. 3.2.2 Generalized least-squares The procedure of the generalized least squares considers the case when the additive noise is a colored noise i.e. the mathematical model of the process is described by A(q-i)y(t) = B( -i)„(t) - Biq-^uty g + J ^ L (3.46) + wit). (3.47) + d q- , (3.48) where A(q~ ) and B(q~ ) are as before and l l Dig- ) = 1 l + d - A... l nD iq nD and w[t) = £>(g-i) is the colored noise. w(t) can also be expressed as: w(t) = -diw{t - 1) - ... - d w{t - n ) + e(t). nD D (3.49) Multiplying both sides of equation (3.46) by D(q~ ) yields 1 A(q- )y {t) l f = B{q- )u (t) l f + e{t), (3.50) Chapter 3. An Extended off-line Least-Squares Method and Time Delay Estimation 27 where y (t) = Diq-^yit) u (t) = f f and (3.51) D{q- )u(t). (3.52) l If the autocorrelations of w(t) are known, then, define the covariance matrix W of the colored noise w(t) as: E{w(t)w(t - n )} E[w(t)w(t)] W =R W = D E[w(t - l)w(t - E[w(t)w(t - 1)] E[w(t)w(t — n )\ D n )} D (3.53) ... E[w(t — n )w{t — n )) D D The inverse of this matrix can be expressed as W~ = V V, l T where V being triangular matrix. The estimate 9 of the vector 9 is determined byminimizing the following loss function V{8) = (3.54) (Ve) (Ve). T In this case it can be shown [64] that e= {{V(i> yv<t> )-\v<i> YVY . N N N N • If the polynomial D(q~ ) is known then 9 can be determined using the least-squares 1 method based on filtered data yj and uj of equations (3.51)-(3.52) above. • On the other hand when prior knowledge of the noise is lacking, the polynomial D(q~ ) has to be estimated first. 1 Chapter 3. An Extended off-line Least-Squares Method and Time Delay Estimation 28 Clarke [15] suggested the following off line estimate D given by D = -(E E)- E H T 1 (3.55) T where e(N) = y'N)-(p(N)§ S = and [e(n + 1) e(n D (3.56) + 2) ... e(N)} D (3.57) and e(nc) ... e(l) E (3.58) e(N-l) ... e{n -N) D The filtered sequences no Y, iy( - ) (3.59) Uf(t) = u(t) + E diu(t — i) i=l (3.60) yf(t) = y(t) + d i=i t i are used for calculation of the estimate 9 = [Myf> f)<l>N{yf,Uf)] u (f) (yf,u ) N f Y (y ,u ) N f f (3.61) which can be considered as the first approximation of 6 in equation (3.56). This procedure is repeated to improve the estimation of polynomials A, B and D. It is noted that the calculation of 6 gives A and B while D is calculated separately. Therefore, the estimation of D depends on that of A and B and vice-versa. A recursive solution exists and can be found in [64]. Chapter 3. 3.2.3 An Extended off-line Least-Squares Method and Time Delay Estimation 29 Extended least-squares In the extended least-squares method the mathematical model considered is given by A(q-')y(t) = B(q-i)u(t) +^^e(t). (3.62) For simplicity we consider D(q~ ) = 1, therefore the model considered is 1 Aiq-^yit) = B{q- )u(t) + C{q- )e{t), l (3-63) l where A, B and C are polynomials in the backward shift operator q~ . The estimation of x 9 can be done as before except the changes in the arrays 0and <p which are now expressed as follows: <p(t) = [—y(t — 1) ... — y(t — TIA) 9 = [ai ... a u(t — 1) ... u(t — n ) e(t — 1) ... e(t — n )] and B bi ... b nA ng c ... t c c ]. T nc The calculation of the error estimate e(t) by means of past estimates of parameters 9 which is the same, by means of the following relation e(t) = y(t) - <p(t)9. (3.64) The recursive algorithm for this method is similar to that of a recursive least squares except the calculation of the prediction error. It is given in the following. §(N + 1) = 9(N) + K [y(n+l)-(p {N T N K P^ + l)9{N)} N = P ^ ( N + l)[l + (^ (A^+l)Piv¥'(^ + l ) ] " + 1 = T e(N) N [I-K <p (N T N + l)]P = y(N) - (p(N)9(N) N (3.65) 1 (3.66) (3.67) (3.68) The difference between these three methods resides in how the prediction error, which is an estimate of the noise, is taken into account. In the following a modified version of the extended least squares method for a regulatory closed loop system is proposed. Chapter 3. 3.3 An Extended off-line Least-Squares Method and Time Delay Estimation 30 A n A R M A model for a closed-loop system In many cases, a regulatory closed-loop system can be described by an A R M A model. Let us consider the plant model described by the following A R M A X model. Aiq-^yit) = Biq-^uW + Ciq-^eit), (3.69) where Aiq' ) = l + a q- + ... + a q- *, (3.70) B{q~ ) = b^- (3.71) Ciq' ) = l + c q- + ...-rC q- . 1 1 n 1 l + ... + b q~ 1 1 nA and nB nB 1 nc 1 nc (3.72) y(t), u(t) and e(t) are respectively the process output, the process input and the process disturbance which is assumed to be a sequence of independent and identically distributed random variables with zero mean and variance cr . e The process is assumed closed-loop stable under the feedback control law given by «(*) = -f^foW-r), (- ) 3 73 where r is a constant reference. Replacing u(t) in equation (3.69) by the expression given by equation (3.73) yields y(t) = B(q-')R(q-i) A(q-i)S(q-i) + B(q-i)R(q-i) Ciq-^Siq- )A(q-l)S(q-l) + B(q-l)R(q-l) 1 with B(q^)R{q- ) 1 A(q-l)S(q-l) + B(q-l)R(q-l) Chapter 3. An Extended off-line Least-Squares Method and Time Delay Estimation 31 This equation can be rewritten as C(q-i)S(q-i) e(t) A( -i)S(q-i) + B(q-i)R( -i) q (3.74) q (3.75) where y(t) = y(t) — y(t). This equation shows that the closed-loop system described by equations (3.69) and (3.73) can be represented by an autoregressive moving-average model ( A R M A ) defined by equation (3.75). Therefore it is then proposed to estimate the past noise sequence e(t) and use it as available information to estimate the process parameters A and B. The estimate of the past noise is obtained once the parameters of the A R M A model of equation (3.75) are determined. In the ideal case Ciq-^Siq' ) and 1 Aiq-^Siq-^ + (3.76) Biq-^Riq- ). 1 (3.77) A R M A models are widely used in the field of time series analysis because they are a general representation of a linear, stationary, stochastic process. Their parameters can be determined in a number of ways described in Box and Jenkins[13], Ljung[51] or in Astrom and Soderstrom[5]. In this last reference, the maximum likehood method is used to estimate the polynomials ct{q~ ) and P(q~ ) which are defined as follows: l l a(q~ ) = 1 + c^q- + ... + a q- Piq- ) = l+p q- + n(i l 1 1 1 1 na na and ...+P q- i>. n The cost function to be minimized is given by (3.78) Chapter 3. An Extended off-line Least-Squares Method and Time Delay Estimation 32 where the residual e(t) is a function of the observations y(l), y(2), y(t) defined by a(? )e(*) = kQ~ )y(t) _1 (3-79) l where &{q~ ) l = •3{q- ) = 1 l + + ... + a q~ and n& n& l + jhq- + ... + 1 Pn q- ». n fi d ( 9 ) and /?(g ) are assumed to be stable. Astrom and Soderstrom [5] showed that a, _1 _1 /3 which minimize V(a,(3) satisfy <% ) aiq' ) _1 1 (3.80) and al = E{e {t))>E{e\t)) = ol 2 (3.81) where E(.) denotes the expectation function, given that ria > n and a The equality o — o is obtained if n„ = n > rip. 2 2 a (3.82) and rip = rip. It is also shown that the residual e(t) is a good estimate of e(t) satisfying E(e(t)e(t -k)) = 0 for k > 1. (3.83) The estimation of e(t) can be very accurate. In the following, this accuracy is illustrated by an example. Example 3.3.1 ]52] The process considered is described by _ V { t ) _ 0.2q-'u{t) (1 - 0.89") 1 (l-0.5q-' + + 0.6q- )e(t), (1 - 0.8(ri) 2 1 ' Chapter 3. An Extended off-line Least-Squares Method and Time Delay Estimation 33 where the process disturbance e(t) is a sequence of independent and identically distributed random variables with zero mean and variance o~ = 1. The control law is given by e ^ = ^ 6 - 0 ^ 0.2 - 0.2?- w 1 ( 3 g 5 ) v " y K K 1 The estimated A R M A model for this process is given by ^ p r j where a'q- ) = 1 - 1.3955g + 0.3835g~ + 0.4585?- - 0.7784?- + 0.4394g~ - 0.0997g~ p(q-i) = l - 1.1387?- - 0.3336?- + 0.8466g" + 0.2867g- . 1 -1 2 1 3 2 4 3 5 6 4 The residual e(t) which is an estimate of e(t) is obtained by e(t) = - !)• (- ) 3 The variance obtained for e(t) is of = 1.0186 which is very close to cr = 1. 2 •y tj* 5- <b © o o e> © o 86 Its © Figure 3.3: Covariance of the residual. autocorrelation is shown in figure (3.3). It can be seen that the condition of equation (3.83) is satisfied. Figure (3.4) shows the near perfect match between the residual (dashed Chapter 3. An Extended off-line Least-Squares Method and Time Delay Estimation 34 line) and the white noise (solid line). Clearly the estimate e(t) oie(t) can be very accurate, which confirms what is already known from the references cited above. For these reasons, in the following, an off-line extended least-squares method to estimate the parameters of a process under closed-loop is proposed. The proposed technique is a batch method, and is based upon the estimate e(t) of e{t) as a first step, followed by estimates of the process parameters. 2 .O I - 2 .S 1 eoo 1 1 1 1 eoo eio eio 020 1 1 1 1 1 1 e2S Tim* 1 1 eao r 1 eas — 1 84o 1 e-as 1 eoo Figure 3.4: White noise and its estimate. 3.4 An off-line extended least-squares method As mentioned above, it is proposed to estimate first the past noise which is obtained by estimating the residual e(t), and then to use it as prior knowledge to determine the process parameters. This method is restricted to processes that can be described by the following equations: y® = f f e l T ^ ) + (3-87) Chapter 3. An Extended off-line Least-Squares Method and Time Delay Estimation 35 <t) = -|^l(y(t)-r), (3.88) where r is a constant reference (regulatory control), e(t) is a white noise and d is the process time delay. In new the following it is first considered the case of d = 0 and r = 0. B y defining a variable yd(t) = y{t) — e(t) the following set of equations is deduced from the above equations fjjrry«(*) Vd{t) = u{t) (3-89) = -j^j(yd(t) + e(t)). (3.90) These last two equations represent a closed loop system where the reference signal is the white noise —e(t) and the output signal is yd(t). Because —e(t) is a white noise we are guaranteed to satisfy the closed-loop identifiability condition for this loop. There we can approximate yd(t) as y(t) — e(t), where e(i) is an approximation of e(t) obtained above. Therefore, it appears that there is no need for an external signal to be applied to determine the transfer function which relates the input u to the output y . In this case d both processes represented by (3.87)-(3.88) and (3.89)-(3.90) have the same closed loop transfer function. The estimation of the polynomials A(q~ ) and B(q~ ) is done as will l 1 be described in the following. Equation (3.89) is also equivalent to y (t) - ip(t)d (3.91) d where ip(t) = [-y (t-l)...-y (t-n )u(t-l)...u{t-n )} (3.92) 9 = [ai ... a (3.93) d d A B h ... b ) . T nA nB Chapter 3. An Extended off-line Least-Squares Method and Time Delay Estimation 36 Given the measurements y ( l ) , y(t) and the estimates e(l), e(t) of e(l), e(t) the variable yd(t) is reconstructed from the measurements of y(t) and the estimates e(t) of e(t) as: y (t)=y(t)-e(t), (3.94) d then the extended least-squares estimate 9 of the vector 9 is determined by minimizing the following loss function i N (3.95) L t=i where e(t) is the equation error defined by e(t) = y (t) - <p(t)9. d (3.96) The estimate 9 is given by 9= (^V)- ^ 1 (3.97) if (j) (j) is invertible. This is an excitation condition see Astrom and Wittenmark[3], where T Vd(l) and (j) = Vd(N) (3.98) ^(iV) Some statistical properties of 9 are given by the following result. Proposition 3.4.1 Consider the estimate 9 defined by equation (3.97). Assume that e(t) is white noise and e is an estimate of e which satisfies equation (3.83) with zero mean. Then the following properties hold: 1. 9 is the unique minimum point of V(9). Chapter 3. An Extended off-line Least-Squares Method and Time Delay Estimation 37 2. E(§) = 9. 3. cov{6) = o-l_ {4?4>)-\ e) Proof: To prove this result, first write V{9) as v{e) = (3.99) hy -<t e] [Y -ct>e], T d ) d then setting the gradient of the loss function V{9) to zero yields dV (3.100) O = - ^ = - F / 0 + 0 (0 0), T T therefore 9 is obtained as in equation (3.97) which is also equivalent to (3.101) § = (<j> ft- <l> (<f>9 + T where 1 T " e(l) ' e= e-e) ' e(l) ' and e = (3.102) e(AQ e(N) Therefore 0 = 0 + f> 0)-V (e-e), T T (3.103) which implies that E(6) = 9, because E(e — e) = 0. Equation (3.103) shows that if e is a good estimate of e then (f/> <^) 0 (e — e) —• 0 and 9^9. T _1 T Rearranging Equation (3.103) yields 9-9 = ((j) <f))- (t> (e - e) T l T (3.104) Chapter 3. An Extended off-line Least-Squares Method and Time Delay Estimation 38 Then taking the expectation of both sides of the last equation yields E[(6 - 6)(6 - 6) ] = ol_ {<t> <t>r . T T (3.105) l e • The estimate 6 depends on the nature of the matrix (<j) <f)). T This matrix should be positive definite [4, 28] which is not always the case specially for process that is under linear feedback of sufficiently low order. For such processes there exist some linear dependencies among the columns of the matrix (b which means that the parameters can not be determined uniquely [4]. The problem with lack of identifiability due to feedback disappears if a linear feedback of sufficiently high order is used. Then the columns of the matrix (b are no longer linearly dependent [4]. This, however requires the knowledge of the true plant order, a practical impossibility. Another possibility is to consider an external input signal that is persistently exciting. In Forssell and Ljung[28] it is shown that this method leads to (<b (b) being a positive definite matrix. Furthermore, according to the T same reference[28]:"the general condition for (4> 4>) to be positive definite is that there T should not be a linear, time-invariant, and noise-free relationship between u and y. With an external reference signal this is automatically satisfied ". Therefore, for the process described by equation (3.89)-(3.90) the resulting matrix {(b (b) is positive definite. In T deed, it can be viewed as a closed loop system where the reference signal is a white noise and also because in the expression of u(t) = — g(g-i) A e(t)) there is e(t). The following section uses the method introduced here to estimate the parameters of a process with time delay. Chapter 3. An Extended off-line Least-Squares Method and Time Delay Estimation 39 Off-line time delay estimation 3.5 The process considered is described by the following transfer function where the notations and assumptions are those of the previous sections. V(t) = | | ^ 9 - w ( t ) + e(t), (3.106) d d is the process time delay for which a maximum and a minimum value are known i.e. d g [d in dmax]. A(q~ ), B(q~ ) and yd(t) are as defined in the last section. The control 1 1 m law is given in equation (3.73). Equation (3.106) is equivalent to = tp(t, d)6 y (t) d where <p(t, d) = [-yd(t - 1) ... - yd(t - n ) u(t - d) ...u(t-d A - n )} B and 9 = [ai...a bi...b, T nA Assuming the same conditions as before except for the presence of d, the loss function is V(e,d) = lJ2e (t), 2 z (3.107) i=i where e(t) is the equation errors defined by e(t) = y -<p(t,d)0. d (3.108) The estimates 9 and d are obtained from the following minimization problem [9,d\= min This problem is solved in two stages. V(9,d). (3.109) First solve equation (3.109) for a given d g [d i d ax] to obtain an estimate 9 defined by m n m d 9 = {<t> \d) j {d)y f{d)Yd 1 d l ( ) (3.110) Chapter 3. An Extended off-line Least-Squares Method and Time Delay Estimation 40 if (b (d)(b{d) is invertible, where Yd and (b(d) are defined by T <p{l,d) VdiX) and (b(d) Y = d (3.111) <p{N,d) y (N) d then the estimates 9 and d are those which minimize V(9,d) for all d e [d i d \. m n max Proposition 3.5.1 Consider the estimated process time delay d and the estimated 6 defined by equation (3.110). Assume that e(t) is white noise and s(t) is an estimate of e(t) which satisfies equation (3.83) with zero mean. Then [d, 9] is the unique minimum point of V(9,d). Proof: For a fixed d, it is shown in the previous section that 9d defined in equation (3.110) is the unique minimum point of V(6d). d can only take values between d i and m dmax, therefore d which minimizes V(9 ) exists and is included in [d d min n d \. • max In the following, the results of this chapter are illustrated by the following example. Example 3.5.1 [52] The process considered is described by (3.112) y(t) where the process disturbance e(rj) is white noise of zero mean and variance a = 1. The e control law is given by the following proportional controller: u(t) = l-y(t). (3.113) The first step consists in determining the closed loop A R M A model or the filter H(q ) : which will be used to estimate the noise sequences. The estimated A R M A model for this Chapter 3. An Extended off-line Least-Squares Method and Time Delay Estimation 41 process using A R M A X command of M A T L A B is given by H(q *) = ^ p r j where a{q- ) = 1 + 0.0843?- - 0.24939 fitq- ) = 1 + 0.1132?- - 0.2359g - 0.0504g- + 0.0208?- 1 1 1 1 -2 and -2 3 4 - 0 . 1 4 6 9 g - 0.1466g - 0.1179g . -5 -6 -7 Figure 3.5: Estimation value of the coefficient a. Using 500 different noise sequences, the estimates value of the parameters a and b are shown respectively in Figure (3.5) and Figure (3.6). The estimated values of the time delay are plotted in Figure (3.7). Figure (3.8) shows the values of the loss function V(6,d). The average estimated model is given by: (0.1842)?y(t) = (3.114) It can be seen from comparing this model to the one described by equation (3.112) that the time delay is exactly identified and the process model reasonably approximated. Chapter 3. An Extended off-line Least-Squares Method and Time Delay Estimation 42 0.2 0.1051- 0.165' 1 SO O 1 IOO 1 1 SO 1 1 1 200 2SO 300 N u m b e r o1 & l m u l a 1 i o r m 1 35a 1 400 1 -J50 1 SCO Figure 3.6: Estimation value of the coefficient b. e SB - 5JS - fc*5.4 - £ a 5 9 .a •w E -».B SU* 4.4 - 42 O 50 IOO ISO 200 2SO 3 OO 350 400 4SO N u m b * r o1 & l m u l s ' lo n& Figure 3.7: Estimation value of the time delay d. BOO Chapter 3. An Extended off-line Least-Squares Method and Time Delay Estimation 43 o.ozs — 0.0£4 - 0.023 - O SO TOO I SO 2CO ZSO 300 N u m b e r o l » l m u l a l l o rm 3SO -*00 -ISO BOO Figure 3.8: Loss function V . Applying ordinary least-squares to the example above, and by considering 500 different noise sequences, the estimates value of the parameters a and b are shown respectively in Figure (3.9) and Figure (3.10). T h estimated value of the parameter a is completly wrong. The estimated values of the time delay are plotted in Figure (3.11). The average estimated model is given by: V(t) = ( f l ^ ^ ^ + e^), (3.115) It can be seen from comparing this model to the one described by equation (3.112) that the ordinary least-squares method has difficulty to approximate the process model. The following table summarizes the obtained results from using the proposed method and ordinary least-squares. Real values Parameters Extended least-squares Ordinary least-squares a=-0.8 a 4 -0.775 3.85 IO" -0.24 5.68 10~ b=0.2 b 5 & 0.1842 0.1842 Table 3.1: Results. d=5 d oi 1.99 10~ 3.2 IO" 4 5 5 5 0 10.51 Chapter 3. A n Extended off-line Least-Squares M e t h o d and Time Delay Estimation - 0 . 2 2 \- —0.2T7 ' O " SO " IOO ' 1 SO 1 1 1 ZOO 2SO 3QO Number oi ilmuta'iorm 1 3SO « 400 I 450 I 500 Figure 3.9: Estimation value of the coefficient a by ordinary least-squares. 0.2B 0.24 \- 0.14' O 1 SO 1 IOO 1 1 SO 1 " 1 — ZOO ZSO 300 N u m b s r o l « i m u l s l l o rm 1 350 1 400 1 A SO 1 SCO Figure 3.10: Estimation value of the coefficient b by ordinary least-squares. 44 Chapter 3. An Extended off-line Least-Squares Method and Time Delay Estimation 45 350 -*00 -*SO SCO Figure 3.11: Estimation value of the delay d by ordinary least-squares. 3.6 Case of colored noise The case where the disturbance is colored is more complex to solve, as knowledge of the disturbance dynamics is required to compute yd = y — w. Consider the following process T7?^9~V*) = + (- ) 3 116 (3.117) where r = 0 is a constant reference (regulatory control), w(t) is a colored noise i.e. w(t) = H i (q~ )e(t) l no se and d is the process time delay. Equation (3.116) and (3.117) are also equivalent to y (t) d = y(t)-w(t) u{t) = = ^^q- u(t) d -j^(y ,(t)-w(t)). d (3.118) (3.119) These last two equations are similar to equation (3.89)-(3.90) except that w(t) is colored noise. To apply the previous analysis, it is needed to estimate e(t) as described above and determine the filter Hnm^q" ) then reconstruct the two signals: 1 w(t) = H (q- )e(t) l noise and (3.120) Chapter 3. An Extended off-line Least-Squares Method and Time Delay Estimation 46 y (t) d = y(t)-w(t). However, if no external signal is considered, the determination of the filter (3.121) H (q~ ) l noise can only be obtained if the process considered operates for a while in an open loop manner. Therefore, the general case of non white disturbance is not solved by the proposed method. 3.7 Conclusion The method presented here considered the case of a regulatory process where the polynomial C(q~ ) = 1. For that special case it was shown that the process is identifiable 1 without an external output signal. The general case though briefly discussed was not solved. The following chapter addresses some problems related to the use of minimum variance control as benchmark in assessing the performance of controllers. Chapter 4 Realistic Performance Assessment in the SISO Case 4.1 Introduction Harris's approach [35] for control loop performance assessment is simple and effective. The method requires routine operating data and little prior knowledge of the process. When using this method to assess control loop performance, problems may arise, some of which are addressed in this chapter. The organization of this chapter is as follows. In section 2, two existing performance index definitions are used to illustrate some of the problems related to the use of this method. In section 3, satisfactory versus unsatisfactory performance is investigated and the transition between them is established as the point at which the process starts to oscillate. In section 4, a tuning procedure to bring the output variance to an acceptable level is described and shown to be a trade-off between other criteria of interest. 4.2 Problem formulation The comparison of the output variance to the minimum output variance a^ v can be carried out by using one of two performance index definitions encountered in the literature [18, 35]. Here the performance index is defined as: (y-r) 1 < I < oo r is the setpoint. p mv 47 (4.122) Chapter 4. Realistic Performance Assessment in the SISO Case 48 When the actual performance is close to the best achievable performance, I fa 1. When p the actual performance is significantly worse than optimum, I p >> 1. To bound the performance measure between zero and one, a performance index can also be defined as V In the definition used here, only the extreme situations are well defined: good variance performance when the performance index is equal to one and poor performance when the index is large. In practice, the index takes on values between these extremes, and it is necessary to determine whether the actual controller performance is acceptable, acceptable with a small tuning correction, or so unacceptable that redesign of the controller is needed. As discussed in the introduction, the minimum variance control strategy has some limitations; therefore it is not always useful to implement. However, a good control strategy is often obtained when compromising between output variance as a performance criterion and other measures such as steady-state rejection of disturbances. To illustrate this idea, consider the following example. Example 4.2.1 The process is described by Figure 4.12 where the variance, cr , of the 2 white noise disturbance is equal to 0.36. Figure 4.13 shows I and the steady state error S p se ( defined as E[^f^} in % ) versus the proportional controller gain K. It can be seen from that figure that for small gains the output variance is relatively low but the steady state error is very high. For K > 1.5, the output variance becomes high and the steady state error is reduced significantly. A compromise between these two criteria may be needed if a proportional controller is the only available design. Thus, although I has advantages as a performance index, a better p overall design is obtained when combining it with other performance criteria. Chapter 4. Realistic Performance Assessment in the SISO Case Ref ) Controller k w l-0.67q J Figure 4.12: Process control. Figure 4.13: I and steady state error versus P-gain. p 49 Chapter 4. Realistic Performance Assessment in the SISO Case 50 Recognizing that minimum variance performance is not always desirable, the next section discusses how good operating performance can be distinguished from unsatisfactory operating performance. 4.3 Satisfactory and unsatisfactory performance From the performance index definition, poor performance is obtained when the performance index takes large values. This extreme situation corresponds to a nearly unstable loop. In reality, poor performance will be obtained long before this extreme situation occurs and, for a process with no integration, only a small range of values of I actually P represents realistic performance alternatives. In practice, the control signal inside the loop is bounded between u i m n and u max and will saturate if values outside this range are called for. Now, consider the single loop system, without integration shown in Figure 4.14. The nonlinear component (N. L.) is a saturated amplifier with a linear gain equal to 1 represented by its equivalent gain ^ y . The remainder of the process is modeled by a linear component represented by KP(jui). The saturation characteristic is shown in Figure(4.15). The equivalent gain is defined[17] as Ni(e) = ——— e = — / ire Jo f(ecos(u)t))(cos(cot)+jsm(ut))d(ujt). The stability of such a system can be investigated in the Nyquist plane, Figure 4.16, by using the Describing Function method. When no intersection between KP(ju) (curve 1) and - ^ y occurs, the loop is stable and the performance can be assessed by use of the output variance criteria. Increasing the gain K, the curve KP(ju) moves to the left to intersect - ^ y (curve 2) for the first time, causing the loop to oscillate at frequency and amplitude (cJn, e ) determined by solving P(JUJ) = Increasing the gain K causes the loop to oscillate at a constant frequency, but with increasing amplitude (curve 3). The conclusions to draw here are: Chapter 4. 51 Realistic Performance Assessment in the SISO Case e(t) N(E) lit) p> KP(jffi) Figure 4.14: A single nonlinear loop. Figure 4.15: Saturation characteristic. vlt) Chapter 4. Realistic Performance Assessment in the SISO Case Figure 4.17: Partition of Ip curve. Chapter 4. Realistic Performance Assessment in the SISO Case 53 1. in this case there are two situations which correspond to either an oscillating or a non-oscillating loop, 2. the performance may or may not be considered acceptable up to the point I o, at p which the loop starts to oscillate, 3. performance is definitely poor from Ipo up to IpUmit which is obtained when replacing a saturated amplifier by a relayor k —• oo. This analysis is not restricted to a saturated amplifier but can be extended to other types of nonlinearity. These conclusions are summarized by Figure 4.17. This approach provides a link between oscillation and performance. It recognizes that whenever oscillation is present in a process, the output variance represents poor performance. It is also noted that the beginning of unsatisfactory performance is obtained at the latest for I = I o which is finite << co and corresponds to an oscillatory system. p p The loops causing these oscillations should be localized in order to remedy the problem. In the following chapter the problem of localizing the oscillating loops will be addressed. IpO plays a major role in some tuning methods because it defines a limit between acceptable performance and unacceptable performance. The Ziegler-Nichols method, for example, uses (indirectly) Ipo to determine the controller tuning parameters as functions of the gain k and the period P . u u The ultimate gain k and period P are obtained u u by setting integral and derivative control action to zero, and increasing the proportional gain k until the measured variable oscillates at period P for a certain gain k (Shinskey c u u [59]). Figure 4.18 from Fu and Dumont [30], shows an oscillating loop. Tuning was performed at t = 27 minute following which the oscillation disappeared and the performance index dropped to a lower value. Chapter 4. Realistic Performance Assessment in the SISO Case 54 Mill Trial Result: A Dryer Pressure output Controller Tuning was performed based o n the evaluation performance Figure 4.18: Performance analysis. 4.4 Control chart The curve of figure (4.17) can be used for supervisory purposes. Given IpQ, determined during the tuning procedure and the estimate cr^„ of the minimum achievable output variance, it is proposed to rescale the calculated performance index as given in the following. I (%) P ^-^100. I pO (4.123) This new performance index is expressed in percentage (% ). This establishes a new scale for I , which can be interpreted as how far the actual controller performance is from that p of a minimum variance controller and how close it is to the worst performance. The alarm should be given as soon as the performance reaches the range 85 ~ 100. The performance index within the non-oscillatory region (figure 4.17) is satisfactory only if it results from a useful trade-off between all criteria considered. In the next section, a tuning procedure is presented. \ Chapter 4. 4.5 Realistic Performance Assessment in the SISO Case 55 Tuning procedure Before compromising between the criteria it is necessary to estimate the useful portion of I curve, corresponding to I G [1 Ipo], in terms of the controller parameters. p p Note that I o can be estimated by the relay method. p Although I itself may be estimated without knowledge of the process (except any p delay), the variation of I with controller parameters can only be predicted once process p dynamics have been estimated. To estimate the I curve, a design method is proposed in which: p 1. The process is approximated, Figure 4.19, by one of the following models o-T s d or M ( s ) = M ^ s ) = f ^ , M (s) = ~ s + as + b s(s + as + b)' 1 + rs' 9 6 2 T d S 3 2 2 ge •TdS l+TS ge •T s d s + as +b 2 ge- TdS s(s ^as+b) Noise Ref C ontroller Process Figure 4.19: Process approximation. This choice covers the most common process behaviours encountered in process industries [44]. The model chosen is that minimizing the error between its output and that of the process, i.e., minimizes the cost function: J{Mi) 1 N il2 Chapter 4. Realistic Performance 56 Assessment in the SISO Case The models can be obtained using the method of chapter 3. They can also be obtained by estimating the time delay and dynamics separately: • the time delay Td is estimated by one of the two techniques described in Elnaggar and Dumont [22] or in Zheng et al.[71]. • g, T, a and b are estimated for each model using a least squares method assuming sufficient excitation. 2. a noise model is then estimated by one of the two following techniques: Time series analysis[35] or Laguerre network[52]. 3. the controller is adjusted in order to obtain the I curve and corresponding curves p for other performance criteria of interest. e(t) l-OV l-O.CTcf 1 Figure 4.20: Process control. In considering these curves one can: • establish the compromise controller tuning parameters, • be guided by these curves to define "satisfactory" performance, and • decide the direction in which the controller should be adjusted. Chapter 4. Realistic Performance Assessment in the SISO Case 57 The controller which results from this procedure is, therefore, a trade-off between all criteria considered. The following section illustrates this method by means of an example. 4.6 Simulation Example The process considered is shown in Figure 4 . 2 1 . A proportional controller is taken in the first test, and a proportional plus integral (P.I) in the second test. The estimation of the I curve for the interval [1 7 ] is carried out following the approach described in section p p o 4. The following model provides the closest representation: 0.232g- y(t) = y\ > , , 4 0 . 6 2 3 4 ? " u(t) 1 - 1 1 ncoo^ H T -i 0.5339- w -r , . 1 — e(t). 7 1 _ 0.62349- 1 w 0.33(1-0.6^^-* Ref Controller (1-0.6q-)(\-0.2ci) A - l l Figure 4 . 2 1 : Process control. F i r s t test: Cbntroller=P-gain i.e. proportional only. The two criteria, output variance and steady state error S er are considered (S is defined as E[ ' ^] r ef er in %). Figure 4 . 2 2 show the simulation results. The solid curves represent the true process, and dashed curves are obtained using the method described above. It can be seen that the estimation of I p curve and steady states error on the basis of the simplified process model is Chapter 4. Realistic Performance Assessment in the SISO Case 58 very close to true values, especially at low gains. From these curves, the best choice of tuning parameter can be determined as a compromise between steady state response and disturbance rejection. Figure 4.22: Comparison of true and approximate model. Second test: Controller^ P i.e. Proportional plus Integral control. Interest is in the process step response, and in output variance. I , tr (tr is defined as the time needed p for output y(t) to reach 90% of y(oo)) versus controller gains for the actual process, and its approximation is shown in Figure 4.23 and Figure 4.24. From these curves it can be seen that: • The estimated curves are close to the true ones and in this case lead to satisfactory controller design. • Combining the curves Figure 4.23-4.24 one can establish the best compromise controller tuning parameters. These curves also indicate the direction in which the Chapter 4. Realistic Performance Assessment in the SISO Case 59 controller parameters should be changed to improve a certain performance criterion. • Similar curves may be obtained for any controller candidate with a small number of adjustable parameters. Figure 4.23: Estimated curves versus Controller parameters . Remark: This technique uses the minimum variance o mv = (1 + ft + ... + fx -i)4 d to normalize curves of the output variance. Note that although incorrect delay estimation may lead to one or more terms, such as f^-i^l of o , mv being added or missing in the estimation the effect is to move the I curve up or down. This estimation error, thus, need not p invalidate the method since the choice of the tuning parameters is always a compromise among all the criteria considered, and the actual value of Ipo is usually less important than its variability. Chapter 4. 60 Realistic Performance Assessment in the SISO Case O O .3 1 1 -S 2 O O .3 1 1 .3 2 Figure 4.24: Estimated curves versus Controller parameters . 4.7 Conclusions In this chapter it was linked oscillations to poor performance and shown that by broadening Harris's performance index as a measure of output variance to consider other criteria a better overall assessment of control loop performance can be obtained. A n approach to controller tuning based on a combination of performance criteria and output variance was also proposed. Because oscillations amputate the control loop performance, oscillating loops should be localized and controllers should be tuned. In the following chapter the problem of localizing the oscillating loops is investigated. Chapter 5 Monitoring Oscillations in a Multiloop System 5.1 Introduction Poor performance of process control in the pulp and paper industry may be attributed to three causes: 1. The first one is due to nonlinear behavior of control valve and sensors [9]. For example at the Mead Fine Paper mill in Chillicothe, O H , about one loop in five cycles at steady state[23] because of nonlinearity in control valves. 2. The second cause is simply use of the control strategy or poor process design. 3. The third cause results from the fact that most mills have a great number of control loops and a very limited staff of instrument technicians and maintain to check the adequate functioning of each loop. A cycling loop affects the performance of others as it spreads the oscillations to its surrounding loops; therefore, it is critical for overall loop performance to detect and eliminate oscillations. To illustrate this fact, consider a loop where the process is described by (1 - 0.74 - )y(t) = 0.26q- u(t) 1 A 9 + (1 - O.bq' + O.Qq~ )e(t) 1 2 and the control law is given by 0.5 - 0.37?- 1 u(t) = w 7(r(t) - 0.26-0.13?- - 0 . 1 3 g 1 61 4 V y ) y(t)), y y > h Chapter 5. Monitoring Oscillations in a Multiloop System 62 where r = 5 is the setpoint or reference signal, y(t) is the process output variable and e(t) is a white noise sequence (^(t) = 0-81) • Assume that this loop is oscillating at frequency OJ = 0.314 rad/sec as result of its interaction with an oscillating loop. This situation 0 is simulated by adding a sine wave disturbance signal equal to 2.5sm(0.6287r£) at the output of the loop. The output signal y(t) can be expressed as: y(t) = y + % ) e ( £ ) + A sin{u t _ 1 0 0 (5.124) + <f> ) 0 where ft(<? ) is the transfer function which relates e{t) to y(t). A n estimated A = 3.207 -1 a is obtained using equation (5.127)-(5.130). The minimum achievable output variance is °~mv = Ip{v{t)) 1-145. The actual output variance is o ^ = 7.84. The performance index is 2 = 6.877. Based on this value of I , the performance of the actual controller p may be judged poor. However, this loop is not generating the oscillation. Then, by subtracting the contribution of the frequency u> to the output variance, which may be 0 A 2 estimated as -f = 5.1424, the estimated output variance without the frequency u would 0 be A 2 _ 2 °~(y(t)-oscillation) — y(t) a O _— ^~ r, n Z D -- Therefore, a more realistic performance index for this loop would be I i = 2.3, which is rea significantly less than the one calculated before. This example shows how the variability of a non-oscillating loop is increased when interacting with other oscillating loops. It also shows that not taking into account oscillations may yield a wrong performance index calculation. Therefore, it is of interest to localize oscillating loops and, using an appropriate controller tuning method, to reduce the amplitudes or eliminate the oscillations completely. Knowledge of the presence of oscillations in loops may be obtained by techniques such as Discrete Fast Fourier Transform ( D F F T ) or Singular Value Decomposition (SVD)[47]. Chapter 5. Monitoring Oscillations in a Multiloop System 63 The issue of localizing the oscillating loops in a cascade multiloop system is addressed in this chapter. The organization of this chapter is as follows. In section 2, background on the concept of equivalent gain for nonlinear systems is given. In section 3, an oscillation index to help localize the oscillating loops, in a cascade multi-loop system, is define. In section 4, the accuracy of using the describing function method is discussed, followed by a simulation example, and then a conclusion. Note that the analysis in this chapter is performed in the continuous time domain. 5.2 Background A n industrial process consists of a combination of linear and nonlinear components. When it is possible to separate the process into a nonlinear gain and linear dynamics, then such process can be represented by Figure (5.1) where e(t) is a white or colored noise. The behavior of such systems has been extensively discussed in the literature Chapter 5. Monitoring Oscillations in a Multiloop 64 System [17, 10]. This investigation discusses, for any loop with (or without) nonlinear components, conditions under which oscillations are present and, if so, the form of their characteristics. One method commonly used to address these questions, for systems of the form of Figure 5.1, is the describing function (DF)[17], introduced during the late 1940s. The D F is an approximation procedure that recognizes, for the simple feedback loop of Figure 5.1, that when limit cycles occur the input to the nonlinearity N is almost sinusoidal due to the low pass filtering action of a typical control system plant, G(s). The nonlinearity can therefore be approximated by an equivalent gain or describing function, N(e), defined as the complex ratio of the fundamental of the output of the nonlinearity to the amplitude of a sinusoidal input at that frequency [10]. Assuming that the nonlinearity input is e(t) = ecos(cui), its output is thus approximated by u(t) ~ uiCOs(tvt) + U2sin(ut). If u = f(e) represents the nonlinearity input/output relation, the equivalent gain is given by Ni(e) = e = — / ns Jo f (e cos(ut))(cos(ut) + jsm(ut))d(ut). For the autonomous case (r(t) = 0), the prediction of the limit cycle rests upon finding a solution to the equation l + N(e)G(s)\ „ J = 0, (5.125) which is equivalent to a unit loop gain for the fundamental component of the signals e(t) and y(t) denoted respectively Ef(t) = efsm(ut Vfit) = yfsm(wt + <l)y ) £f f where (f> — <f) = —180°, resulting in yf + <fi ) ef = -Vf(t)- Chapter 5. 5.3 Monitoring Oscillations in a Multiloop System 65 Monitoring the oscillations The reverse problem, that of limit cycle prediction, involves recognizing which loop in a multiloop system is causing oscillation at a particular frequency. For that, an oscillation index is defined as follows: Definition 5.3.1 For the linear feedback loop of Figure 5.2 where an oscillation of frequency LU is detected, define the oscillation index, at frequency LU, by / 0 S C H = |1-||, (5.126) where ^ is the gain of the considered loop at frequency LU. This gain is determined by extracting the magnitude and the phase of the oscillation LU from both signals y(t) and e(t) using the following set of equations: y = y\ = y/yi+y*, (5.127) j,J (5.128) y{t) cos (ut)dt, 2 2/2 <f> y [T = j; J y(t) sin (Lut)dt, = arctan(^). Q and (5.129) (5.130) e is obtained from the following set of equations: e = \Jel + sl (5.131) £ l = 2 fT — f (t) cos (ut)dt, (5.132) e = 2 r — / e(t)sm(ut)dt, T Jo £ T 2 <f> £ T Jo = arctan(—). and (5.133) (5.134) £2 For a cascade multiloop system, where loops are defined by sets of a closed paths according to Mason's loop rule, the following result is proposed: Chapter 5. Monitoring Oscillations in a Multiloop 66 System e(tj r(t) rL[5)t<\ f r M IJ\5 } Figure 5.2: Process control. P r o p o s i t i o n 5.3.1 For oscillation at frequency to at the output of loop U in a multiloop system, the loops associated with the oscillation u are those with /^ (cu) ~ 0. c Proof: To prove this result, first consider the single loop of Figure 5.2. This loop is oscillating at frequency u if ju> is a solution to the following characteristic equation 1 + C( )G( )U = 0. 5 (5.135) S This last equation is in fact a set of two equations which are =ju) s lC(s)G(s)\ , s=JU — 1 gain equation, (5.136) = - 1 8 0 ° phase equation. (5.137) The frequency UJ was detected at output of that loop; therefore, it remains to extract the magnitude and the phase of that frequency from both signals e(t) and y(t) to find out whether equations (5.136)-(5.137) are met. In fact this is what is done by equation (5.126) and equations (5.127) to (5.134). Chapter 5. Monitoring Oscillations in a Multiloop System 67 The equation (5.135) is a necessary and sufficient condition for any loop to generate oscillations at frequency UJ; therefore, for the loop k taking from the considered multiloop system, to generate oscillations at frequency UJ equation (5.135) should hold which implies that Iii (UJ) osc ~ 0. • The extension of these results to a nonlinear multiloop system is proposed by the following result. Proposition 5.3.2 For oscillation at frequency UJ at the output of a nonlinear cascade multiloop system, the loops associated with the oscillation UJ are those with /o (cu) — 0. l Proof: sc To prove this result, first consider the oscillating nonlinear feedback loop of Figure 5.1, where the input r(t) = 0, and then consider the case of r(t) = RQ+RI sin(a;it). The extension of these results to a multiloop system will be discussed in the following. The autonomous case, r(t)=0 From the above, any limit cycle satisfies Equation 5.125. To obtain the equiva- (UJ,O) lent gain for the nonlinear component, the D F approach of [10] is followed. Write the approximation at the output of the nonlinear component as U(t) Ui cos (cut) + « 2 u sin (ut + (b ) u sin (ujt) (5.138) (5.139) where (5.140) 1 Jo arctan (—) u 2 (5.141) (5.142) Chapter 5. Monitoring Oscillations in a Multiloop System 68 then, u is given by u= y/ul + v%. (5.143) For the nonlinearity input a similar approximation gives e(t) ~ £i cos(ut) + £ sm(cot) (5.144) = esm(ut +(j) ) (5.145) = —t (t) cos(ut)dt (5.146) = | ; f e(t)sm{tot)dt (5.147) arctg{—) (5.148) 2 £ where £ l e 2 (j> = e T Jo £ T 1 Jo £2 and £ is given by £ = yjel + el (5.149) The estimated equivalent gain N for (cu, a) is given by \N\= ^^1^, Z0 -0eu (5.150) For the linear component, because it satisfies the superposition theorem, the gain at co can be obtained by applying the same reasoning to the fundamental component: y(t) = yi cos(ut) + y sm(ut) + terms (5.151) = ysm(cot + (j) ) + terms, (5.152) 2 y where y\ and y are obtained by the same technique as above, y is given by 2 y= \fyJ+?2- (5.153) Chapter 5. Monitoring Oscillations in a Multiloop System 69 The estimated linear gain at u> is given by \G(u)\ = ^ A L^-K. (5-154) It should be noted that the integrations considered above are affected by the presence of noise, n(t), introducing the following integral terms / n(t)cos(wt)dt, Jo [ n(t) sm(wt)dt Jo (5.155) (5.156) in the expressions above. To minimize the contribution of the band limited noise terms, the integral in equation (5.155)-(5.156) should be taken over a number of periods, kT, instead of one period. Replacing G and ./V by their estimated value in equation (5.125) yields, for this simple case: s I (u) osc = 1 (5.157) = | 1 - | | = 0. (5.158) Equation (5.125) is a set of two equations: the magnitude equation which is used for equation (5.157) and the phase equation: </> -<f> = - 1 8 0 ° . y e (5.159) Equation (5.159) must be satisfied for a uniformly maintained oscillation OJ, both the total gain must be equal to one, and there must be a total of 180° of dynamic phase lag in the loop. R e s u l t 1: If the oscillation to is generated inside the loop of Figure 5.1, then Io8c{u) = 0. Chapter 5. Monitoring Oscillations in a Multiloop System 70 The above analysis considered the case of an autonomous system. To generalize this result, consider the cases of an external input, and of an external input, combined with auto-oscillation. Forced system The input into the system of Figure 5.1 is now assumed: r(t) = R + Ri sm(wt). (5.160) 0 The most obvious possibility is that all signals in the feedback loop can be approximated by (5.161) e(t) ~ eo + eis'm(ut + (j) ), £ and y(t) by y(t) ~ G(0)N e 0 + e i | G f » | W ( e i ) sm(wt + </>'). 0 Since e{t) = r(t)-y(t), it follows that e(t) (5.162) ~ £ + £isin(ut ~ R^ + R^m^-G^Noeo-exlGijLJ^Nie^smicot+ ~ Ro - e G(0)No + yjR\ + [e^G^N{e^ 0 + <j> ) = r(t) - y(t) e (/)') (5.163) sm(ut + 0,) 0 (5.164) and by equating coefficients results R 0 = e (l + G(0)N ) o o and R\ = (l-[|G0u;)|Ar( )] )4 2 £l (5.165) (5.166) Chapter 5. Since R\ Monitoring ^ 0 implies Oscillations in a Multiloop 71 System ^ 1, the following result is established. [\G(JUJ)\N(EI)] 2 Result 2: If the oscillation UJ is not generated inside the loop then The estimated total gain is obtained from data by the procedure given \G(JUJ)N(EI)\ above. Now, consider the case of an external input r(t) at frequency UJ into the loop of Figure 5.1 combined with internal oscillation at frequency UJ . Then the input into the system C of Figure 5.1 is r(t) = R + (5.167) R sm(ujt). 1 0 The most obvious possibility is that all signals in the feedback loop can be approximated by e(t) ~ €Q + e\sin(ujt + <j>i) + E sin(uj t + (b ) c c c and y(t) by y(t) = G{0)N e +\G{juj )\N{e )e sm{ujt 0 o 1 1 +\G(juj )\N(e )e sm(ujt c c + f> ) 1 ( 1 + <f> ). c c Rearranging this equation and equating coefficients as before results in: Ro = R\ = G(0)N ) (5.168) ( l - O G W I ^ e O H e ? (5.169) e (l o + o (l-[\G(juj )\N(e )] )el 0 = 2 c c Repeating the same analysis leads to the following conclusion: I c(u) ^ 0 and Iosc(u ) = 0. OS c (5.170) Chapter 5. Monitoring Oscillations in a Multiloop 72 System Thus the evaluation of the oscillation index for each detected oscillation frequency at the output of a loop will distinguish between oscillation generated inside the loop and that coming from outside into the loop. Now return to the multiloop system, the condition for limit cycles at frequency UJ and amplitude a to occur in a loop (i) is that the fundamental component of the signals y (t) l and e (t) l Ef(t) = efsm(ut + <f) ) ef y){t) = y}sm(ujt + <j) ) yf where <p yf — (b £f verify E){t) = -y}(t). (5.171) The reasoning in the above analysis completes the proof of the proposition (5.3.2). • The approximation considered here can be applied still to a signal containing more than one frequency [17]. Note that when a loop is externally excited, auto-oscillation may be hidden. This behavior does not cause trouble once the remaining detected frequencies have been assigned to loops. B y tuning oscillating loops in sequence, the hidden oscillations will appear and the tuning procedure may be repeated until all loops are tuned. In the following section the effectiveness of the describing function method will be discussed. 5.4 Accuracy of the Describing Function Approach The number of terms that should be included in the oscillatory signal component, xap(t) can be determined using the normalized Mean Square (MS) error given by E{\x(t) - x (t)\ } E{\x(tW} • 2 2 6 ap ( 5 - 1 7 2 ) Chapter 5. Monitoring Oscillations in a Multiloop 73 System Taking e to have a value less than 0.05 implies that the approximation accounts for 95% of normalized MS variation of x(t). 5.5 Some Limitations The proposed method is limited to cascade loops. Note that a given loop may oscillate at frequency LO either if ju> is solution to the equation (5.135) or as result of a periodic disturbances at frequency LO acting on that loop. If periodic disturbances may be present these two cases can not be distinguished directly with the present method. 5.6 Simulation example disturbances i»,(t) mi) NI.2 N.LI sfl Figure 5.3: Process control. To illustrated the proposed method of oscillating loop detection it is considered in the following two examples. Example 5.6.1 Consider the cascade control system of Figure 5.3 where the disturbances are white noise (aj isturbances G(s) = 0.4), = 1 and s + 2s + 4s 3 2 (5.173) Chapter 5. Monitoring Oscillations in a Multiloop System 74 (5.174) The nonlinearity N . L . 1 is described by if ini(t) > 5 5 outi(t) k ini(t) if - 5 < ini(t) < 5 p -5 if ini(t) < - 5 while the nonlinearity N . L . 2 is described by 5 out (t) 2 ifin (t)>5 2 k in (t) if - 5 < in (t) < 5 2 2 -5 ifin (t)<-5 2 The process considered here, is formed from two loops. Therefore, if an oscillation is detected at the output y(t), and assuming that it is not resulting from an external periodic disturbance, which is not handled in this study, then this could not happen unless at least one of the two loops is oscillating. In such a case, and provided the scheme of the two loops, the oscillation index of the outer loop at the detected oscillation frequency will be equal to one and that of the inner loop at the same frequency will be equal to one only if it is causing the oscillation. Therefore, only the oscillation index for the inner loop at the detected frequency will be calculated in the following two test. Test 1: First, set linear gain k = 4 and k = 10. This corresponds to an oscillation p frequency at to = 0.589 rad/sec and a period of T = 10.667 sec. The oscillation is caused by the outer loop. The output of the outer loop is shown in Figure 5.4. The output of the inner loop is shown in Figure 5.5. The output of the nonlinearities N . L . 1 and N . L . 2 are shown in figure(5.6) and figure(5.7) respectively. From these figures it can be seen that the oscillation drives the system into saturation regions for N . L . 1 and N . L . 2. The estimated total gain from the input i n to the output y at oscillation to = 0.589 rad/sec, 2 2 Chapter 5. Monitoring Oscillations in a Multiloop System 75 is total gain (0.589) = 0.67. (5.175) Figure 5.4: The output y(t) of the outer loop, example(5.6.1)-test 1. Therefore the oscillation index for the inner loop is 7 (0.589) = |1 - 0.67] = 0.33 f 0. OSC (5.176) It can be concluded that the inner loop is not generating the oscillation at frequency to = 0.589 rad/sec. It is generated by the outer loop. Test 2: Now consider k = 10 and k = 0.5. This corresponds to an oscillation frequency p at LO = 2.47 rad/sec period of T = 2.53 sec. The oscillation is caused by the inner loop. The output of the outer loop is shown in Figure 5.8. The output of the inner loop is shown in Figure 5.9. The output of the nonlinearities N.L. 1 and N . L . 2 are shown in figure(5.10) and figure(5.11) respectively. From these figures it can be seen that the Chapter 5. Monitoring Oscillations in a Multiloop System 76 Figure 5.5: The output yi(t) of the inner loop, example(5.6.1)-test 1. " nnnnn ~ • .n n" i i —— i — "n r ~1 - JUL tOO 150 200 2SO Tlrm(Mc) 300 350 Figure 5.6: The output outi(t), example(5.6.1)-test 1. <00 Chapter 5. Monitoring - - Oscillations in a Multiloop —r - - - - System - 77 —r- — , r - - -1 - - - _ 61 IOO 1 _ _ 1 5 0 - - _ 1 2 0 0 - - 1 2 S O Tlrrm(S9c) i_1 - J 3 0 0 - - - - I I 3 5 0 I O O Figure 5.7: The output out (t), example(5.6.1)-test 1. 2 oscillation drives the system into saturation regions for N . L . 2. The estimated total gain from the input m to the output y 2 2A7rad/sec, 2 at frequency u = is The oscillation index for the inner loop is I (2.47) = |1 - 1.0097| = 0.0003 ~ 0, osc (5.177) therefore, the oscillation is caused by the inner loop. E x a m p l e 5.6.2 The process considered is defined in the example (5.6.1) where the nonlinearity N . L . I is as defined above where k p = 1 and the nonlinearity for N . L . 2 is described by the following equation. out (t) = sign(in (t))[|in (t)| + 3]. 2 2 2 (5.178) The discontinuity offset of 3 at zero models coulomb friction and the linear gain of one models viscous friction. Chapter 5. Monitoring Oscillations in a Multiloop System Figure 5.9: The output yi(t) of the inner loop, example(5.6.1)-test 2. 78 Chapter 5. Monitoring 0.2 0 Oscillations in a Multiloop System 79 ' 100 " 1 SO ' 200 ' 250 ' 300 1 350 Figure 5.10: The output outi(t), example(5.6.1)-test 2. -400 Chapter 5. Monitoring Oscillations in a Multiloop System 80 The simulation results of this new process showed an oscillation frequency at co = 2.25 rad/sec and a period of T = 2.78 sec. The output of the outer loop is shown in Figure 5.12. The output of the inner loop is shown in Figure 5.13. The output of N . L . 1 and N . L , 2 are shown respectively in figure (5.14) and figure (5.15). Figure (5.16) shows in2(t) versus out2(t) which indicates that the inner loop is driven to the non linear region of its non linear component. The estimated total gain from the input i n to the output 2 3 - 2.5 - _, I , IOO . ISO 200 . 250 Tlm0(MC) , 300 , 3SO 1 AOO Figure 5.12: The output y(t) of the outer loop, example(5.6.2). y2 at frequency LO = 2A7rad/sec, is total gain = 1.039. (5.179) The oscillation index for the inner loop is 7 (2.25) = |1 - 1.039| = 0.0039 ~ 0, osc therefore, the oscillation is caused by the inner loop. (5.180) Chapter 5. Monitoring Oscillations in a Multiloop System Figure 5.13: The output y\(t) of the inner loop, example(5.6.2). Figure 5.14: The output outi(t), example(5.6.2). 81 Chapter 5. Monitoring Oscillations in a Multiloop System Figure 5.16: The output out (t) versus in-2(t), example(5.6.2). 2 82 Chapter 5. Monitoring Oscillations in a Multiloop System 83 What can be learned from this study, is that oscillations are transmitted between interacting loops, cause excess variability, and may yield wrong performance index calculation. 5.7 Conclusion In this chapter, it was shown that oscillations transmitted between interacting loops, cause excess variability, and may yield wrong performance index calculation. Therefore, for cascade loops a unified oscillation index that can be used for testing linear and nonlinear loops to help localize any loops oscillating at a given frequency is introduced. The procedure does not require knowledge of the transfer functions; however, it requires knowledge of the frequencies involved. Knowledge of the limit cycle conditions can be used to tune the controller (for example using Ziegler and Nichols method [59]) and so eliminate oscillation and guarantee better performance. Up to now only SISO systems were considered; however, many processes are M I M O systems. The assessment of process control of such multivariable systems will be presented in the next chapter. Chapter 6 Control Loop Performance Monitoring in the M I M O case This chapter is concerned with the extension of Harris's control loop performance assessment technique for SISO systems to Multi-Input Multi-Output (MIMO) systems. This extension was considered in Harris et al.[36] and in Huang et al.[40]. The approaches developed in these two references are quite similar as both require the knowledge of the interactor matrix to perform the assessment. This knowledge may be considered the main drawback for these two techniques as the interactor is not unique and requires significant knowledge of the system for its determination. The approach presented in this chapter is based on defining an absolute lower bound on the achievable value for each output variance, thus avoiding the need to identify the interactor matrix. The organization of this chapter is as follows. In the first section, strategies of control loop performance assessment based on multi-input multi-output (MIMO) minimum variance control (MVC) are discussed. This is followed by the definition of ultimate output variance bounds. These bounds are used in section 3 to define a performance index associated with each process output. The accuracy of the proposed method is compared to the existing ones and is investigated by means of an example. In section 4, a procedure for estimating the performance index is given, followed by a simulation example, and, finally, some conclusions are drawn. 84 Chapter 6.1 6. Control Loop Performance Monitoring in the MIMO 85 case Assessment strategies based on M I M O minimum variance control ( M V C ) In this section M I M O minimum variance control strategies are discussed together with methods of assessing performance based on M I M O M V C . Consider a M I M O linear discrete-time system described by the following transfer function y(t) = T{q- )u(t) l + Niq-^eit) (6.181) where y(t), u(t) are m x l output and r x 1 input vectors, respectively. e(t) denotes the m x 1 innovations white noise sequence. T(q~ ) and N(q~ ) are rational matrices in the 1 1 unit delay operator q" . 1 In recent years, many algorithms have been described for the adaptive control of linear stochastic systems. For more discussion of these algorithms see, for instance, Goodwin et al.[31] (1984). Early work on (MVC) by Astr6m[l] focused on Single-Input/Single-Output (SISO) systems. Borisson [12] (1979) considered the Multi-Input/Multi-Output (MIMO) M V C where a common delay was associated with all system outputs. In Goodwin et al.[32], a different delay was associated with each output. A more complicated M I M O time delay representation was proposed by Elliott et al. [21] (1982), who showed that the notion of delay in the M I M O case is related to the system interactor matrix £(g) introduced by Wolovich et al.[70] (1976) which satisfies the following two properties: 1. l i m ^ o o £(q)T(q~ ) l 2- = K for any non singular matrix. = #(9)diag(<7 , ...,q ) where H(q) is a unimodular matrix. dl dn Chapter 6. Control Loop Performance Monitoring in the MIMO case 86 where H(q) is a unimodular matrix of the form 1 0 0 h21 1 0 H(z) 0 0 — h m—l 1 m and hij(q) is either divisible by q or is zero. B y considering the Markov parameters representation and by explaining the meaning of property (2), Shah et al. [58] gave an interesting interpretation of the interactor matrix as follows: TV ) 1 = Y.Giq •t-i i=0 = foq* + Ziq*- + ... + U-iq 1 1 where d is the order of f (q) lim T(q)t(q) = £ G - i + ^G 0 d d 2 + ••• + td-iG 0 = K. K represents, then, the linear combination of the rows of the first d Markov parameters. That means that the interactor matrix satisfies the following algebraic equations: CoGo = + ioGd-i + £,iGd-2 + ... + £d-iG 0 = 0 = 0 Solving these equations gives another way of evaluating the interactor matrix. The first method to evaluate £(q) is the one given originally by Wolovich and Falb [70] and used in [31]. Chapter 6. Control Loop Performance Monitoring in the MIMO case 87 W i t h the above interactor, which is in lower triangular form, the first output variable is decoupled from the other outputs. The spectral unitary interactor matrix introduced by Bittanti et al.[ll] (1994) is an interactor matrix preserving the spectral properties of the underlying system. It satisfies t(q)Z{q- ) 1 T = i, (6-182) and also satisfies property (1) above. Regardless of the interactor matrix adopted, it can be shown that the minimum achievable output variance, for a non minimum phase system, is a finite moving average polynomial of the form y(t)-y*(t)=C\q)F(q)e(t), (6.183) where {y*(t)} is a given reference trajectory and £(q) is the time delay representation considered. F(q) and R(q~ ) are obtained from the following equation 1 Z(q)N(q- )e(t) = F(q)e(t) + R(q- )e(t), 1 where F(q) = F q d F(q)e(t) d 1 (6.184) + ... + F q and i ? ( g ) is a proper transfer function. The terms _1 x and R(q~ )e(t) are respectively, future and past noise sequences. 1 Note that the expression for this tracking error (equation 6.183) does not depend on controller parameters. As for SISO assessment, one might be tempted to use the predicted performance of a M I M O minimum variance controller as an absolute lower bound against which current performance can be assessed. A bound is obtained when taking the variance of both sides of equation (6.183). However, this lower bound will depend on the type of interactor used. In the next section a new procedure for performance assessment without knowledge of the interactor matrix is proposed. Chapter 6. 6.2 Control Loop Performance Monitoring in the MIMO 88 case Ultimate output variance bound In the last section it was shown that the minimum achievable output variance depends on the type of time delay matrix used. This non-uniqueness, together with the fact that the structure of an "optimal" interactor matrix is not yet known, leads to the following technique in which an absolute bound is independently established for each output. Proposition 6.2.1 Consider a minimum phase system. The minimum achievable output variance for output y.i(t) is given by: K)L = EWCFjEf, (6.185) k=l where F is a vector containing the coefficients of the polynomial F (q) obtained from l l the following Diophantine equation q where d l N , the i m th Proof: i n ) N \ F l and R l d L n N ^ q - 1 ) = F \ q ) + R ^ q ' 1 (6.186) ) , are respectively, the minimum delay in row i, the i th row of row of F and that of R. First consider the output ordering y(t) = {yi(t), y2(t), ...,y (t)) m and consider the lower triangular interactor matrix as a M I M O time delay representation. Dugard et al.[19] (1984) showed that an absolute bound or an unconditional minimum output variance is achieved for the first component of the process output, i.e. for yi(t). For the other components, the optimization is achieved subject to constraints. To determine the lower bound achieved by output (yi(t) — y*(t)), return to the major steps in the derivation of the M V C strategy. • First note that the first row of the interactor matrix is 6(9) = ..,o), where e ^ is the minimum delay in the first row of in T ( q ~ 1 ) . Chapter 6. Control Loop Performance Monitoring in the MIMO case 89 • Then, the first row of equation (6.184) is ^ where F ( q ) 1 = q F \ + L ... n + j q V d ^ 1 - 1 ) , m i n F i l F = ^ q ) + R and N (q), x ^ R q ' ) 1 ^ 1 , 1 ) are respectively the first min row of F ( q ) , the first row of and that of N(q~ ) l R ( q ~ 1 ) . • Under the M V control and in steady state the tracking error is a finite moving average polynomial given by equation (6.183). Therefore, yi{t)-y{{t) = q - d l ^ F \ q ) e ( t ) , A n d the minimum achievable output variance (<r^) „ for output yi(t) m is d . 1 - E{(vi(t) 2 F}Q(F}) T (6-187) yi(t), consider the following output ordering Now, for an arbitrary output yif) yKt)) } = E = (yt(t), y2(t),...,yi(t),.,y (t)), m then replace the subscript 1 by the subscript i in the above equations to complete the proof. • It is shown in Dugard et al. variance ( C y . ) m v (1984) [19] that the unconditional minimum output is achieved for all outputs if and only if the system is minimum phase and the interactor matrix is diagonal. The test for a diagonal interactor is established by the following result. Proposition 6.2.2 Define matrix T ( q " 1 ) . d m i = mini<j< dkj, i.e. the minimum delay on row i of the r n Define the matrix D D ( q ~ ^ 1 ) 1 ) = as { k i j q - d - } Chapter 6. Control Loop Performance Monitoring in the MIMO case D(q~ ) is the first non-zero Markov parameter of T(q ). 1 1 90 Define D associated with D(q~ ) as D = {kijSij} where 1 6^ = lifd ij = = dl (6.188) nin 0 if not. (6.189) T(q~ ) has a diagonal interactor matrix 1 £(q) = diag(q ™™, ...,q ™™) d (6.190) d if and only if the associated matrix D has full rank. Proof: To prove this result let us express the transfer function as: T(q~ ) = {hjq-^h^q- )}, l (6.191) 1 where kij and dij are constant and delay between the input j and the output i. , / -IA n i j [ q _ i +M - +- + M ~ 1 m >~ l + a^-1 + ... + M ^ is a transfer function without delay between the input j and the output i. £{q) = diag(q ,q ) ai am is an interactor for T(q~ ) if and only if 1 limefaOTfar ) 1 = Um {kijq'^^hijiq- )} 1 =D D is nonsingular only if a* = min, d^. Otherwise if cn, < min^ d^ then the i (6.192) th row of D is a zero row therefore D is singular. If instead OJ; > min^ d^ then at least one infinite zero is present in D; therefore, it is nonsingular. • The next section deals with performance indices. Chapter 6. Control Loop Performance Monitoring in the MIMO case 6.3 91 M I M O Performance Indices Proposition 6.2.1 defines a unique absolute bound associated with the output yiit). If the process is minimum phase then this bound is achievable by using the lower triangle interactor and considering the output ordering (yi(t), ...,y (t)). n Furthermore, if the pro- cess is minimum phase with a diagonal interactor matrix, all outputs can achieve their lower bounds simultaneously. Therefore an absolute measure of controller performance may be built based on these bounds by defining the performance index associated with output yiit) as: 2 W ) ) = 2 ^ ( • (6.193) y yi{t))mv a Where CT .^ is the actual output variance of the output yi(t), and (o- (t)) 2 yi mv is the min- imum achievable output variance bound defined by equation (6.185). Note that this unique bound is achievable by at least one output. For an estimation of the overall performance, it is proposed the following reference bound which is unique and can only be achieved if and only if the interactor is diagonal. m (° y(t))ref = Y,K(t))rnv. 2 i=l (6.194) The M I M O performance index is defined as W ) ) = J^l2 ? M( • (6-195) It is clear that the procedure of assessment will be mainly based on equation (6.193). In many cases, not all outputs have the same importance, and it is thus useful to know what would be the achievable limit for an individual output. The benefits of the proposed performance indices over existing ones will be shown by means of the following example. Chapter 6. Control Loop Performance Monitoring in the MIMO case 92 Example 6.3.1 [19] (Dugard et al. (1984)) The process considered is described by ,-1 yi(t) q- 1 - 2q -2 0 Ui(t) .-3 u (t) 1 + + T T 2 l 0 ei(t) 1 e (t) 0 2 where £[ei(*) ei(*)] = of and E[e {t) e {t)} q~ Cl 2 = o\. 2 In the following, three different minimum variance control strategies, MVCi, and MVC$, MVC 2 obtained for three different interactor matrices, are considered. MVC\ 3 For the output ordering y(t) = (yi(t) y2(t)), the interactor is 0 q -q + 2q 2 3 q 3 Based on this interactor matrix, the output variances are E{bn(t)-ym?} = oi, E{[ (t)-y* (t)} } y2 and 2 and is obtained when a unitary MVC interactor matrix is considered. MVC : 2 are obtained respectively when considering the lower triangular interactor matrix with the output order (yi(t), y2(t)) and (y2(t), yi(t)). MVC\'. MVC 2 = (l + ^l)al 2 + o 2 2 E{[y(t) - y*(t)] [y(t) - y*(t)]} = (2 + 5cj)a + o\. 2 T If instead y(t) = (y (t) Vi(t)), then the interactor matrix is 2 0 q 6(?) = -(q 3 + 2q ) 2 q 3 Based on this interactor matrix the output variances are E{[yi(t)-yKt)} } 2 Emt) and = (i + cl)o , 2 - y* (t)?} = 4 E{[y(t) - 2 y*(t)f[y(t) - y*(t)}} = (1 + c\)a\ + o\. Chapter 6. MVCz' Control Loop Performance Monitoring in the MIMO case 93 The unitary interactor is determined as [ 0.35? + 0.35? - 0.35q 2 -0.35q 0.35q + 0.35g 3 - 0.35<? 3 -0.35g + 0.35g + 0.35g 3 2 3 The M V C based on this interactor matrix (equation 6.183) leads to E{[yi(t)-ym} } = (l + l4)ol 2 EiM^-ymn^A+o-i and E{[y(t) - y*(t)] [y(t) - y*(t)}} = (1 + \c\)ol T As can be seen, this control strategy achieves better overall performance than M V C l and M V C 2 . However, it does not provide an absolute minimum output variance for any of the individual outputs. Further, the unconditional minimum variance control is achieved for (t) Vl by M V C l and for y {t) by M V C 2 . 2 This example shows that the results of the M V C strategies differ depending on the type of interactor matrix considered. Now, consider them as three different controllers and compare them using the performance index as introduced in this chapter and also that proposed in [40, 36] (Huang (1997); Harris et al. (1996)). If these controllers ( M V C l , M V C 2 and M V C 3 ) were implemented on the process considered, the output variances would be those determined above. Therefore, according to the proposed approach these output variances should be compared to the following performance reference bounds: (°~yi )ref {vyzYref = = °\ f ° output yi and r 2 a f o r output y. 2 The overall performance will be compared to (o~ (t))f f = o~\ + o~\. y e The control performance of M V C l is first analyzed. Using the previous performance + cx . 2 Chapter 6. Control Loop Performance Monitoring in the MIMO case 94 index definition, the following results are obtained: \°~y\)ref = h\W)) = 7 2 \ y2)ref WW = (2 + 5c )a + o\ , J/2, °2 a 2 f O T (2 + 5c?)<r + a\ a„2 , J °2 2 2 N2 W r e / l+ - for an overall performance assessment. When applying the approach described in [36, 40], the reference output variances are those obtained when M V C 3 is implemented which are: 5 (°yi)ref = (! + g i > i for output y i , = ^ i°"i+ 2 for output y . . c c a 2 The overall performance is then compared to (cr (t))f f y (1 + f i ) ° " i + °2- Based on that = c e the performance indices would be cr 7—T2~ 5 2 ^P(?/I(*)) = = + 1 \ yi)ref a r^ w\ IpW)) r ( aw v2 , \ 1 output for (1 + 5c )a + a\ a = > QCI ° 2 2 = v — i 2 o y u 2 , o — t for output 2/2, (2 + 5c )a + <r a 2 2 2 2 The same analysis is carried out for the two remaining controllers (MVC2 and M V C 3 ) . Results are shown in Tables 6.2, where the notation used is: • U . I.: the performance indices based on the use of unitary interactor. • U . B . : the performance indices as defined in equation 6.193 and 6.195. Note that in Table (6.2) the values of I are expected to be greater or equal to one, but p in column 2 and 4 of this table, I (yi(t)) p and I (y2(t)) show values less than one. This p Chapter 6. Control Loop Performance ipMt)) i (yi(t)) Ip MVd MVC MVC P 2 3 Approach 0.66 1.2 1 U . I. Monitoring 1 1.8 1.5 U. B. 5.5 0.9 1 U . I. in the MIMO 95 case i (y(t)) P 6 1 1.1 U. B. 2.7 1.07 1 U . I. 3.5 1.4 1.3 U. B. Table 6.2: d = 0.9 o = \ a = 1 x 2 occurs because the M V C reference used is based on a unitary interactor which is not optimal when applied to individual outputs. The overall performance of MVC3 is better than that of MVC\ and MVC 2 and is very close to the proposed performance references (compare column 6 with column 7). This example shows that assessment based on the bounds defined in Section 3 is accurate. In the next section the steps for estimating the control loop performance based on previous bounds are given. 6.4 Estimation of the lower bounds Proposition 6.4.1 The bounds defined by equation (6.185) are feedback invariant. They can be estimated from routine closed-loop operating data. Proof: Express T(q~ ) as T(q~ ) = [2 ...T ...T ] and i V f a ) = [Wi.../y ..JV ] where 1 -1 x 1 Ti and Ni are the i th y(t) = (yi(t),y (t), 2 i m i m row respectively of T and N. Define y(t) = (yi(t), ...,y (t)) m ...,y (i)). m It is easy to see that y(t) = My(t) and where the (1, i) component of M is equal to one (i.e M\ = 1), M ^ i = 1 and ti 0 0 0 .. 1 0 M 0 1 0 .. 0 0 0 0 0 .. 0 .. then let the feedback control law be u(t) = —Ky(t). Substituting this equation into Chapter 6. Control Loop Performance Monitoring in the MIMO 96 case equation 6.181 yields y(t) = ((/ + TKy'NYq-^eit). (6.196) Multiplying equation (6.181) by M yields y(t) = My(t) = MTiq-^uit) + MN^q-^eit). (6.197) The new output ordering is y(t) = (yi(t), y2(t), ...,y (t)). The new transfer function is MT{q~ ) 1 m l = [T ...T ...T ] i 1 and the new noise term is MN^ ) m = [N ...N ...N ]. i 1 The m lower triangular interactor matrix £ (q) is determined for the new transfer function MT. m Multiplying equation (6.197) by £ (q) yields m U{q)My{t) = U(.q)MT(q- )u(t) l + U(q)MN{q- )e(t)- (6-198) l The noise term is expressed as U(Q)MN(a- ) 1 = F(q) + R(q- ). (6.199) 1 Substituting equation (6.196) and equation (6.199) into equation (6.198) yields U{q)My{t) = [R- UMTK{I + T K ) " i V ] ( g - ) ( t ) + F(q)e(t). 1 e Because the matrix H(q~ ) = (R-£ MTK(I x (6.200) 1 + TK)- N)(q~ ) 1 m 1 is proper, the first term on the right hand side of equation (6.200), which depends on the controller, contains only present and past values of e(t). The second term, which does not depend on the controller, contains only future values of e(t). Therefore, the two terms are independent. The first row of equation (6.200) is given by dl .(t) q iny = [R- UMTK(I where [R- £ MTK(I l + F (q)e(t), (6.201) 1 + T K ) ' ^ is the first row of the matrix H = R - £ MTK(I 1 m TK)~ N + TKy'NUq-^eit) + m and F = F (q) is defined in equation 6.186. Equation 6.201 is equivalent to 1 l yi{t) = (F* + ... + Fk . <THe(t) + (a mm di mm . q~<^ +1 + ....) (t). e (6.202) Chapter 6. Control Loop Performance Monitoring in the MIMO case 97 From the right hand side of equation (6.202) the first term, which is the i th process output under minimum variance control, is not affected by any feedback control. Therefore, modelling the process output y(t) by a multivariate M A process y(t) = $ a(t) + $ q- a(t-l) +... = # ( g 1 2 1 where a(t) is zero mean white noise, the i th Vi(t) = ((ft)* + - + . )*q- ^)a(t) d _ 1 K (6.203) process output is + . mm + 1 ) <-- 1 g + ...)a(t). (6.204) min Thus, the first terms on the right hand side of equation (6.204) can be estimated as (<j>iy = E[y (t)al }Q- , i j 1 where Q = E[a a[], then the lower bound can be determined as t «( u = (tiro-mr+...+(4. t) ^ v ' YQMS mtn . mm rr (6.205) m Based on this, the procedure of control loop performance assessment consists of the following steps. Algorithm for M I M O conrol loop assessment: 1. Determine the delays between process inputs and outputs. 2. Determine the process output M A or A R model to obtain the estimated sequence a(t). (equation (6.203)) 3. Determine the estimated lower output variance bounds of different outputs by use of equation (6.205). 4. Determine the performance indices associated with different outputs (equations (6.193) and (6.195)). Chapter 6. Control Loop Performance Monitoring in the MIMO case 98 To know whether those bounds are all achievable it is necessary to determine the matrix D and determine whether the interactor is diagonal. A l l steps in this procedure can be performed on line with a small amount of computation. In the following, two simulation examples are considered. 6.5 Simulation examples 6.5.1 First example This first example is from Harris et al.[36]. The process considered is a catalytic packedbed reactor having the following transfer function 5.158<r -6.643<r +1.8<r l+0.3412 - -1.04g5 4 2 G{q- ) X = 3 0.5384gl+0.5203g- -1.262q- 1 2 9 -3 0.1143<7- +0.1592<?1-0.6143?6 -1.012?- - 3.398?- - 0.7498g" 5 4 1 5 1 The disturbance process N(t) is a multivariable Integrated-Moving Average process of order (1,1): -l N(t) = 1 + 0.126?0.212?- 1 1 -0.1249?- 1 1 + 0.126?- 1 - q' 0 0 1 - q- 1 1 e(t), 1 where e(t) is a Gaussian noise sequence with covariance matrix given by: Q= 0.134 0.043 0.043 0.2 Though this process is nonminimum phase, the assessment described above can still be used for comparison purposes [36]. The objective of this example is to compare the bounds obtained when the unitary interactor matrix is used to those obtained when using equation (6.185). Chapter 6. Control Loop Performance Monitoring in the MIMO 99 case The unitary interactor associated with this process is given by: 0.355 - 0.355 q~ 1 0 q l 3 -0.4166 0.8521 + 0.1479 q' 0.4167 q 1 4 q 4 0.852071 q + 0.147928 q -0.355 q + 0.355 q -0.35497 q + 0.35507 q 0.1479 q + 0.8521 q 3 4 3 3 4 4 3 4 then F(q), i?(<? ) defined in equation 6.184 and ( £ F ) ( < 7 ) defined in equation 6.183 _1 -1 -1 are given by 1.126g + 1.126<7 + 1.093q + 0.148<? -0.125q - 0.125q - 0.157q + 0.355q 2 F(q) = 4 3 2 0.211 q + 0.211 q + 0.225 q + 0.355 q 2 3 4 1.1265 + RiQ- ) = ^ 4 0.609 q + 0.609 q + 0.62 q + 0.852 q 2 3 4 0.1249 -0.1249 ' 1 0.2113+ 3 0.6096 + (-1+9) 0.6096 (-1+9) and 1 + 0.241 q~ + 1.126 g 3 - 2 + 1.126g" 1 0.58 q~ + 0.211 q~ + 0.211 q3 2 1 0.197q~ - 0.124q~ - 0.125?3 1 2 1 + 0.475 q~ + 0.609 q~ + 0.60993 2 Therefore, the output covariance matrix is 0.551 0.127 E{[y(t)-y*(t)}[y(t)- *(t)] } y T 0.232 0.420 Based on this strategy the output variances to take as references against which the control loop performance will be compared are: KUf = E{( (t) = 0.551, yi - yl(t)) } 2 (6.206) (6.207) 1 Chapter 6. Control Loop Performance « W Monitoring = E{(y (t) in the MIMO - 2 case 100 (6.208) y* (t)f) 2 0.420 (6.209) and KW 2 (6.210) E{[y{t)-y*{t)) [y{t)-y\t))} = T = 0.551 + 0.420 (6.211) = 0.971. (6.212) When applying the approach introduced above to determine the output variance references, the minimum delay in each row is three; therefore, equation 6.186 for each output is q'N^q) = q N (q) 3 = 2 q + 1.126? + 1.126? -0.125? - 0.125? + 0.211 1-g- 0.211 1-g 0.212? + 0.212? ? + 0.6096? + 0.6096? + 0.211 0.6096 1-q- l-q- 2 3 2 2 3 2 1 1 - 1 1 Therefore the lower output variance bound for yi(t) is K)ref = E{( (t) = 0.4558, = E{(y (t) = 0.3827. Vl - y{(t)) } 2 (6.213) (6.214) and, for y (t) 2 (<4W 2 -y* (t)) } 2 2 (6.215) (6.216) The overall output variance bound would be KW 2 = = E{(y(t)-y*(t))*} 0.8385. (6.217) (6.218) Chapter 6. Control Loop Performance Monitoring in the MIMO case 101 The results of this example are summarized in Table 6.3. It can be seen that the references given in equation (6.207)-(6.209) are bigger than those established in equation (6.214)(6.216). That is because the M V C based on the unitary interactor does not achieve an ultimate minimum variance for the individual output. Note, also that these results are very close to those obtained using the approach introduced here where the knowledge of the interactor matrix was not used. Table 6.3: Reference values for example 2. references for yi(t) references for y (t) references for y(t) 0.551 0.4558 0.420 0.3827 0.971 0.8385 Unitary Int. U . bound Unitary Int. U . bound Unitary Int. U . bound 2 6.5.2 Second example For this example the process considered is described by -0.2(?i -0.12( —yj.J-t.il +0.12g- +0.2g1 0.06(7- 3 l-0.86(j- -|-0.12q- +0.2g- 2 ? 2 l-0.86g- +0.12 - +0.2g1 2 3 1 2 l-0.86(j- -f0.12q- +0.2(j- 3 1 g 1-0.369l-0.86q- +0.12g- +0.2(j- 2 ui(t) 3 J 3 u {t) 2 1 + 1 2 3 (t) ei u 1-OAq- 1 0 l-0.86g- +0.12 - +0.2 1 2 g 3 g where E ei(t) (t) ei )\=Q e (t) 2 By choosing the output ordering (yi(t), 0.42 0 0 0.38 = y (t)), the interactor matrix is 2 q -2.5g - 0.5? 3 0 2 <t 3 Chapter 6. Control Loop Performance Monitoring in the MIMO case 102 The M V C based on this interactor leads to the output tracking error yi(t)-yl(t) 1 0 y2(t)- m -1.250- + 1.02g" , (6.219) 1 y 1 + 0.469- + 0.27g- 2 1 2 C2(*) and the output variances are E{[yi(t)-yl(t)} } = 0.42 and 2 y»(t)?} E{\jfy(t) - = 1-59- To apply the proposed control loop performance assessment procedure, it is needed to know the minimum delay for each output and a M A or A R model for process outputs. In this case the minimum delay for both outputs is one, and equation (6.219) is the M A model needed. Therefore, according to equation (6.185), the minimum achievable output variances for outputs y\ and y are 2 d . 1 E{(m(t) E{(y (t) 2 - yl(t)) } = 2 - y* (t)) } 2 2 mtn J2 jQ( j) F F = 0-4229 and T j=i d? . = £" F Q(F ) 2 = 0.3802. 2 T j=i Therefore, the performance index for y\ is I (yi(t)) = 1 and that of y is I (y (t)) p If instead y(t) = (y (t) 2 2 p 2 = 4.19. yi(t)) then, the interactor matrix is o = -0.40? + 0.080g 3 2 The M V C based on this interactor leads to the output tracking error 1 + O.bq' + O.Slq- -0.18Q- + 0.07q- 0 1 1 vi(t)-yi(t) y2(t)-y* (t) 2 2 1 2 (6.220) e (*) 2 Chapter 6. Control Loop Performance Monitoring in the MIMO case 103 and the output variances are E{[yi(t) - y*M } = 0.5838 and E{[y2(t)-y* (t)?} = 0.3802. 2 2 In this case the minimum delay for both outputs is still one and equation (6.220) is the M A model to be used in the same way as before to obtain that the minimum achievable output variances for outputs yi and y are 2 d . 1 min E{( (t) - yl(t)) } = 2 yi F}Q(F}) £ = 0.4229 and T d . 2 E{(y2(t) - y* (t)) } = 2 2 2:^(^ = 0.3802. Therefore, the performance index for y\ is I (yi(t)) = 1.39 and that of y is I (y (t)) p 2 p 2 = 1. The process described above is now simulated with the following controller. " ) 1 0.04a~ k , C(q 1 2 -0,0344q~ 3 +0.0048Q~ 4 4-0.0080a~ - 0 . 0 6 O ~ 5 fc (0.02+0.019<l-2)(l_q-l) 1 2 + 0 . 0 5 1 6 Q ~ 3 - 0 (0.02+0.019q-2 0 0 7 2 q - 4 - 0 . 0 1 2 q ~ 5 ) ( 1 - q - 1 ) = . 0.5 — O . O a q 2 - 1 -0.284q~ +0.148q~ 2 3 +O.Q8q~ (0.02+0.019q-2)(l-q-l) 4 —0 . 2 + 0 . 0 5 2 q ~ . 2 1 +0.0792q~ 2 (0.02+0.019q-^ — 0 . 0 5 4 4 q ) ( l - q - l ) 3 - 0 . 0 2 4 q ~ 4 - The design of this controller was done in two stages. The first stage dealt with decoupling; then a proportional plus integral controller was designed for each loop. The simulation results are obtained when applying the procedure of assessment given above. Delays are considered known. The A R filter is determined using "arx" command of Matlab. Once the parameters of the filter are known, its recurrent form is used to obtain the residuals on line then updating the minimum output variances and output variances are updated considering a forgetting factor in the range of [0.9 0.99]. These simulation results are shown in Figure 6.17. It can be seen from combining all the results obtained above with different control strategies applied to the considered process that at least one output can achieve a performance index of one (first two controllers). The performance index for all outputs can Chapter 6. Control Loop Performance Monitoring in the MIMO case 104 be close to one for a good controller such as the above last one. The proposed procedure is able to track on-line the control loop performance. -i i 1 _l I 1_ i r - 1 ~I I 200 I 300 ! 1 I 4oo eoo eoo Figure 6.17: On line performance index estimation kl=0.08 and k2=0.05. 6.6 Conclusions It was shown that the proposed multivariable procedure for control loop performance assessment is an improvement over existing techniques. No prior knowledge of the interactor matrix is needed, and only knowledge of the minimum delay for each output is needed. The method is then a natural extension of SISO control loop performance assessment to M I M O control loop performance assessment. This procedure can be used to identify the loops showing a large output variance and so in need of special attention, from those having an output variance close to the M V C benchmark. A significant advantage of the bounds derived here is that they are independent of the interactor matrix. In the following chapter the proposed method is evaluated by application to the loops Chapter 6. Control Loop Performance Monitoring controlling a lime kiln in a Kraft pulp mill. in the MIMO case 105 Chapter 7 Performance Assessment of Control Loops on a Lime Kiln In this chapter, the method proposed in the last chapter is evaluated by application to the loops controlling a lime kiln in a Kraft pulp mill. The data was obtained in collaboration with engineering staff at Canfor, Prince George and P A P R I C A N , Vancouver Laboratory. The organization of this chapter is as follows. In the first section the lime cycle is described. This is followed by the development of a dynamic model of the lime kiln. The results of the assessment procedure are then given in section 4, followed by a conclusion in the last section. 7.1 Process description The lime cycle is an important part of the kraft recovery and recausticizing system in a pulp mill. The purpose of the cycle is to regenerate the alkaline white liquor from the spent green liquor recovered from the pulping process, see Figure (7.18) from Dumont [20]. The white liquor, which is a solution of sodium hydroxide (NaOH) and sodium sulfate (Na2S), is used to break down the lignin in the wood. The green liquor is a solution of sodium carbonate (Na2COz), sodium hydroxide (NaOH), and sodium sulfite (Na S). 2 The first step in the lime cycle is to add lime (CaO) to the green liquor where it reacts with water in solution to produce calcium hydroxide. This reaction is called slaking and 106 Chapter 7. Performance Assessment of Control Loops on a Lime SLMUU CL A RIME It LIHE fUEL 107 WHITE WHITE LIQUOR C A U 5 T I C I Z CR L14 U 4 ft Kiln LIQUOR HUD WASHER HUD Figure 7.18: Lime cycle. can be written as: CaO + H 0 — • Ca(OH) 2 2 + heat. The calcium hydroxide then reacts with the sodium carbonate in the green liquor to produce sodium hydroxide and calcium carbonate. This reaction is called causticizing and can be written as Na C0 2 The white liquor (NaOH) 3 + Ca(OH) 2 — • 2NaOH + CaC0 . 3 is separated from the calcium carbonate (lime mud) by passing the solution through a clarifier, where the calcium carbonate settles out as an insoluble mud. To complete the cycle, the calcium carbonate or lime mud is washed and filtered, then dried and heated in large rotary kilns. The heat causes the calcium carbonate to break down, giving off carbon dioxide gas and calcium oxide or lime. This reaction can be written as CaC0 3 + heat —• CaO + C0 2 T. Chapter 7. Performance Assessment of Control Loops on a Lime Kiln 108 The lime produced is then used in slaking. From the above, it is evident that the lime kiln is a crucial element of the lime cycle. 7.2 Lime kiln process FUME S ID FAN LIME MUD HO T END Figure 7.19: Lime Kiln. The lime kiln, depicted in Figure (7.19), is a long rotating cylinder, sloping typically at four centimetres per meter, and rotating at one revolution per minute. The lower end of the kiln is referred to as the hot end or the feed end. The fuel and air for combustion enter at this end. The lime mud enters into the kiln at the feed end or the cold end. It is then dried and converted into hot lime as it flows down the kiln toward the hot end, which is maintained at a temperature of about 1000°C. Detailed descriptions of the lime kiln can be found in [61]. Clearly, lime kilns are multivariable distributed systems. The objective of a lime kiln control system is to produce lime of good reactivity while minimizing the energy Chapter 7. Performance Assessment of Control Loops on a Lime Kiln 109 consumption. Underburned or overburned lime from an ineffectively controlled kiln can be very costly in terms of reduced pulping efficiency and increased energy costs. A review of the principles and current practice of lime kiln control is given in Dumont [20]. The manipulated variables may include mud flow rate, fuel rate, rotational speed of the kiln, primary air flow (FD Fan), and draft (ID Fan). However, setting of the mud flow rate depends on the availability of stored mud and the desired lime production rate and so may not be adjustable for control purposes. Typically, the controlled variables include the feed end and the back end temperatures. The quality of the lime produced depends upon the temperature profile in the longitudinal direction of the kiln. Another control concern is the oxygen content of the fume gasses, which indicates how efficiently fuel has been burned. For the K i l n under consideration, the controlled variables (outputs) are: feed end temperature (FET), hot end temperature (HET), and the excess of oxygen (XS 02)-The rotational speed of the kiln is adjusted according to the load. The manipulated variables (inputs) are the fuel rate and air flow (ID fan). The lime kiln is considered here as a multivariable system with two inputs and three outputs. The models are assumed to be first-order with time delay determined from the two bump tests shown i n figure (7.20) and figure (7.21). The resulting transfer function model is shown below: 0.52e87s+l HET FET Xs0 2 = 67s -0.002elOs+1 0.062e10.52s+l 0.228e- 23s+l 13s -0.0088e-°2.5s+l 7s FUEL 1 7 s FAN 5 67s 12s 0.018e-° 1.75s+l (7.221) where all time constants and delays are expressed in minutes. Dynamic operating response data were collected and compared to the response of the model determined above. The model delays show that the excess of oxygen can be adjusted with less than a minute of delay. However, a delay of at least 12 minutes from input to any change in the hot end temperature output and 5.7 minutes for the feed end temperature. In Figure (7.22), the Chapter 7. Performance Assessment of Control Loops on a Lime E Q. < 10 20 30 40 Tim e ( m n ) Figure 7.20: First bump test. Kiln 110 Chapter 7. Performance Assessment of Control Loops on a Lime Kiln 9 d LL E Q. z < 150 —1 1000 980 960 100 Tim e (mn) Figure 7.21: Second bump test. 200 111 Chapter 7. Performance Assessment of Control Loops on a Lime Kiln 112 true response curves from mill data in solid line, while those obtained from the model are in dashed line. It can be seen that the model provides a reasonably good approximation of the kiln's response curves. In the next section the control performance of this process is investigated. 1 100 III \ 500 380 1000 ~* / 1500 - 2000 2500 3000 3500 4000 i i i i 1 i i I 1 2500 3000 3500 O>360 3-340 £±320 300 500 i 1000 i 1500 2000 Tim 4000 e ( m n ) Figure 7.22: Validation of the model. 7.3 Loop performance assessment for the lime kiln The control loops of this process are assessed using the approach described in the last chapter. It consists of the following steps: 1. Determine the delays between process inputs and outputs (from the model above). Chapter 7. Performance Assessment of Control Loops on a Lime Kiln 113 2. Determine the multivariate M A time series model for the process outputs (y - y)(t) = $ ( t ) + ^q^ait - 1) + ... = H(q- )a , 1 lfl t where a(t) is a zero mean white noise disturbance, q is the unit delay operator and y is the mean of y. The i process output is th (lU ~ Vi)(t) = ((<!>{)* + ... + . Tq- ^)a(t) + ((cb . d (tit min m i n di ' + 1 ) q - « - - » + ...)a(fi).222) 3. Determine the estimated lower output variance bounds of different outputs by use of the following equation ("J-*)™ = WG[($)T + + 3 1 W m . -i)*QW i. -,Y) , mm l T d (7.223) mm where Q = E[a(t)a(tY] the covariance matrix of the filtered data and ((b))* can be estimated by ((b)Y = E[(y -y ) al }Q-\ i i t j 4. Determine the performance indices associated with different outputs using W ) ) = jH r , (7-224) where 1 < I with higher value indicating poor control. p These steps can easily be performed on-line. Neither the nature of the control strategy nor the process model is needed, but the delays between different input output pairs of the process are required. In this case these delays, expressed in minutes, are obtained from the model above. The delays required by the assessment procedure are the smallest ones from each row of the above model. The on-line output variance is computed using the following equation K)N = A(tT )jv_i 2 + (i - x)( (N) Vl - m)?, Chapter 7. Performance Assessment of Control Loops on a Lime 114 Kiln where A e [0.98 1] is the forgetting factor. The resulting on-line performance index for each output is shown in figure (7.23). This preliminary study shows that the performance 41 I 1500 30 I 1 i 2000 1 1 I 1 5 0 0 2 0 0 0 1 1 i 2500 3000 i 1 1 I 2 5 0 0 I 3 0 0 0 1 r i 3500 i 4000 1 1— I 3 5 0 0 I I 4 0 0 0 Tim e (m n) Figure 7.23: On-line performance indices. index of the hot end temperature is in the range of [1.8 3.2] , that of the feed end temperature is within [4 8] and that of excess of 0 2 is within [10 20]. The conclusion to draw at this point is that the hot end and the feed end temperatures are relatively well controlled while the control of the excess of oxygen may not be. Note however that there will be a contribution to the performance index arising from the frequent setpoint change. Based on this study it was decided to design a new controller based on MPC (model predictive controller) to replace the one in place which is a combination of PI and Smith predictor. Chapter 7.4 7. Performance Assessment of Control Loops on a Lime Kiln 115 Conclusion The application of control loop performance monitoring presented here shows the important role that the monitoring tools can play in maintaining high process performance. By comparing the minimum variance level that each output can achieve with the actual variance, a decision can be made to accept, retune, or replace the control algorithm. The assessment was performed using routine operating data and little prior knowledge of the process. The approach proposed in [25] follows the spirit of control loop performance assessment initiated by Harris [35] for single-input single-output systems. Chapter 8 Summary 8.1 Contributions of this thesis The purpose of this thesis is to develop tools and techniques for monitoring and diagnosing the performance of control loops based primarily on operating data. The contributions of this thesis to this topic are: • In chapter 3, the development of a Least-Squares method for parameter identification and time delay estimation for processes where the acting disturbances are white noise, is presented in an attempt to overcome the difficulties of the existing method for time delay estimation. • In chapter 4, Harris's performance index was broadened as a measure of output variance to consider other criteria, so obtaining a better overall assessment of control loop performance. It was shown that oscillations result in poor performance and proposed an approach to controller tuning based on a combination of performance criteria and output variance. • In chapter 5, a unified oscillation index is proposed that can be used for linear and nonlinear cascade loops to help localize the loops causing oscillation at a given frequency. The procedure does not require knowledge of the transfer functions involved. Analysis of responses is used to determine the set of frequencies at which further invistigation is needed. 116 Chapter 8. Summary 117 • In chapter 6, a multivariable procedure for control loop performance assessment is proposed. This method leads to good results. No prior knowledge of the interactor matrix is needed and only knowledge of the minimum delay for each output is required. 8.2 Recommendation for future work The issues addressed in this thesis can be viewed as background for several interesting research directions. In chapter 3, the closed-loop method developed needs further study specially when the acting disturbances are described by colored noise. In chapter 4, by considering the presence of a nonlinear element in a given loop, it was possible to narrow the useful range for the performance index and to show that the performance index under normal operation is finite. The need to develop a framework for nonlinear systems is still an open problem. In Chapter 5, an oscillation index was introduced for SISO systems (cascade loops only); therefore, an extension of that study to multivariable systems is needed. In chapter 6, the M I M O case was addressed, but there remains a need to investigate the weight assigned to each performance index for practical use. The topic of control loop performance monitoring is relatively new and is an active area of research. Some other areas for future work would be to study the effect of the sampling rate on the estimation of the minimum achievable output variance, and to develop a confidence interval for performance index in both cases SISO and M I M O . This is needed for on-line performance assessment to compensate for varying time delay or parameter uncertainties and therefore, for robustness. Another area of special interest to paper industry is the machine cross directional control and it would interesting to extend the M I M O performance index defined in this thesis to cover such processes as well. Bibliography Astrom K . J . "Computer Control of a Paper Machine an Application of Linear Stochastic Control Theory". IBM J. Res. Dev., 11:389-405, 1967. Astrom K . J . "Assessment of Achievable Performance of Simple Feedback Loops". Int. J. Adapt. Cont, 5:3-19, 1991. Astrom K . J . and B . Wittenmark. "Computer Control System, Theory and Design". Prentice-Hall, Inc, 1986. Astrom K . J . and B . Wittenmark. "Adaptive Control". Addison-Wesley, 1995. Astrom K . J . and T. Soderstrom. "Uniqueness of the Maximum Likelihood Estimates of the Parameters of an A R M A Model". IEEE Trans. Automat. Contr., vol. AC-19, No. 6:769-773, 1974. Astrom K . J., C. C. Hang, P. Persson and W . K . Ho. "Towards Intelligent P I D Control". Automatica, 28:1-9, 1992. Basseville M . "Detection changes in signals and systems - a survey ". 24(3):309-326, 1988. Basseville M . and I. Nikiforov. Automatica, "Detection of abrupt changes: Theory and Applica- tions". Prentice-Hall, Inc, 1993. Bialkowski W . L . "Process Variability, Control, Standards and Competitive Program". Proc. ISA/92 Canada, Toronto, Canada, pages 395-418, April 1992. Billings S. A., J.O. Gray and D.H. Owens . "Nonlinear system design". I E E Control Engineering Series 25, 1984. Bittanti S., P. Colaneri and M . F. Mongiovi. "The Spectral Interactor matrix for the Singular Ricati Equation". Proc. 33rd CDC, pages 2165-2169, 1994. Borisson U . "Self-Tuning Regulators for a Class of Multivariable Systems". matica, 15:209-215, 1979. Box G . E . P. and G . W . Jenkins. "Time Series Analysis forcasting and Auto- Control". Holden-Day, 1976. Carter G . C , A . H . Nuttall, and P. G . Cable. "The Smoothed Coherence Transform" . Proc. of the IEEE, Oct. 1973. 118 Bibliography 119 Clarke D. W . "Generalized-least-squares estimation of the parameters of a dynamic model". In Identification in Automatic Control Systems. IFAC Symposium, Prague, 1967. Clarke D. W . "Sensor, Actuator, and Loop Validation". IEEE Control Systems Magazine, pages 39-45, Aug 1995. Cook P. A . "Nonlinear dynamical systems". Prentice-Hall, Inc, 1986. Desborough L . and T. Harris. "Performance Assessment Measures for Univariate Feedforward/Feedback Control". Can. J. Chem. Eng., 71:605-616, 1993. Dugard L . , G . C . Goodwin and X . Xianya. "The Role of Interactor Matrix in Multivariable Stochastic Adaptive Control". Automatica, 20:701-709, 1984. Dumont G. A . " A Review of the Principles and Current Practice of Lime K i l n Control". Course notes. Elliott H . and W . A . Wolovich. " A Parameter Adaptive Control Structure for Linear Multivariable Systems". IEEE Trans. Aut. Cont., AC-27:340-352, 1982. Elnaggar A . , G . A . Dumont and A . L . Elshafei. "Adaptive Control with Direct Delay Estimation". Control Systems '92: Dream vs. Reality, Whistler, B.C., pages 13-16, Sept. 1992. Ender D. B . "Implementation of Deadband Reset Scheduling for the Elimination of Stick-Slip Cycling in Control Valves". In Tappi Process Control, Electrical, & Information Conference, pages 83-88, 1997. Ettaleb L . , G. A . Dumont and , M . S. Davies. " A n Extended off-line Least-Squares Method for Parameter Identification and Time Delay Estimation". Procedings of the 37th IEEE Conference on Decision and Control, Dec. 16-18, 1998. Ettaleb L., G. A . Dumont and M . S. Davies. "Performance Assessment of Single and Multivariable Feedback Controllers: a Unified Approach". submitted to Automatica. Ettaleb L . , M . S. Davies, G. A . Dumont and E . E . Kwok. "Monitoring Oscillations in a Multiloop Systems ". Procedings of the 1996 IEEE International Conference on Control Applications, pages 859-863, 1996. Ettaleb L., M . S. Davies, G. A . Dumont and E . E . Kwok. "Realistic Performance Assessment in the SISO Case". Procedings of Mediterranean Conference on Electronics and Automatic Control, pages 114-117, Sep. 17-19 1998. Forssell U . and L . Ljung. 35:1215-1241, 1999. "Closed loop identification revisited". Automatica, Bibliography 120 Fu Y . and G. A . Dumont. "Optimum Laguerre Time scale and its On line Estimation". Int. J. adpt. Cont. Sig. Proc, 5:313-333, 1991. Fu Y . and G. A . Dumont. "On-line Evaluation of Control Loop Performance". Procedings of 4th IEEE cca, pages 655-656, 1995. Goodwin G. C. and K . S. Sin. "Adaptive Filtering Prediction and Control". Prentice- Hall, Inc, 1984. Goodwin G. C , P. J. Ramadge and P. E. Caines. "Discrete time multivariable adaptive control". IEEE Trans. Aut. Cont, vol. AC-25:449-456, 1980. Gunnarsson S., and B . Wahlberg. "Some Asymptotic results in Recursive Identification Using Laguerre Models.". Int. J. adpt. Cont. Sig. Proc, 5:313-333, 1991. Gustavsson I., L . Ljung and T. Soderstrom. "Survey Paper: Identification of Processes in Closed Loop -Identifiability and Accuracy Aspects". Automatica, 13:59-75, 1977. Harris T. J. "Assessment of Control Loop Performance". 67:856-861, 1989. Can. J. Chem. Eng., Harris T. J., F . Boudreau and J. F. MacGregor. "Performance Assessment of Multivariable Feedback Controllers". Automatica, 32(11):1505-1518, 1996. Harris T. J., J.F Macgregor and J. D. Wright. " A n Overview of Discrete Stochastic Controllers: Generalized PID Algorithms with Dead-Time Compensation". Can. J. Chem. Eng., 60:425-432, 1981. Horch A . and A . J. Isaksson. " A methode for detecting of stiction in control valves". In IFAC Workshop on On-line Fault Detection Process Industries, Lyon France, 1998. and Supervision in the Chemical Horch A . and A . J. Isaksson. "A modified index for control performance assessment". In American Huang B . Control Conference, Philadelphia, "Multivariate statistical PA, pages 3430-3434, 1998. methods for control loop performance assess- ment". Ph.D thesis, Dep. of Chemical Eng., University of Alberta, 1997. Huang B . and S. L . Shah. "Practical Issues in Multivariable Feedback Control Performance Assessment". In IFAC ADCHEM, Banff, Canada, pages 429-434, 1997. Huang B . , S. L . Shah and E . E. Kwok. "Good, Bad or Optimal? Performance Assessment of M I M O processes". Automatica, 33(11):1175-1183, 1997. Huang B . , S. L . Shah and H . Fujii. "Identification of the Time delay /Interactor matrix for M I M O systems using Closed Loop data". IFAC World Congress, 1996. Bibliography 121 [44] Isermann R. "Results on the Simplification of Dynamic Process Models". Int. J. Cont, 19(1):149-159, 1974. Isermann R. "Digital Control System". Springer-Verlag, 1981. Isermann R. "Fault diagnosis of machines via parameter estimation and knowledge processing - tutorial paper". Automatica, 29(4):815-835, 1993. Kanjilal P. P. and S. Palit. " On Multiple Pattern Extraction Using Singular Value Decomposition". IEEE Trans, on Automatic Cont, 43:1536-1540, 1995. Knapp H . C. and G. C. Carter. "The Generalized Correlation Method for Estimation of Time Delay". IEEE Trans. Acoust. Speech Sig. Process, 27:320-327, 1976. Kozuband D. J. and C . E . Garcia. "Monitoring and diagnosis of automated controllers in the chemical process industries ". In AIChE Meeting. St. Louise, MO, 1993. Landau I. D. and A . Karimi. " A n Output Error Recursive Algorithm for Unbaised Identification in Closed Loop". Automatica, 33(5):933-938, 1997. Ljung L . . "System Identification". Lynch C. and G. A . Dumont. Prentice-Hall, 1998. "Control Loop Performance Monitoring". IEEE Trans. Cont. Sys. Tech., 4,2:185-192, 1992. Mah R. S. H . . "Design and Analyse of Process Performance Monitoring Systems". In Chemical process control, Proceedings of the Engineering January 18-23, 1981. Matsubara M . "On the Equivalent Dead Time". IEEE pages 464-466, oct. 1965. Fondation Conference, Trans, on Automatic Cont., Nomikos P. and J. F . MacGregor. "Multivariate S P C charts for monitoring batch processes". Technometrics, 37(l):41-59, 1995. Perrier M . and A . Roche. "Towards Mill-Wide Evaluation of Control Loop Performance". Control Systems, Whistler, B.C., pages 51-60, Sept. 1992. Pryor C. "Autocovariance and Power Spectrum Analysis- Derive New Information from Process Data". Control Eng. 2, pages 103-106, Oct. 1982. Shah S., C. Mohtadi and D. W . Clarke. "Multivariable Adaptive Control without a prior knowledge of the Delay Matrix". Syst. Cont. Let, pages 295-306, 1987. Shinskey F. G . "Process Control System". McGraw-Hill, New York, 1988. Bibliography 122 [60] Smith K . W . " A Methodology for Applying Time Series Analysis Techniques". In Proc. ISA/92 Canada, Toronto, Canada, pages 419-435, 1992. Smook G. A . "Handbook for Pulp & Paper Technologists". Soderstrom T. and P. Stoica. "System Identification". T A P P I , C P P A , 1986. Prentice-Hall, 1989. Staton B . "Characteristics and Control of Process Disturbances". Process, 63:51-60, 1984. Strejc V . "Least Squares Parameter Estimation". Automatica, Hydrocarbon 16:535-550, 1980. Taha O., G . A . Dumont and M . S. Davies. "Detection and diagnosis of oscillation in cotrol loops". In 35th IEEE Conference on Decision and Control Kobe, Japan, pages 2432-2437, 1996. Tyler M . and M . Morari. "Performance assessment for unstable and nonminimum phase systems". In IFAC Workshop on On-line Fault Detection and Supervision the Chemical Process Industries, Newcastle upon Tyne, England, 1995. in Van den Hof P. M . and R. J. P. Schrama. "An Indirect Method for Transfer Function Estimation from Closed Loop Data". Automatica, 29(6):1523-1527, 1993. Weldon A. D. "Control valves- The key component in managing process variability". Tappi journal, 80(8):71-74, 1997. Wise B . M . and N . B . Gallagher. "The process chemometrics approach to process monitoring and fault detection". Journal, of Process Control, 6(6):329-348, 1996. Wolovich W . A . and P. L . Falb. "Invariants and Canonical Forms Under Dynamic Compensation". SIAMJ Cont. and Opt., 14:996-1008, Nov. 1976. Zheng W . X . and C. B . Feng. "Identification of Stochastic Time Lag Systems in the Presence of Colored Noise". Automatica, 26:769-779, 1990.
- Library Home /
- Search Collections /
- Open Collections /
- Browse Collections /
- UBC Theses and Dissertations /
- Control loop performance assessment and oscillation...
Open Collections
UBC Theses and Dissertations
Featured Collection
UBC Theses and Dissertations
Control loop performance assessment and oscillation detection Ettaleb, Lahoucine 1999
pdf
Page Metadata
Item Metadata
Title | Control loop performance assessment and oscillation detection |
Creator |
Ettaleb, Lahoucine |
Date Issued | 1999 |
Description | In this thesis an analysis of process control loop performance based on the Harris's performance index is given. Good and poor performances will be defined and it will be shown that poor performance may be obtained when some loops are cycling. Then a procedure for localizing the oscillating loops in a cascade multiloop system will be introduced. The procedure uses data collected under normal operation and assumes knowledge of the different frequencies of the periodic component of the output signal. An oscillation index is introduced to characterize the oscillating loops. The estimation of time delay is a major step in assessing control loop performance. To overcome some difficulties of the available methods for its estimation, a new off line extended Least-Squares technique for parameter identification and time delay estimation will be introduced for processes where acting disturbances are white noise (i.e. C(q⁻¹) = 1). The Harris performance index for single-input single-output (SISO) systems is a very useful tool, therefore its extension to multivariable systems becomes necessary. The multi-input multi-output (MIMO) performance indices available in the literature involve use of the MIMO time delay matrix, also called the interactor matrix, to perform the assessment. The method proposed in this thesis does not require knowledge of the interactor matrix. The method uses knowledge of the delays between different input/output pairs of the process and defines an absolute lower bound on the achievable output variance for each output. This bound is then used to define a performance index associated with each output. The sum of these bounds is used to characterize the overall control loop performance. The results are independent of the matrix time delay representation, or interactor. This method is evaluated by application to the loops controlling a lime kiln in a Kraft pulp mill. |
Extent | 4961909 bytes |
Genre |
Thesis/Dissertation |
Type |
Text |
File Format | application/pdf |
Language | eng |
Date Available | 2009-07-02 |
Provider | Vancouver : University of British Columbia Library |
Rights | For non-commercial purposes only, such as research, private study and education. Additional conditions apply, see Terms of Use https://open.library.ubc.ca/terms_of_use. |
DOI | 10.14288/1.0065169 |
URI | http://hdl.handle.net/2429/9963 |
Degree |
Doctor of Philosophy - PhD |
Program |
Electrical and Computer Engineering |
Affiliation |
Applied Science, Faculty of Electrical and Computer Engineering, Department of |
Degree Grantor | University of British Columbia |
Graduation Date | 1999-11 |
Campus |
UBCV |
Scholarly Level | Graduate |
Aggregated Source Repository | DSpace |
Download
- Media
- 831-ubc_1999-463419.pdf [ 4.73MB ]
- Metadata
- JSON: 831-1.0065169.json
- JSON-LD: 831-1.0065169-ld.json
- RDF/XML (Pretty): 831-1.0065169-rdf.xml
- RDF/JSON: 831-1.0065169-rdf.json
- Turtle: 831-1.0065169-turtle.txt
- N-Triples: 831-1.0065169-rdf-ntriples.txt
- Original Record: 831-1.0065169-source.json
- Full Text
- 831-1.0065169-fulltext.txt
- Citation
- 831-1.0065169.ris
Full Text
Cite
Citation Scheme:
Usage Statistics
Share
Embed
Customize your widget with the following options, then copy and paste the code below into the HTML
of your page to embed this item in your website.
<div id="ubcOpenCollectionsWidgetDisplay">
<script id="ubcOpenCollectionsWidget"
src="{[{embed.src}]}"
data-item="{[{embed.item}]}"
data-collection="{[{embed.collection}]}"
data-metadata="{[{embed.showMetadata}]}"
data-width="{[{embed.width}]}"
async >
</script>
</div>
Our image viewer uses the IIIF 2.0 standard.
To load this item in other compatible viewers, use this url:
http://iiif.library.ubc.ca/presentation/dsp.831.1-0065169/manifest