Control Loop Performance MonitoringChristopher B. LynchB.A.Sc. University of Waterloo, Waterloo, 1990A THESIS SUBMITTED IN PARTIAL FULFILLMENT OFTHE REQUIREMENTS FOR THE DEGREE OFMASTER OF APPLIED SCIENCEinTHE FACULTY OF GRADUATE STUDIESDEPARTMENT OF ELECTRICAL ENGINEERINGWe accept this thesis as conformingto the required standardTHE UNIVERSITY OF BRITISH COLUMBIADecember, 1992© Christopher Lynch, 1992In presenting this thesis in partial fulfilment of the requirements for an advanceddegree at the University of British Columbia, I agree that the Library shall make itfreely available for reference and study. I further agree that permission for extensivecopying of this thesis for scholarly purposes may be granted by the head of mydepartment or by his or her representatives. It is understood that copying orpublication of this thesis for financial gain shall not be allowed without my writtenpermission.(Signature)Department of ELEGTR^ea&m.E.eRIA.i6) The University of British ColumbiaVancouver, CanadaDate^ 1-1F:P--R 1 41 352-DE-6 (2/88)AbstractControl loop monitoring provides the opportunity to maintain high process performance. In thisthesis, the regulation performance of a process is monitored by examining the time delay, static input-output relation, and the minimum achievable output variance of the process. The time delay estimateis obtained from the Fixed Model Variable Regressor Estimator and the static input-output relation isestimated using the Adaptive Nonlinear Modeller. The minimum achievable output variance of theprocess is determined by modelling the closed loop noise transfer function by a Laguerre network andfinding the Laguerre gains of the network using Recursive Extended Least Squares estimation. Anextensive simulation study was performed to examine the operation, usefulness and limitations of themonitoring tools. These simulations confirmed the operation and benefit of the monitoring methods,and defined their limitations. The monitoring tools were applied to data generated by industrialprocesses. The successful application to industrial processes highlighted the benefits obtainable bycontrol loop monitoring. Due to the ease of applying the monitoring tools and the valuable insightthey provide, these monitoring tools should be applied to any process where regulation performanceis a concern.iiTable of ContentsAbstract ^List of Tables List of Figures ^Acknowledgments iiviviixiii1 Introduction 11.1 Background ^ 11.2 Motivation 21.3 Contributions ^ 21.4 Thesis Outline 32 Control Loop Monitoring 42.1 The Discrete Linear Model ^ 42.2 Monitoring Closed Loop Processes ^ 62.3 Purpose of Monitoring ^ 72.3.1^Time Delay 82.3.2^Static Characteristics ^ 92.3.3 Minimum Achievable Output Variance ^ 103 Time Delay Estimation 123.1 Fixed Model Variable Regression Estimator ^ 123.1.1^Derivation of FMVRE ^ 123.1.2 FMVRE Algorithm 143.1.3^Examples ^ 143.2 Estimation of Delay in Control Loops ^ 173.2.1^Effect of Stochastic Disturbances 173.2.2^Closed Loop Estimation ^ 19iii4 Static Characteristics 214.1 Adaptive Nonlinear Modeller ^ 214.1.1 ANM Algorithm 234.1.2 ANM Example ^ 244.2 Application to Stochastic Feedback Processes ^ 265 Output Variance 275.1 Minimum Achievable Output Variance ^ 275.2 Method of Finding Minimum Achievable Output Variance ^ 305.3 Estimation of Minimum Achievable Output Variance 335.3.1^Discrete Laguerre Model ^ 335.3.1.1 Laguerre Network State-Space Representation ^ 345.3.2^Recursive Extended Least Squares Estimation ^ 365.3.2.1 Recursive Extended Least Squares Algorithm ^ 375.3.2.2 Convergence ^ 385.3.2.3 Estimation of Laguerre Gains ^ 395.3.3 Minimum Achievable Output Variance Estimation Algorithm ^ 405.3.4 Example ^ 416 Simulation Study 456.1 Time Delay Estimation ^ 456.1.1^Unbiased Estimation 466.1.2^Estimated Delay of One Sampling Interval ^ 476.1.3^Arbitrary Bias ^ 476.1.4^Insufficient Excitation ^ 486.2 Causes and Solutions for Incorrect Delay Estimation ^ 496.2.1^Biased Delay ^ 496.2.2^Insufficient Excitation ^ 53iv6.3 ANM Estimation of Stochastic Feedback Processes ^ 576.4 Estimating Nonlinearities ^ 616.4.1 Saturation ^ 616.4.2 Deadzone 646.5 Effect of Process Changes on Minimum Achievable Output Variance ^ 666.6 Choice of Laguerre Time Scale and Number of Filters ^ 766.6.1 Real Closed Loop Noise Transfer Function Denominator Roots ^ 776.6.2 Complex Closed Loop Noise Transfer Function Denominator Roots ^ 826.6.3 Closed Loop Noise Transfer Function Numerator Not SPR ^ 876.7 Conclusions ^ 927 Industrial Data 957.1 Pressure Screen Reject Flow ^ 957.2 Reject Refiner Motor Load 1007.3 Kamyr Digester Chip Level ^ 1048 Conclusion^ 111References 114vList of Tables6.1 ANM results for estimation performed under closed loop operation ^ 586.2 ANM results for stochastic process using various low pass filters 606.3 ANM results for estimation of saturation ^ 636.4 ANM results for estimation of deadzone 657.5 The results of minimum achievable output variance estimation ^ 1087.6 The results of minimum achievable output variance estimation on the AGPC control . . ^ 109viList of Figures2.1 Block diagram representation of the discrete linear model ^ 52.2 Block diagram representation of closed loop process 62.3 Block diagram of the monitoring scheme ^ 73.4 Results of FMVRE time delay estimation with (a) presenting the estimation over timewhile (b) shows the E1 function at sample one hundred ^ 153.5 Delay estimation with varying forgetting factor, with (a) A=1, (b) A=0.99, (c) A=0.975,(d) A=0.95 ^ 163.6 The stochastic process used for the example ^ 183.7 Results of FMVRE time delay estimation with (a) presenting the estimation over timewhile (b) shows the E1 function at sample one thousand ^ 193.8 E1 functions for (a) the input and plant output and (b) the input and disturbance output . ^ 194.9 Block diagram of the Adaptive Nonlinear Modeller ^ 214.10 Results of the ANM estimation showing (a) the input used, (b) the resulting output, (c)the table value estimation over time and (d) a plot of the final estimated table values . . . 254.11 Results of using (a) the input and the ANM to (b) reconstruct the output ^ 255.12 Block diagram of Laguerre filter network ^ 345.13 Process used in the example ^ 415.14 The Laguerre gains estimated during minimum achievable output variance estimation ^ 425.15 The (a) estimate of the F(c1 1 ) polynomial and (b) the estimated noise autocorrelation ^ 425.16 The (a) FMVRE time delay estimation and the (b) estimated minimum achievableoutput variance and process output variance ^ 43vii5.17 Distribution of the results of five hundred simulations ^ 446.18 Process which results in correct delay estimate 466.19 E1 function from FMVRE estimation resulting in correct delay of seven ^ 466.20 Process which results in a biased delay estimate ^ 476.21 E1 function from FMVRE estimation resulting in a bias at one ^ 476.22 Process which results in a biased delay estimate ^ 486.23 E1 function from FMVRE estimation resulting in a bias at four ^ 486.24 PID controlled process which results in a biased delay estimate 496.25 E1 function from FMVRE estimation resulting in a bias at one ^ 496.26 The disturbance output autocorrelation from the process (a) resulting in a delay bias,Figure 6.20, and (b) resulting in a correct delay estimate, Figure 6.18 ^ 506.27 E1 functions for (a) the input and plant output and (b) the input and disturbance outputresulting in a bias at one ^ 506.28 E1 functions for (a) the input and plant output and (b) the input and disturbance outputresulting in a correct estimate ^ 516.29 E1 functions for (a) the input and plant output and (b) the input and disturbance outputresulting in a bias at four ^ 526.30 E 1 functions for (a) the filtered input and filtered plant output and (b) the filtered inputand filtered disturbance output resulting in correct delay estimate ^ 536.31 E1 function from FMVRE estimation using filtered signals, resulting in a correct delayestimate ^ 536.32 E1 functions for (a) the input and plant output and (b) the input and disturbance outputresulting in a bias at one ^ 546.33 Power spectral density function for the input signal from the PID controller ^ 55viii6.34 Power spectral density function for the input signal from the Deadbeat controller ^ 556.35 Magnitude response of the PID controller ^ 566.36 Process used to test ANM ^ 576.37 Signal used to estimate the static input-output relationship ^ 586.38 The (a) process input and (b) process output from closed loop operation ^ 59^6.39 The (a) process input and (b) process output from stochastic closed loop operation 596.40 The (a) filtered process input and (b) filtered process output used by the ANM ^ 606.41 Saturation nonlinearity to be estimated ^ 616.42 Process with saturation nonlinearity 626.43 The (a) process input and (b) process output from the saturation only ^ 626.44 The (a) process input and (b) process output from closed loop operation with Dahlincontroller ^ 636.45 The (a) process input and (b) process output from closed loop operation with gaincontroller ^ 646.46 Deadzone nonlinearity to be estimated ^ 646.47 The (a) process input and (b) process output from the deadzone only ^ 656.48 The (a) process input and (b) process output from closed loop operation with Dahlincontroller ^ 666.49 Original process and associated properties ^ 676.50 The (a) estimated noise autocorrelation function, (b) the FMVRE delay estimation, (c)the estimated F(q-1 ) coefficients and (d)^and a;, for the original process ^ 686.51 Distribution of the results of five hundred simulations for the original process ^ 686.52 The (a) estimated noise autocorrelation function, (b) the FMVRE delay estimation, (c)the estimated F(q-1 ) coefficients and (d) 6-,2m, and cr- 2 for a change in delay ^ 69ix6.53 The (a) estimated noise autocorrelation function, (b) the FMVRE delay estimation, (c)the estimated F(cf 1 ) coefficients and (d) 6.12nv and o-v2 for a change in plant ^ 716.54 Distribution of the results of five hundred simulations for the change in the processplant ^ 716.55 The (a) estimated noise autocorrelation function, (b) the FMVRE delay estimation, (c)the estimated F(cf1 ) coefficients and (d) (3-11Y and (7,2 for a change in disturbance filter . . . 736.56 Distribution of the results of five hundred simulations for a change in disturbance filter . . . 736.57 The (a) estimated noise autocorrelation function, (b) the FMVRE delay estimation, (c)the estimated F(q -1 ) coefficients and (d) i3"!iv and o-v2 for a change in controller ^ 756.58 Distribution of the results of five hundred simulations for a change in the controller . . . . 756.59 Gain controlled process ^ 776.60 The (a) estimated noise autocorrelation function, (b) the FMVRE delay estimation, (c)the estimated F(q -1 ) coefficients and (d) &!,, and o for process with real poles ^ 786.61 Distribution of the results of five hundred simulations for process with real poles ^ 786.62 Effect of changing the time scale and filter number on 6-1„v for process with real poles . . ^ 796.63 Dahlin controlled process with real poles ^ 806.64 The (a) estimated noise autocorrelation function, (b) the FMVRE delay estimation, (c)the estimated F(q -1 ) coefficients and (d)^and o- for Dahlin controlled process withreal poles ^ 806.65 Effect of changing the time scale and filter number on 6 -,2„v for Dahlin controlledprocess with real poles ^ 816.66 The (a) estimated noise autocorrelation function, (b) the FMVRE delay estimation, (c)the estimated F(q -1 ) coefficients and (d) (3-11v and o-2 for process with complex poles . . . 826.67 Effect of changing the time scale and filter number on 611v for process with complexpoles ^ 836.68 Effect of reducing filter number on (a) noise autocorrelation function and (b) theestimated F(q -1 ) coefficients ^ 846.69 Process controlled by Dahlin controller with complex poles ^ 846.70 The (a) estimated noise autocorrelation function, (b) the FMVRE delay estimation, (c)the estimated F(q-1 ) coefficients and (d)^and^for Dahlin controlled process withcomplex poles ^ 856.71 Effect of changing the time scale and filter number on ^for Dahlin controlledprocess with complex poles ^ 866.72 Effect of reducing filter number on (a) noise autocorrelation function and (b) theestimated F(q-1 ) coefficients ^ 866.73 Process controlled by Dahlin controller without SPR R(q -1 ) ^ 876.74 The (a) estimated noise autocorrelation function, (b) the FMVRE delay estimation, (c)the estimated F(q-1 ) coefficients and (d)^and o for Dahlin controlled processwithout SPR R(q 1 ) ^ 886.75 Distribution of the results of five hundred simulations for Dahlin controlled processwithout SPR R(q 1 ) ^ 896.76 Effect of changing the time scale and filter number on ^for Dahlin controlledprocess without SPR R(q-1 ) ^ 896.77 The (a) estimated noise autocorrelation function, (b) the FMVRE delay estimation, (c)the estimated F(q -1 ) coefficients and (d)^and o-2 for process without SPR R(q -1 ) . . . 906.78 Effect of changing the time scale and filter number on i3-.2„ for process without SPRR(q 1 ) ^ 917.79 The pressure screen reject flow input and output ^ 967.80 The output autocorrelation function and input power spectral density function ^ 97xi7.81 The FMVRE (a) estimation of the delay and E1 function for (b) all samples, (c)samples 300 to 700 and (d) samples 800 to 1200 ^ 977.82 The estimated noise sequence and F(q -1 ) coefficient for (a&b) ten filters and (c&d) sixfilters ^ 997.83 The estimated minimum achievable output variance using (a) ten filters and (b) sixfilters ^ 997.84 A double-disc refiner ^ 1017.85 The input and output signals of the reject refiner ^ 1027.86 The filtered input and output signals of the reject refiner 1037.87 Various results from the ANM estimator ^ 1047.88 A typical Kamyr digester ^ 1057.89 The output signals of the digester for (a) P control, (b) AGPC, (c) CMV control and(d) more AGPC ^ 1067.90 The autocorrelation functions of the output signals of the Kamyr digester for (a) Pcontrol, (b) AGPC, (c) CMV control and (d) more AGPC ^ 107xiiAcknowledgmentsI would like to express my thanks to my supervisor Dr. Guy A. Dumont for his guidance andsupervision during my study at the University of British Columbia. I also appreciate the advice ofmy colleagues in the Control Group at the Pulp and Paper Centre, especially Ashraf Elnaggar whoprovided a patient ear and was always willing to discuss my ideas. I wish to thank Bruce Allisonfor providing the industrial data used within this work. I would like to thank Dr. Dumont and theWoodpulp Network of Centres of Excellence for providing financial support and excellent computerfacilities. And to the people who are always there, my family and Katherine.Chapter 1: IntroductionChapter 1Introduction1.1 BackgroundA serious problem faced upon the completion of any engineering research is its acceptanceand use by industry. This problem is very pronounced in the field of control engineering [7].The reasons are numerous and the blame lies with both industry and the research community.The bottom line in industry is production, which dictates that processes and therefore controllersoperate continuously, with as few disruptions as possible. This bottom line will allow promisingcontrol methods to pass by untried for fear of disturbing operation, even though process performancecould be improved. Universities, however, must realize the needs and restrictions of industry andgauge research accordingly. Only once industry and universities work together will new research beintegrated into industrial processes.Control loop monitoring is a way to bring both sides closer together. Monitoring in its mostbasic form involves observing process operation, recording behavior and commenting on any unusualoccurrence. It is very passive, in the sense that no changes are required to the setup of the loop,so process operation is not disturbed. Meanwhile, monitoring offers the research community anopportunity to use various identification and modelling techniques in an industrial environment. Thusindustry is introduced to alternative methods and new research, while universities are given theopportunity to interact with industry. Successful results will demonstrate the usefulness of newmethods and should help encourage its acceptance and possible use in control algorithms. Anyunsuccessful results can be examined by researchers to determine problems, all without adverselyaffecting the process operation in any manner. Such cooperation would ultimately result in thefulfillment of the desires of both industry and the research community.A pulp and paper mill can have thousands of control loops [8] which typically are the responsi-bility of one or two process control engineers. With this large number of loops, it is next to impossibleto ensure that they all are operating satisfactorily and it is time consuming to investigate each loop1Chapter 1: Introductionindividually. This suggests that there is a requirement for methods that could detect which loopsneed the attention of the engineer [23]. The monitoring methods presented in this thesis can be usedto flag loops whose performance has deteriorated. The cost of implementing control loop monitoringis small, as it requires only the addition of the software code to manipulate the process signals. It isthis ease of implementation and the need of such methods in an industrial setting, that makes controlloop monitoring an attractive way to improve process performance.1.2 MotivationAn industry wide problem in the control of pulp and paper processes is that in practice manycontrol principles are ignored. As a result, many loop controllers including low level single-inputsingle-output controllers, are often poorly chosen and tuned [4]. If these low level processes operatepoorly due to the control, the overall operation of the process will suffer. One way of improving theoperation of these controllers is through monitoring the process.The goal of control loop monitoring is to provide information about a process during closedloop operation. This information can be used to confirm if the process is operating as expected andalso highlight any changes in the process. If the process is not operating as expected, then changesto the way it is controlled may be necessary to improve performance. By highlighting changes inthe process over time, changes which often signal deteriorating operation, maintenance of the loopcomponents can be undertaken to improve performance.1.3 ContributionsThree properties of single-input single-output control loops are monitored. The time delay of theprocess is monitored using the Fixed Model Variable Regressor Estimator proposed by Elnaggar [12].The static input-output relationship of the process is determined using the Adaptive NonlinearModeller proposed by Asti -Om [3] but the relationship is estimated using the Recursive Least Squaresmethod. Based on the work of Harris [17], the minimum achievable output variance of the controlloop is monitored. In this work, the delay is not assumed known and Laguerre filters are used withRecursive Extended Least Squares estimation for process modelling.2Chapter 1: Introduction1.4 Thesis OutlineThe types of processes considered and the notation to be used are outlined in Chapter 2. Aswell, the purpose of monitoring the three properties chosen is discussed and the problems associatedwith closed loop monitoring are outlined. Chapter 3 describes the theory behind the Fixed ModelVariable Regressor Estimator, lists its algorithm and presents some examples of its operation. Theproblems encountered when using the estimator on closed loop stochastic processes are discussed.The Adaptive Nonlinear Modeller used for the estimation of the static input-output relationship isstudied in Chapter 4. The various components required for determining the minimum achievableoutput variance are presented in Chapter 5. This includes discussion of Laguerre filters, RecursiveExtended Least Squares estimation and the method used to calculate the minimum achievable outputvariance. Chapter 6 contains numerous computer simulations of the various techniques studied.Practical aspects of the techniques are examined via computer simulation in order to gain furtherinsight on their operation. Then, Chapter 7 presents the results when monitoring is used on industrialprocess data. Finally, in Chapter 8, conclusions are drawn and summarized based on previous results.3Chapter 2: Control Loop MonitoringChapter 2Control Loop MonitoringIn this chapter, the foundation for the research is set, providing a starting point for the remainingchapters. First, the model used to describe an industrial process is given and so then the notation tobe used is defined. Next, the assumed operating setup of the process is given, which is followed by adiscussion on the importance of monitoring and the reason why the specific properties are monitored.Finally, this discussion includes benefits to be achieved from monitoring the chosen process properties.The importance of monitoring and the approach taken in this research will then have been presented.2.1 The Discrete Linear ModelThe discrete linear model which is used to describe processes consists of two terms, one whichrepresents the plant dynamics and another which accounts for the influence of the surroundings onthe plant. These effects are denoted separately by a single term and combined together because themodel used is linear. A block diagram form of such a model appears in Figure 2.1 and representsthe plant and disturbance of a process. The equation:A iB((€1-1C1-1))^kyo (t) = ^u(t) (2.1)completely describes the plant and its properties. The plant output is yo(t), its input u(t), and thepolynomials Ai (q-1 ) and B(q-1 ) contain information about its dynamics. The backward shift operatorq-1 is used and is defined as u(t — 1) = q—l u(t). The polynomial AI (q -1 ) is monic and of degree n,while the polynomial 8(q -1 ) is of degree m, where m n. The plant dynamics are determined bythe value of the roots of Ai (q -1 ) and 8(q-1 ), known as the poles and zeroes of a plant respectively.The last piece of information about the plant is contained in k, which is the delay. For sampled-dataprocesses, this delay is always greater than or equal to one, to ensure process causality.The disturbance term represents the effect of the environment in which the plant exists and cantake on many forms. It could be purely deterministic, such as randomly occurring steps changes, orit could be stochastic in nature. The Autoregressive—Integrated—Moving—Average (ARIMA) model4DisturbanceN(0,cse )^e(t)Plant^ ye(t)Yo(t)u(t) y(t)Chapter 2: Control Loop MonitoringFigure 2.1. Block diagram representation of the discrete linear modelcan be used to represent many different types of typical industrial disturbances. It makes up thesecond part of Figure 2.1 and the equation:ye (t) = C(q1) e(t)A2(cri)Vdcontains information about the properties of the disturbances. Its output is ye(t) while the input e(t)is a sequence of uncorrelated random variables with a zero mean and constant variance, which canbe denoted as N(00-,2). This input signal is then a white noise sequence with a normal distribution.The polynomials C(q-I ) and A2(q-1 ) have degree q and p respectively and like the plant polynomials,the roots of these polynomials determine the dynamics of the disturbance. It is also assumed that theroots of the A2(q -1 ) polynomial lie inside the unit circle. The operator, V, is defined as V = 1 —and the term Vd permits the mean of the disturbance output to change over time so that it canexhibit some nonstationary behavior. For example, when d=1 the disturbance will exhibit a randomwalk behavior where its mean varies over time. The specific order of the disturbance is denoted byARIMA(p,d,q) where p, d, and q are as defined above and thus set the properties of the disturbance.These are the two terms, when combined, that can be used to represent most industrial processes.The output of a process which combines the two terms is given by:B(q-1 ) k^CW1) y(t) — ^u(t)-1-^e(t)Ai(c1-1 ) A2((1-1)Vd (2.3)(2.2)5DisturbanceChapter 2: Control Loop Monitoringwhich is made up of the effects of the plant and disturbance. The monitoring to be performed isnot done on the process as it appears in Figure 2.1. This is because the processes to be monitoredwill be under feedback control.2.2 Monitoring Closed Loop ProcessesA goal of this process monitoring is to be able to extract information from the process successfullywhile it is operating under feedback control. This means that all estimation is performed under closedloop operation, which can cause some difficulties when it comes to estimation. These difficultiesmust be addressed to ensure successful monitoring and in future chapters, they will be taken intoaccount as they become a concern.Closed loop operation of a process can be represented by the block diagram in Figure 2.2. Usingthe definitions above, the plant transfer function is:+Gp (q—i)^Yo(t) = B(q 1 ) _k bp u(t)^Ai(q-1 ) cl — 1+ aic1-1+ ... ancru q^(2.4)and the disturbance transfer function is:-I- cicri^...cqq—qGd(q— =1)^YYe(t) — C(c11)^1^(2.5)e(t)^VtIA2(q-1 )^1 + aCcrl + ...a' +d ci-v-dpThe final block to be defined is the controller, it has the form:Gc (q—i)^u(t)^E(c1-1 ) _ eo eicr i c(t)^H(q-1 )^ho h1q-1and is chosen so that the plant output response performs in some desired manner.(2.6)Figure 2.2. Block diagram representation of closed loop process6ProcessOutputReferenceSignalChapter 2: Control Loop MonitoringProcessInput StaticInput-OutputEstimatorProcessOutputNoiseProcessInput^-JProcessOutputProcess IDelayMinimum ^AchievableOutput VarianceProcessOutputProcessInputProcessDelayFigure 2.3. Block diagram of the monitoring scheme2.3 Purpose of MonitoringThe successful operation of a process requires the inclusion of some sort of performancemonitoring and/or fault detection. There are many insights that can be obtained by monitoringthe process, all which will aid in choosing a control scheme for the process or provide a betterunderstanding of the operation of the process. Monitoring will help determine if the process isoperating as expected, which is important for successful control. The detection of changes in theprocess over time will be noticed by monitoring, which will allow preventive measures to be takenbefore the changes become too serious. The cumulation of the insights provided by monitoring willbe used to indicate if and how performance can be improved.There are different properties of a process that can be monitored to gain knowledge of processperformance. When the following three attributes of a process are monitored, namely the timedelay, static characteristics, and output variance, many insights will be provided which will aidin the successful operation of a process. These properties can be used separately or together toserve this purpose. These process properties were chosen because of an interest in monitoring theoperation of processes in regulation mode. Examining the process output variance and determining7Chapter 2: Control Loop Monitoringthe minimum achievable output variance becomes very important for regulation loops. Knowledgeof the time delay and static characteristic will help in determining and interpreting the minimumachievable output variance level. They can also stand alone as monitoring tools. The block diagramof the monitoring scheme appears in Figure 2.3, it showing the way information is extracted from aprocess. It is obtained using only the available input and output signals. The benefits obtained frommonitoring these properties will now be discussed, emphasizing the value of successful monitoring.2.3.1 Time DelayMonitoring the time delay of a process serves several purposes. Most importantly, knowledgeof the true time delay provided by monitoring will ensure that the controller can compensate for theeffect of the delay. Also, monitoring can indicate if changes in the process are required for improvedperformance. Lastly, as changes in the delay can be noted as they occur, the reason for the changecan be determined and corrected if possible.The knowledge of the correct time delay of a process is very important for successful controlof any industrial process. Many processes have long and often varying time delays, which makesproper control difficult. This is because the delay adds a phase contribution to the process frequencyresponse which slows the output response and must be compensated for by a controller. Monitoringthe process to obtain the delay can allow the controller to be designed correctly, so the performanceof the process can be improved.Typically, processes with delay are controlled by one of two approaches. One approach uses atype of controller, such as a Dahlin, a Deadbeat, or a Minimum Variance controller, that we shall calldead-time controllers, whose design requires the knowledge of the exact delay. If the controller usesan incorrect delay value, the resulting process performance is poor or even unstable. Monitoring thetime delay of the process allows confirmation of correct controller design and allows the controllerto be changed on-line if the delay varies, avoiding the associated possible poor performance.The other approach involves the use of Proportional-plus-Integral (PI) and Proportional-plus-Integral-plus-Derivative (PID) controllers. These controllers ignore the exact time delay in theirdesign. The performance of these controllers is therefore not as sensitive to variation in delay as8Chapter 2: Control Loop Monitoringthe dead-time controllers, so they are often used to control processes with unknown time delay.The reason for using these controllers lies in the fact that the performance of dead-time controllersdeteriorates or even becomes unstable if the controller does not use the correct delay. However, theperformance of PI and PID controllers is worse than the dead-time controllers when the time delayis known. By monitoring the process and obtaining the correct time delay, dead-time controllers canbe correctly designed, thus replacing PI and PID controllers, in these situations. The process canthen enjoy the superior performance of the dead-time controllers.There are situations where the control scheme is not able to improve performance to a satisfactorylevel due to limitations imposed by the process itself. The information obtained by monitoring theprocess may reveal that changes in the process are necessary to further improve performance. If, forexample, monitoring the time delay reveals that the time delay changes frequently or is very long,changes in the process may be required to correct this so performance can be improved. It may benecessary to consider moving sensors to shorten delay or perform other measures which affect thelength or variability of the time delay. If the results of these changes make the delay shorter or moreconstant, then the process performance will be improved.A final benefit which arises from monitoring is that changes in the delay will be noticed. Oncea change is noted, the process operation can be investigated to determine the cause of the change.Knowledge about the cause of the change or the degradation in process behavior will provide insightinto the important aspects of the process. If these aspects are watched and kept from changing, theoverall operation of the process can be improved. This discussion has shown a few advantages ofmonitoring the process time delay.2.3.2 Static CharacteristicsA major assumption made when controllers are designed, is that the process to be controlled islinear. Because no process is truly linear, it would be informative to know the degree of nonlinearitypresent in the process. By monitoring the static input-output relationship, the degree of nonlinearitywill be found and can be examined. If it differs a great deal from what was expected then steps canbe taken to take accommodate this difference and improve control performance.9Chapter 2: Control Loop MonitoringAlso, nonlinearities can develop over time, often as a result of deterioration in actuator operation,such as control valves. If this nonlinearity is detected, then maintenance can be performed to returnthe actuator operation to normal. Often there will be nonlinearities present that were not expectedwhen a control strategy was chosen. If this occurs, then the control strategy must be reworked toaccount for the nonlinearity. It must also be realized that achievable performance will be limitedby any nonlinearities.Finally, control performance will deteriorate if sensors or actuators are poorly sized. For example,an undersized actuator will most likely exhibit saturation which in turn will limit the process responsetime to changes in reference or disturbance upsets. Conversely, an oversized actuator may not becapable of making fine adjustments required to keep the output very close to the reference. Similarproblems will be experienced with sensors as well. Poorly chosen sensors may be driven outsidetheir linear region, or may not have a resolution fine enough to detect small, yet undesirable changes.Problems such as these can be determined from estimating the static input-output relationship.2.3.3 Minimum Achievable Output VarianceIt has been shown [7] that the variability of a control loop output can often be lower when theprocess is in open loop or manual operation as opposed to feedback control. Sometimes, this situationresults from poorly maintained control valves and sensors, but often it arises from poorly tuned PIcontrollers which cause much of the variability themselves. Bialkowski [7] estimates that the outputvariance of processes in a pulp and paper mill can be reduced by as much as thirty to fifty percentby improved loop tuning, controller selection and equipment maintenance. In industries such as pulpand paper, where a small variation in the final product is an important goal, this reduced variabilitycan result in major economic savings.It has been determined by Astrom [2] that the minimum variance in an output is achievablethrough the use of a minimum variance controller. The output variance when this minimum variancecontroller is used is the minimum achievable output variance obtainable by feedback control. Amethod of calculating this minimum achievable output variance level for a process under any linear10Chapter 2: Control Loop Monitoringcontrol was determined by Harris [17]. There are several reasons why this minimum achievableoutput variance level should be estimated and monitored.The ratio of output variance to minimum achievable output variance can be used to assess theperformance of the control algorithm which is employed. If this ratio is large, it suggests thatalternative controller parameters be chosen. This retuning of the controller may result in a smallerratio. It may be found that no choice of parameters for the present controller will reduce the outputvariance to the desired level. A new control algorithm may then be required to further reduce theratio. In the case of processes with a large time delay, a dead-time controller should be employed.It is just as likely that the minimum achievable output variance determined is not satisfactory.No feedback controller can reduce the output variance below this level, so some other steps must betaken to reduce the output variance. Modelling of the disturbance can be done so that feedforwardcontrol can be employed and this control scheme will reduce the minimum achievable output variance.Another way of reducing the minimum achievable output variance is by decreasing the time delay.Thus the cause of the delay must be found so that one can be determined if it can be reduced. Theseare a few uses of the minimum achievable output variance once it is known.There are several more process properties which must be estimated in order to determine theminimum achievable output variance. It is only natural that these properties are monitored, becausethe cause of any change in the minimum achievable output variance level can be determined byexamining the estimation of these properties. Simulations will be used to study different changes ina process and the effect these changes have on the calculated minimum achievable output variance.In practice, many of these changes will suggest that maintenance be performed on the process, whichwill help reduce the output variability and operation stoppages.Now that the process properties to be monitored have been chosen, the tools required to determinethem will be examined. In turn, the following three chapters examine the steps necessary to determinethe time delay, static relationship and the minimum achievable output variance.11Chapter 3: Time Delay EstimationChapter 3Time Delay EstimationThere is a great advantage in estimating the time delay of the process and there have been variousalgorithms suggested recently to provide that estimate [14][12]. The chosen method, Fixed ModelVariable Regressor Estimation, is described in this chapter and the steps required to implement thescheme and estimate the delay are listed. Examples are then given to demonstrate the operation ofthe estimator. This is followed by an examination of the estimator to the operating conditions underwhich monitoring is to be performed.3.1 Fixed Model Variable Regression EstimatorThe Fixed Model Variable Regression Estimator (FMVRE) was proposed by Elnaggar [12] as asuitable method for estimating the delay of industrial processes. The FMVRE has advantages overother delay estimation techniques as it is very straightforward to implement and it converges veryquickly. The FMVRE is decoupled from the estimation of any other process parameters and this aidsin its speed of convergence and ease of implementation.In the derivation of the FMVRE, a fixed auxiliary model is used to represent the process. Thisauxiliary model is chosen as a first-order plus delay, since most industrial processes are overdampedand can be approximated by such a model. If the process is of a higher order, it is shown in [12]that this auxiliary model is still suitable to estimate the delay.3.1.1 Derivation of FMVREFor a plant represented by a sampled-data linear model described earlier, the output is:B(q-1 )^k^bincr^y(t) = crk A (q-1 ) u(t, ) = q—^u(t)1 + aiq-1^... aiiq—u (3.7)where the delay is k and is always larger than or equal to one. When this delay is known only to liewithin a range [krnin/kmax], it can be correctly estimated, k, if the cost function:J(1;.) = E{ E 2 (1:)}^ (3.8)12Chapter 3: Time Delay Estimationis minimized, where E is the expectation function. The prediction error E is defined as:E(k) = y(t) — k(t, fc.)^ (3.9)where Sr(t, k) is the output determined for every delay in the range of interest, as given by a chosenauxiliary model. This model is chosen to be first-order, for the reason given above, so its output is:^y(t, 1c) = ay(t — 1) bu(t —^ (3.10)and the cost function can be minimized once a and b are chosen. This minimization is not easy to doyet, so the squared prediction error is evaluated as it will simplify finding the delay that minimizesthe cost function. The squared prediction error is given as:E2 (1-0= [y(t) — ay(t — 1) — bu(t — k)1 2= y2(t) a2y2(t — 1 ) + b2u2(t — fi.) — 2ay(t)y(t — 1)— ^(3.11)2b [y(t)u(t — fc) — ay(t — 1)u(t — k)]Before taking the expectation of this expression, which is required to evaluate the cost function, thedefinitions of autocorrelation and cross-correlation are introduced:ry (r). Ely(t)y(t f r)}ruy (r) = Efu(t )y(t r)}(3.12)= E{y(t)u(t — r)}=^( —7-)Using these definitions, the expectation of the squared prediction error is rewritten as:E{ E.2(i)}^T.:v (0) a2 ry (0) b2 rum — 2ary (1) — 2b [ruy (k) — ar„y (k — 1)]^(3.13)^E0^_ 2bEiUpon examination of the terms that make up Eo and Eli it can be seen that Eo does not depend uponthe delay, whereas E1 does. Thus to minimize 4), E1 must be maximized if b > 0, or minimized ifb < 0. This means that for processes with a positive gain E1 is to be maximized and then for processeswith a negative gain E1 is to be minimized to determine the correct delay. Now E1 is given as:E1^ruy (k) — aruy (k — 1)^ (3.14)13e — TdsG(s) — TS + 1 (3.18)Chapter 3: Time Delay Estimationwhich still requires that the auxiliary pole a is chosen. It was shown in Elnaggar [12] that the choiceof a=1 gives unbiased estimation of delay for low order processes and good delay approximation forhigher order processes. Thus:E1 = ruy (10 — ruy (1: — 1)Ef[y(t) — y(t — 1)]u(t —11.)}=^Ay(t)u(t — 1.'0}ruAy (k)An advantage of choosing a=1 becomes apparent by the appearance of Ay. It means delay estimationwill not be affected by step changes in either the reference signal or disturbances. This is how theprocess time delay is determined and now the FMVRE algorithm will be outlined followed by twoexamples demonstrating its operation.3.1.2 FMVRE AlgorithmThe FMVRE algorithm is very simple to implement and is outlined below. These steps includethe important equations required to implement the algorithm:1. Update Ei for every delay in the interval of interest. A forgetting factor, A, can be used as shownto track time-varying time delay.Ei (t, ki)^AEi (t — 1, ki) + [y(t) — y(t — 1)]u(t — ki) V kiE[k„,in , k„,„„]^(3.16)2. For positive process gain, choose the delay that maximizes El.El (t, k(t)) = max{Ei(t, ki)} V^ki E^kmaxj^(3.17)3. Return to step 1 for the next sampled data.3.1.3 ExamplesThe following two examples demonstrate the operation of FMVRE and confirm its quickconvergence and its ability to quickly track changes in the delay. A first-order process plus delaywas used for these examples:(3.15)14Chapter 3: Time Delay EstimationNow after zero -order hold sampling, with a sample time of T seconds, the process has the followingdiscrete form :kG 1 )= b Cr N + aq-1where:Ta = —e(3.19)b = 1 —^ (3.20),^dK =T—TFor these simulations a white input signal, N(0,1) was used for the estimation of the delay.EXAMPLE 1: The parameters of the process were chosen as follows:Continous Time Discrete Time^(T=lsec)T=2.5 Td=4 a=-0.67 b=0.33 k=5The input and output data were run through the FMVRE and the results are shown in Figure 3.4.It can be seen from Figure 3.4(a) that the estimate converges to the correct delay after about sixsamples and the plot of the E1 function confirms that the delay is definitely five.Figure 3.4: Results of FMVRE time delay estimation with (a) presenting theestimation over time while (b) shows the E1 function at sample one hundred15Chapter 3: Time Delay EstimationEXAMPLE 2: This example illustrates the ability of the FMVRE to adapt to changes in the process delay. For thissimulation, the following parameters were used:Period Continous Time Discrete Time (T=lsec)0- 100 r=2.5 Td =4 a=-0.67 b=0.33 k=5101 -200 T=3 Td=1 a=-0.72 b=0.28 k=2201-300 r=2 Td =3 a=-0.61 b=0.39 k=4and it should be noted that the plant for the first one hundred samples is identical to the plant usedin Example 1. The results of delay estimation are in Figure 3.5 and include four different values ofthe forgetting factor. These results show the change in delay convergence, as it varies from ninetyFigure 3.5: Delay estimation with varying forgetting factor, with (a) A=1, (b) A=0.99, (c) A=0.975, (d) A=0.9516Chapter 3: Time Delay Estimationsamples in Figure 3.5(a) for a forgetting factor of one, down to fifteen samples for a forgetting factorof 0.95 shown in Figure 3.5(d). The reason for this can be understood by examining Figure 3.4(b),which shows the E 1 function at sample one hundred. The difference between the E l function at five,the correct delay at sample one hundred and El at two, the next delay, is quite large and must bemade up before the estimate of the delay will change. This time can be shortened by decreasingthe forgetting factor as it reduces the weighting of previous values of El. This increases the weightof each new sample and will keep the difference between two and five in the El function smaller,so the estimate of the delay will change quicker. Care must be taken so the forgetting factor is notmade so small that the estimation changes too rapidly.3.2 Estimation of Delay in Control LoopsUp until now, the time delay of the plant has been determined while the plant has operatedin open loop and with no disturbances. To be of use for monitoring, the FMVRE must estimatethe delay correctly when the plant is under feedback control and subject to stochastic disturbances.To ensure this is possible, the FMVRE operation is investigated with the addition of a stochasticdisturbance and then with the process under closed loop control. It can then be determined if thedelay estimation becomes biased under any circumstances.3.2.1 Effect of Stochastic DisturbancesBefore treating closed loop operation, the effect of adding a stochastic disturbance will beexamined. The form of such a disturbance has been discussed and a stochastic process was shownin Figure 2.2, where the output was given by:B(q-1 )^k^C(C1-1)y(t) = ^cf-^+ ^ e(t)Ai(cri )^A2(q-1)Vd^ (3.21)= Yu(t) Ye(t)which consists of output due to the white noise as well as output due to the input. The effect of thisextra component on correct delay estimation must be examined.The cross-correlation between the input and output is given by:ruy (t = 1.113, o (t)^ruy jt)^ (3.22)17e(t) Ye(t)Chapter 3: Time Delay Estimationand since the delay is determined by evaluating this cross-correlation, the effect of adding the r use (t)term is important. This cross-correlation can be rewritten once the y e(t) is expressed as a sum ofnoise terms up to time t:ye(t) woe(t) -1- wie(t — 1) + w2e(t — 2) +^(3.23)Then the input-disturbance cross -correlation can be written as:ruye (t) = wor„(t) wirue (t — 1) + w2rue(t — 2) +(3.24)=0This is because in open loop operation, the input signal and noise signal are not correlated. Thismeans that the delay that maximizes the E1 function is unbiased in the presence of process noise ifthe noise is uncorrelated with the process input. This is the case with an open loop process and thisproperty of the FMVRE will be illustrated by an example.EXAMPLE The first-order plant plus delay is used once again, though this time with the addition of a sto-chastic disturbance. Figure 3.6 shows the process used for this example, which has white, N(0,1),uncorrelated sequences as the input and noise signals.FMVRE delay estimation was used to estimate the delay and the results are shown in Figure 3.7.From these results, it can be seen that the delay estimation was not biased by the addition of thestochastic disturbance, as the E1 function has a maximum at four. Since this test was performed ona computer, it is possible to access the signals yo(t) and ye(t) so that their contributions to the E1function can be studied. Figure 3.8(a) shows the E1 function when the signals u(t) and y o(t) are usedand Figure 3.8(b) is the result of using u(t) and ye(t). The E1 function resulting from the signals u(t)u(t)4^Yo(t)^y(t)N(0,1)Figure 3.6: The stochastic process used for the example18Chapter 3: Time Delay EstimationFigure 3.7: Results of FMVRE time delay estimation with (a) presenting theestimation over time while (b) shows the E 1 function at sample one thousandand yo(t) all but matches the E1 function in Figure 3.7(b) while the E1 function for the signals u(t)and yo(t) is essentially zero. This example then highlights the fact that when the process input andnoise signal are uncorrelated, as in open loop operation, the delay estimation is unbiased.3.2.2 Closed Loop EstimationAttention is now turned to the process configuration of Figure 2.2, which has the process underfeedback control. FMVRE delay estimation will be performed on this setup and the possibility ofbiased delay estimation will be investigated. It will be shown that the possibility of biased delaydoes exist and examples to demonstrate its seriousness will be given in Chapter 6.Figure 3.8: E 1 functions for (a) the input and plant output and (b) the input and disturbance output19Chapter 3: Time Delay EstimationIn the majority of cases the process will be assumed to be operating with a constant referencesignal, so any effect of changing the reference signal will be ignored. Now return to the cross-correlation of the input and output which was found previously:ruy (t) = ru3, 0 (t) + ruy. (t)(3.25)ruy. (t) worue(t) wirue(t — 1) + w2rue (t — 2) +By adding feedback to the process, the input and noise signal will no longer be independent and thusthe E1 function will have contribution from ruse (t). The amount of this contribution will depend onthe noise filter and controller and may result in biased delay estimation of the delay. This problemwill be examined further in Chapter 6 to determine its extent and to explore possible solutions.So, now that the FMVRE algorithm has been introduced and its operation demonstrated, attentionwill be turned to the problem of estimating the static input-output relationship and the minimumachievable output variance. Further simulations to investigate the operation of the FMVRE are saveduntil Chapter 6 so that the other tools can be introduced.20Chapter 4: Static CharacteristicsChapter 4Static CharacteristicsIt is possible to identify the static relationship between the input and output of a process while itis operating in feedback control. The shape of this relationship will provide additional understandingof the process operation. The Adaptive Nonlinear Modeller represents the static relationship by usinga table of values and an interpolation formula. After the operation of the estimator is explained,its algorithm will be outlined. Then an example of the Adaptive Nonlinear Modeller operation willbe given to exhibit its functionality. The application of the estimator to the closed loop regulationproblem is briefly discussed and while simulations are saved for a future chapter.4.1 Adaptive Nonlinear ModellerThe Adaptive Nonlinear Modeller (ANM) was proposed by Astrom [3] as a multiple use toolto model the static input-output relationship of a process. The ANM described here uses RecursiveLeast Squares estimation to identify this relationship avoiding some difficulties of the integrationmethod used by Aström. A static input-output relationship can be described by:y = f(u)^ (4.26)where f(•) is simply a function which relates the process input, u, to the process output, y. The ANMrepresents this function with a table of values and an interpolation formula. The table consists ofFigure 4.9: Block diagram of the Adaptive Nonlinear Modeller21Chapter 4: Static Characteristicsdiscrete values of the input over the interval of interest along with the corresponding output for eachdiscrete input. The role of the interpolation formula is to provide output when the input is betweeninterval values. Estimation is required to identify the output values for each discrete input value inthe table, from sampled data. As an aid in the discussion, Figure 4.9 presents the structure of theANM and the following notation is defined:uk — a discrete input value in the tabley — the process output during estimationSr. — the Adaptive Nonlinear Modeller (ANM) outputfk — the estimated value in the table corresponding to input ukf(uk) — the process output for the input ukUsing this notation, the input interval is [ui,u n] and each entry in the table is simply a pair ofnumbers {uk, f(uk)}. The interpolation formula then uses the table values to determine the output forany arbitrary input which lies between two entries in the table, uk +i?_ u uk. Linear interpolationis a very simple and straightforward formula to use and because of this is used by the ANM forinterpolation. The linear interpolation formula which determines the output for an input is givenby:y f(u) = uk+i — u ff(uk) u — ukUk+1 uk^ f (Uk+1)^(4.27)uk-1-1 —The remaining problem is to estimate the values of f(uk) for the table. These entries will be determinedby using Recursive Least Squares (RLS) estimation.In order to use RLS estimation, the model must be placed into a linear regression form. Theoutput for any input in the interval [ui,u n] can be determined using the table of values along withthe interpolation formula and is given as:y = aif(ui) a2f(u2)^an—if(uu—i) anf(uu)^(4.28)All but two of the ai's will be zero, as the input will lie in only one interval [uk, uk+1]. The linearinterpolation formula will determine the value of the two nonzero ai's as:uk_i —ak = ^Uk+1 — 11ku — ukand ak+1 = ^uk-F1 — uk(4.29)22Chapter 4: Static CharacteristicsFor every sampled data, it is not known which ai's are zero, so in general the regression vector isdefined as:(NO = [^a2 ] T (4.30)and the parameter vector is:0 = [ f (ui) f (u2) f (u3) . • • f (11u) ] Twhich allows the following linear regression form:y = 9T1.(t) v(t)The least squares estimate of the parameter vector is defined as:e = [ ft f2 f3^fIl]Twhich can be defined recursively using:P(t — 1)(t) K(t) = 1 +q,T(t)p(t — 1)4)(t)(4.31)(4.32)(4.33)(4.34)O(t) = O(t — 1) K(t) [y(t) — OT (t — 1) ,T.(t)]^(4.35)P(t) = P(t — 1) — K(t)(1) T (OP(t — 1)^ (4.36)where P(t) is the parameter covariance matrix and K(t) is the gain matrix.4.1.1 ANM AlgorithmThe following stages outline the operation of the ANM. They summarize the previous sectionand highlight the important steps and equations required to implement the ANM.1. Determine the input range of interest and break it into discrete intervals.2. For every sampled data, find the interval in which the input value lies.uk+i ] where Uk < u <^ (4.37)23Chapter 4: Static Characteristics3. Use linear interpolation to determine the values for the regression vector.uk+1 — u^u — ukak —^and ak-hi = ^uk+1 — uk uk+1 — uk1:13 = [ 0 ... 0 ak ak+i 0 ... 0F(4.38)4. Determine the estimate of the table values by using RLS estimation.P(t — 1)(D(t)K(t) =1 +4DT(t)P(t — 1)1, (t)O(t) = O(t — 1) + K(t) {y(t) — O T(t — 1)(D(t)]P(t) = P(t — 1) — K(t)4DT (t)P(t — 1)5. Return to step 2 for the next sampled data.The shape of the f(.) function can be seen by plotting the table at any time during estimation.If the number of computations required by the RLS estimation is a concern, it is possible to reducethis number by realizing the fact that for any data, all but two values in the regression vector arezero. This greatly reduces the number of multiplications required by the matrix operations in theequations of (4.39) to determine the estimates.4.1.2 ANM ExampleThis example was used by AstrOm [3] to exhibit the operation of the ANM. It is repeated hereto show the success of the ANM using RLS estimation. The input and output are related by:y = u2 (4.40)which is definitely a nonlinear relationship. The input and output used for estimation are shownin Figure 4.10(a&b), and the resulting parameter estimation for the table output values is shown inFigure 4.10(c&d). Figure 4.10(c) shows that the estimates arrive quickly to near their correct values,while Figure 4.10(d) shows that the table values, represented by zeroes, do a good job approximatingthe squared relationship. Figure 4.11 shows the ability of the ANM to reconstruct the output usingthe table and interpolation formula for the given signal. This confirms the successful operation ofthe ANM with RLS estimation.(4.39)24150^20050 100 250^300^350^400^450^500Time(c)2 4-2-4Estimated Table ValuesInput(d)2015105-6500350 45040030010050Output150^200 250Time(b)5°02520153025205,.! 18 10Estimated Output Values00000160^180 20014020^40^60^80^100^120Time(a)160 180 200140Test Input20^40^60^80^100^120Time(b)105050- - - Estimated Output 'System Output2520Estimated and System OutputChapter 4: Static CharacteristicsFigure 4.10: Results of the ANM estimation showing (a) the input used, (b) the resultingoutput, (c) the table value estimation over time and (d) a plot of the final estimated table valuesFigure 4.11: Results of using (a) the input and the ANM to (b) reconstruct the output25Chapter 4: Static Characteristics4.2 Application to Stochastic Feedback ProcessesThe ANM described has assumed that the input and output are obtained during steady stateoperation of a process. If all transients are not removed from the signals, then the static input-outputrelationship will not be correctly identified. These transients result from the process dynamics whichare prevalent during reference changes and disturbance upsets. Prior to the installation of the ANM,some precaution must be taken to ensure that the input and output signals are in steady state. Lowpass filtering of the signals will remove much of the undesired transient component of the signalsas well as the effect of stochastic disturbances. Also, the ANM could be turned off during setpointchanges and turned back on once the setpoint remains constant for a period of time. These problemsare very important and will be thoroughly examined by simulations in Chapter 6. Also to be examinedis the problem of setting an appropriate cutoff frequency for the low pass filter.The only remaining task is to introduce the method used to determine the minimum achievableoutput variance of a process. Once introduced, the details of all the estimators used for monitoringpurposes will be studied using simulations.26Chapter 5: Output VarianceChapter 5Output VarianceThe variance in the process output becomes important when a constant output is desired. Thecontroller becomes known as a regulator because it attempts to keep the output close to the reference.This chapter begins by determining the minimum variance achievable by a feedback controller. Thenthe method used to find this minimum level is given, for a process in closed loop regulation. Next, thedetails of estimating the minimum achievable output variance are presented, including an explanationof the Laguerre network used to model the process and an explanation of the Recursive ExtendedLeast Squares method used to estimate the parameters of the Laguerre network. Lastly, the chapterends by outlining the algorithm needed to calculate the minimum achievable output variance anddemonstrating its operation by presenting a simulation result highlighting important aspects.5.1 Minimum Achievable Output VarianceThe following derivation of the minimum achievable output variance follows the method usedby Astrtim [2] to determine minimum variance control and uses the notation introduced in equation(2.3). The goal is to find the lowest possible level of the output variance, E{ y 2(t+k)}. It will thenbe possible to estimate this level, so that there is a comparison value for the present output varianceof the process.The derivation begins by writing the process model in terms of present and future values:B(q- 1 )^CN- 1)y(t k) =^u(t) +^e(t + k)Ai(q- 1 )^A2(c1-1)Vdand the last term on the right can be written as:C(q 1 ) ^(-1]A2(q-1)Vd e(t + k) F(ql)e(t + k) + A2((1-1)Vd e(t)(5.41)(5.42)which is now expressed as known and unknown quantities as of time t. The polynomial F(q -1 ) ismonic and of degree k-1, that is one less than the process delay. Meanwhile, the degree of G(q1)27Chapter 5: Output Variancedepends on q, the degree of C(q-I ), and p, the degree of A2(q -1 ), and is the larger of (p+d-1) or (q-1).Dividing (5.42) through by e(t+k) gives a Diophantine equation:C(q 1 )^F(cri)^G(q-1) __k (5.43)^A2(q-1 )Vd^A2(q-1)Vdqor in a more familiar form:C(q1) F(g 1)A2 (q 1)vd G(q1) q—k^ (5.44)Now substituting (5.42) into (5.41) gives:Bq^ G)y(t + k) = A((cr 1i) ) u(t) + F(q-1 )e(t + k) + A2(c(q-1ri)vd e(t)^(5.45)which shows that the input can affect the future output by removing the present noise. It is notpossible however, for the input to remove the future noise from the future output. Yet the value ofthe present noise is unknown and it is not possible to determine the input needed to remove its effect,so this noise must be replaced with something which is known, namely the present output. This isdone by returning to the process model and isolating e(t):B(q-1 )A2N-1 )Vde(t) = ^ u(t k) + A2 (q1) Vd y(t)^(5.46)Ai(cri)c(q—i)^C(q-1)This can be substituted into (5.45) and the following relation can be found:y(t + k) — A iB((q-1-1c1)) [1^C(q-1 )^C(q-1)^G(ci 1) q—k] u(t)^G(ci 1) y(t) + F(q-1 )e(t + k)^(5.47)Now the Diophantine equation can be manipulated to show that:GN-1)^F(cri)A2(cri.)vd1 ^ q k — ^ (5.48)C(q-1 ) C(q-1)which can be substituted into (5.47) to give:y(t k) = B(q-1)F(q-1)A2(g1)Vd u(t)+ G(q1) y(t) + F(q-1 )e(t + k)^(5.49)^AI(q-1 )C(q ')^C(q-1)Then the variance of y(t+k) is given by:E{y2 (t + k)} = E{ [F(q-1 )e(t + k)] 2 +[B(q-1 )F(q-1 )A2(q-1 )Vd u(t) G (q-1 )^]Ai(q-1 )C(q-1 )^C(q-1) y(t) f (5.50)28Chapter 5: Output VarianceThere are terms missing from the above expression which contain the correlation of the future noisewith past inputs and outputs. As the future noise is independent of the past inputs and outputs, theexpectation of these terms is zero. There are now two terms left and it is now possible to removethe present output by the appropriate choice of u(t). The input cannot affect the future noise and sothe future output will always contain these terms. The minimum achievable output variance is thelevel reached when the present output is removed from the future output and is:E{y2 (t + k)} = Ef [F(q-1 )e(t + k)] 2 }^ (5.51)As can be seen from (5.51), the minimum achievable output variance will depend on F(q -1 ) and thevariance of e(t). If the process delay is one, then the variance will simply be:E{y2 (t)} = E{e2 (t)} = cre2^(5.52)Thus, for a single delay, the output variance will equal the disturbance variance. When the delayincreases, so will the minimum achievable output variance. For a process with a delay k, the minimumachievable output variance is found from (5.51) to be:2^2Grum/ = Cre ( 1^f?^•••^fk2-1) (5.53)The coefficients of F(q-I) must be known along with the time delay to determine the minimumachievable output variance.To be able to achieve the minimum achievable output variance, it is required that the input be:A1(q 1 )G(q-1 )u(t) —131(q-1 )F(q-1 ) A2^ 1)Vd y(t )( t (5.54)The application of this control signal will result in the minimum achievable variance in the outputsignal and is known as the minimum variance control law. This controller is what must be used toobtain the minimum achievable output variance and any other controller will not reach this lower limit.However, this controller can require a very active control signal to obtain this minimum varianceperformance and the process actuators must then be able to follow this signal. Also, in order todetermine the minimum variance control signal, knowledge of the disturbance filter is required. If29Chapter 5: Output Variancethis is not known, then a minimum variance control algorithm cannot be used. It may be decided touse another type of controller, in which case its output variance can be compared to the minimumachievable output variance and its performance judged. Lastly, the autocorrelation function of aprocess under minimum variance control will exhibit an interesting and useful characteristic, whichwill depend on the time delay of the process. The function will have significant correlation only atlags less then the time delay. This number of lags will also correspond to the number of coefficientsin the F(q-1 ) polynomial.5.2 Method of Finding Minimum Achievable Output VarianceNow that the minimum achievable output variance level has been determined, what must be doneto find this level. By knowing the minimum achievable output variance, it can be compared withthe present operating output variance to evaluate performance. Based on the previous discussion, theuse of the autocorrelation function seems to be an obvious choice. Calculating this function wouldshow if the minimum achievable output variance has been achieved, but the use of the autocorrelationfunction has several shortcomings. First, while significant correlation only at lags less than the delayindicates minimum variance performance of the controller, correlation beyond this does not give aconcrete indication of how far from minimum achievable output variance the process is operating.It is then difficult to determine whether changes to the present process to improve output varianceare worthwhile. It would also be difficult to determine changes in the minimum achievable outputvariance by just using the autocorrelation function. Though the autocorrelation would change eachtime it was calculated, it would be difficult to determine if the changes were significant. A bettermethod seems to be required.Another possible method of determining the minimum achievable output variance would consistof estimating the disturbance filter and the variance of the noise and solving the Diophantine equation(5.44) to determine the F(q -1 ) coefficients. Once they are determined, the minimum achievable outputvariance can be calculated and compared to the present output variance. The difference between theminimum achievable output variance and operating output variance will indicate if the performance issatisfactory. The major drawback to this method exists due to the fact that the process is operating in30Chapter 5: Output Variancefeedback control while estimation is performed. Parameter identification may not be possible duringfeedback control, especially if the process input is not rich enough. A way of avoiding this problemwas suggested by Harris [17] and consists of estimating the closed loop noise transfer function andfrom this, determining the minimum achievable output variance.Returning to the closed loop diagram of Figure 2.2, the transfer function between the noisesignal and the output is found as:or:y(t) ^Gd(q1) e(t)^1 + Gc (q-1)Gp (cr1) (5.55)y(t)^R(q1)^Ai(cri)c(cri)e(t)^s(cr i)^A2 (q- 1)vd[Ai (cr 1) Gc ( cr i)B(cri) q—kiNow it was shown in the previous section that the minimum achievable output variance is a functionof the noise variance and coefficients of the F(q -1 ) polynomial. If the above transfer function isestimated, then it is possible to find F(q-1 ) and the noise sequence. It will now be shown how to findthe F(q-1 ) polynomial and we start by rearranging the transfer function as:R(q1)^C(g1)^(q-1) S(q 1 )^A2(q 1 )vd (A1(q-1 ) + Ge(q--1 )B(q 1 )q—k)(5.57)and replacing A2 (q_c(q- )mod with the Diophantine equation (5.43) (and dropping the (q -i rs temporarily):RAi + Gcl3q— k(F A2VdGcrk^Al^)FA1^GAisorkAl + G.,13q—k A2Vd[A1 + Gel3q—k]FG,l3q—k^GAlq—k= F + ^ (5.58)Al + G,13q—k A2Vd[A1 + Gcl3q-1 ]F + [A2VdFG,B + GAi] Crk^= A2Vd[A1 + Gcl3q—k]R(q-1)^Gi(cri)qk^ = F(q-1 ) + ^S(q-1 ) S(q-1)The closed loop transfer function can then be expressed as a Diophantine equation where the F(q -1 )polynomial is the one required to determine the minimum achievable output variance. This is notsurprising as it confirms the fact that the controller transfer function can only compensate for the(5.56)31Chapter 5: Output Variancenoise fed back to it and not compensate for the noise placed directly onto the output. After the R(q-1 )and S(q-1 ) polynomials are found then the F(q -1 ) polynomial can be determined by solving:R( 4-1) F (q—i)s(q—i) G/(cri)q—k (5.59)once the size of the F(q-1 ) polynomial is set. This relies on knowledge of the time delay because theunknown coefficients in the F(q-1 ) polynomial number one less than the delay.It will become necessary to determine the coefficients of the F(q 1 ) polynomial from a state-spacerepresentation of the process. This is because the process will be modelled by Laguerre filters whichuse a state-space representation. The transfer function from the noise to the input has the followingstate-space representation:x(t + 1) = Ax(t) be(t)(5.60)y(t) = cx(t) e(t)where for a single input, single output process the following are defined. The matrix A will bea n by n matrix while b and c are vectors of length n. The state vector is x(t) which also has nelements, while the process output is y(t) and the noise sequence is e(t). The value of n will dependon the process represented. Now if the output is written in terms of the noise signal, the equivalentof equation (5.56), the closed loop transfer function, can be obtained. This is performed by writingthe state equation as:qx(t) = Ax(t) be(t)x(t) = (q — A) —l be(t)and this expression can be substituted into the output equation:y(t) = c(q — A) —l be(t)+ e(t)= [1 c(q — A) —l b]e(t)= [1 + cbq-1 cAbq-2cAk-2 bq- k+1 cAk-1 (€1^1 bq-k]e(t)(5.61)(5.62)The terms cb, cAb, cA2b..., are the Markov parameters of the process. Now it has been shown thatthe terms in the output up to the size of the delay are not affected by any control action, so theminimum achievable output variance is defined as:cre2 [1 + (cb) 2 (cAb) 2^(CAk-2b) 2]^(5.63)32Chapter 5: Output Variancewhich can be found from the state-space representation once the matrix A and vectors b and c areknown.5.3 Estimation of Minimum Achievable Output VarianceThe method for finding the minimum achievable output variance has now been outlined, so itis necessary to determine what must be done to estimate it. The method used by Harris [17] andPerrier [23] deals with the estimation of the time-series model previously described. The method to beused here involves modelling the closed loop noise transfer function with a state-space representationof discrete Laguerre filters. Discrete Laguerre filters are used as they avoid the necessity of choosingthe degree of the closed loop polynomials, as required by time-series models. While there areparameters to be chosen so the Laguerre network can be used, estimation is more robust [25] totheir choice than to the choice of polynomial degree in a time-series models. The parameters ofthe state-space representation are estimated using Recursive Extended Least Squares estimation as itis the least complicated method available to estimate process parameters where unknown variablesexist. Unlike Harris, who assumed the time delay is known, the estimate obtained by the FMVREwill be used to determine the number of terms in the minimum achievable output variance. Finally,the minimum achievable output variance is calculated using the method outlined for a state-spacerepresentation of a process.The following sections describe the discrete Laguerre model and how it is formulated for thissituation as well as outlining Recursive Extended Least Squares estimation. The Recursive ExtendedLeast Squares method must be applied to this problem, which requires some of the parameters ofthe state-space representation be estimated. The complete algorithm to find the minimum achievableoutput variance is outlined, including all steps and important equations necessary to find the correctvalue. Lastly, an example of minimum achievable output variance estimation is given.5.3.1 Discrete Laguerre ModelIn discrete time, the Laguerre set is described by the following functions [16]:— V1 a2 ^1N(t) (N — 1)! de-1 [ N-1(1 ac)N-1 —1 < a < 1^(5.64)q=a33y(t)4^Chapter 5: Output VarianceFigure 5.12: Block diagram of Laguerre filter networkwhere N is the order of the function while the parameter a is the function's time scale. As the Laguerrefunctions are orthonormal and complete in L2{0, 00), they can represent the impulse response of anystable sampled linear process with an infinite expansion:00^h(t) = E gN 1N(t)^ (5.65)N=1where gN is the Nth Laguerre gain. The Z-transform of this set of functions results in the followingLaguerre filters:(1 — aq N-1LN (q) = ^ (5.66)q— a q—awhich can then represent any sampled linear process as :00G(q) = E gN LN(q)^ (5.67)N=1The Laguerre filters can be used to approximate a sampled linear process by using N filters whichleads to the network representation in Figure 5.12.5.3.1.1 Laguerre Network State-Space RepresentationThe most convenient form of a Laguerre model to use when the gains are to be estimated, is astate-space representation. The state-space representation of the discrete Laguerre filters is:1(t + 1) = Al(t) bu(t)(5.68)y(t) = c 1(t)34Chapter 5: Output VarianceTo determine the matrix A and the vectors b and c, the block diagram in Figure 5.12 can be examined.A set of equations is required to define this matrix and vectors, some equations which express thefuture states, 1(t+1), in terms of present states and the input while another equation expresses the plantoutput in terms of the present state values and the input. From this block diagram, the followingequations are found:y(t)^gill (t)+g212(t)= [ gl^g2...H-gN1N(t)grd [WO^12 (t)^1N(t)]T(5.69)11 (t)^a2 (t) + V1 - a2 u(t) (5.70)11(t + 1) = au(t)^q - a12(t)^1 - aq 12(t + 1) = al2(t)^li (t) - a li (t + 1)li(t)^q - a (5.71)12(t + 1) = a 12 (t) + (1 - a 2 )11(t) - a/1 - a2 u(t)13(t)^1 - aq 13(t + 1) = a13(t) +12(0 - a12(t + 1)12(t)^q - a (5.72).•. 13(t + 1) = al3(t) + (1 - a2 )12(t) - a(1 - a2 )11(t) a2 V1 - a2 u(t)1N (t)^1- aq= 1N(t + 1) = alN(t)^1N-i(t) - alN_1(t + 1)IN_1(t) q - a1N(t + 1) = a 1N(t) + (1 - a2 )1N_1(t) - a(1 - a2 )1N_2(t) (5.73)(_a)N-2 ^- a2)11 (t )^(_a)N-101 - a2u(t)From equations (5.70-5.73), the matrix A can be formed as:a^0^0^..^0 -(1 - a2 ) a 0^..^0A = -a(1 - a2 )^(1 - a2 )^a^..^0 (5.74)and the vector(_aN-2) ( 1 - a2)^(_aN-3)b is found to be:(1 - a2)^(-aN-4) (1 - a2)^a_^1 - a2-a V1 - a2b = a2 V1 - a2 (5.75)(-a) N-2^— a235By inspection of the A matrix, its elementsA =Chapter 5: Output Varianceare given by:aii = a^i = jaii = 0 i < j^ (5.76)aii = ( —a) i—j-1 (1 — a2 )^i > jand inspection of the b matrix shows its elements are defined by:bi = (From equation (5.57), the c vector—a)i-1^— a2 ,^i = 1 tois:N (5.77)c= [gi^g2^gN (5.78)and with these definitions, the state-space representation of Laguerre filter model is complete.One slight modification is required in the above state-space representation so that it can be usedfor estimation of the minimum achievable output variance. This change arises from the fact that thereis no delay in the closed loop noise transfer function and the noise appears directly on the output.The state-space representation becomes:1(t + 1) = Al(t) be(t)(5.79)y(t) = c 1(t) + e(t)where the matrix A and the vectors b and c are as previously defined. It should be noted that thematrix A and vector b of the Laguerre model are set once the time scale and number of filters arechosen. This allows them to be calculated prior to estimation and used directly in computationswithout any modifications. Now the estimation algorithm used to estimate the elements of the vectorc as well as the e(t) values will be explained. The method chosen to do this is Recursive ExtendedLeast Squares estimation.5.3.2 Recursive Extended Least Squares EstimationThe Recursive Extended Least Squares (RELS) estimator will be explained using a time seriesmodel and then applied to a state-space representation. A process can be described by the followingmodel:A(q-1)y(t) = B(q-1 )n(t) C(q-1 )e(t)^ (5.80)36Chapter 5: Output Variancewhich is just a general case of equation (2.3) where:A(q-1 ) = 1 + aicri^auccuB (q 1) boq—k bmq—k—uiC(q- 1 ) = 1 + cicr l ...crcrrwhich can be expressed as the following linear regression model:(5.81)y(t) = OT(D(t) e(t)^ (5.82)with regressor:(I)(t) = [—y(t — 1)^—y(t — n) u(t — 1) ...u(t — m) e(t — 1)^e(t — r)]Tand parameter vector:(5.83)0 =- [ai^• all b 1 • • • bin ci •^er T^ (5.84)It would be natural to attempt to estimate the parameter vector using least squares estimation, but thismethod cannot be used because the e(i) terms are not known and so the regressor cannot be formed.The least squares method can be modified slightly to include estimation of unknown variables in411(0 from available data. Such methods, which combine the estimation of a parameter vector andunobserved components in the regressor are known as Pseudolinear Regressions (PLR). A morecommon PLR is the extended least squares (ELS) method.5.3.2.1 Recursive Extended Least Squares AlgorithmThe problem with attempting to use least squares estimation in the above situation is theunobserved variables in the regressor. The unobserved variables could be replaced with an estimateof them, given by rearranging (5.82):e(t) = y(t) — OT1.(t)^ (5.85)With this minor change in the regressor, the recursive extended least squares algorithm is given by:cil(t) = [—y(t — 1) ... — y(t — n) u(t — 1) ...(5.86)u(t — m) E(t — 1) ... E(t — r)]T37Chapter 5: Output VarianceP(t — 1)(1)(t)4T(t)P(t — 1) P(t) = P(t — 1 )^1 + (DT(t)P(t — 1)43(0 (5.87)= O(t — 1) + P(t)(13(t)[y(t) — O T (t — 1)41)(t)]^(5.88)E(t) = y(t) — OT (t)(D(t)^ (5.89)The similarity of the RELS algorithm to that of the RLS is the main attraction of this method. Unlikeother algorithms used to determine the coefficients of unobserved variables, this PLR does not requireany extra filtering or checking that the filter is stable.5.3.2.2 ConvergenceThis PLR is not guaranteed to converge to the true parameter vector, 00, and relies on the positiverealness of certain transfer functions in order to converge. For a process exactly described by thefollowing ARMAX model:A0(q 1 )y(t) = Bo (q-1 )u(t)-i-00(q-1 )e(t)^(5.90)a sufficient condition for convergence of 0(t) to 00 is:Re [ C0(eiw)1 ^11 > 0 V co^ (5.91) 2or equally:IC0 (eiw )—li< 1 V w^ (5.92)Thus the filter co(1,,, )^must be strictly positive real (SPR) to ensure the convergence of theparameter vector to the true values. A necessary condition to assure convergence also exists, whichis:Re [Co (eil] > 0 V w^ (5.93)If this condition does not hold, then there exists a A0, Bo, and u(t) so that the probability of 0(t)approaching 00 is zero.38Chapter 5: Output VarianceBefore continuing and detailing how the Laguerre gains are estimated using RELS estimation,the convergence property will be considered. The transfer function to be modelled so that o canbe found, differs from 5.90 in two ways. There is no Bo(q -1 ) polynomial while Ao(q-1 ) and Co(q-1 )are generally complicated polynomials made up from the plant, controller and disturbance filters ofthe process. This means that the convergence of the RELS estimates depends on whether R(q -1 ), andnot Co(q-1 ) is SPR. In the simulation study the influence of this condition will be further examined.By adopting Laguerre filters to model the closed loop noise transfer function, it is requiredthat the Laguerre gains, which make up the c vector in 5.79, be estimated. The following sectionexplains why RELS estimation is required and how it is applied to this situation, when a state-spacerepresentation is used.5.3.2.3 Estimation of Laguerre GainsThe Laguerre model state-space representation of the closed loop noise transfer function is:1(t + 1) = Al(t) be(t)y(t) = CT 1(t) e(t)^ (5.94)the vector c has N elements, the Laguerre gains, which must be estimated. The parameter vector is:0 = [Ci C2 • • • CNand the regressor is:cp(t) = [ h (t) 12 (t)^1N(t) ]Tso the linear regression equation becomes:y(t) = e(t) — O T4.(t)(5.95)(5.96)(5.97)Now the regressor elements are solved from the state equation:(D(t) = Al(t — 1) + be(t — 1)^ (5.98)39Chapter 5: Output Variancewhich involves the unknown noise sequence, so RELS is required to estimate the elements of thec vector. The above linear regression is then rewritten to give an estimate, E(t), of each elementin the noise sequence:(t) = y(t) — 9T(t)(1.(t)^ (5.99)and then the regressor is defined as:Ib(t) = 1(t) =A1(t — 1) + bE(t — 1)^ (5.100)The RELS algorithm can be carried out with these definitions and so the unknown parameter vectorcontaining the Laguerre gains can be identified.5.3.3 Minimum Achievable Output Variance Estimation AlgorithmNow that all the tools necessary to estimate the minimum achievable output variance have beendescribed, the complete algorithm for its estimation will be given. Following this algorithm willobtain an estimate of the minimum achievable output variance.1. Prior to estimation, the number of filters, N, and the time scale, a, should be chosen. The matrixA and vector b can then be determined.aii = a i = jA= aij = 0ail =i < j— a2 )^i > jbi = ( — a2 , i = 1 to^N2. Next, an estimate of the c vector is obtained along with an estimate of the noise value.1(t) = Al(t — 1) + bE(t — 1)P(t) = P(t — 1)^P(t — 1)1(t)1T (t)P(t — 1)1 + IT(t)P(t — 1)1(t)e(t) = e(t — 1) +P(t)1(t)[y(t) — e T (t — 1)1(t)]E(t) = y(t) — j(t)1(t)3. The estimate of the plant delay is found.Ei(t,ki) = AE1(t — 1,k1)^[y(t) — y(t — 1)111.(t — ki) V kie[kmi„,El (t, k(t)) = max{E(t, ki)} V^ki E [k111111, knia,](5.101)(5.102)(5.103)40e(t)y(t)zeroreferenceChapter 5: Output Variance4. Finally, the estimate of the minimum achievable output variance is calculated using the estimateof the plant time delay, k, the estimate of the Laguerre gains, c, and the estimate of the noisesequence, E(t).[ 1,_1.e,..2w ,_ k2_ , + E (eAi_i b) 2i= 1(5.104) 5. Return to step 2 for the next sample data.Step four does not have to be performed for every new sampled data and to save calculations, theAb, A2 b, up to Al"• b can be computed prior to estimation.5.3.4 ExampleThis example will show the ability of this algorithm to estimate the minimum achievable outputvariance, as well as highlight interesting aspects that result during estimation. The process used forthis example is shown in Figure 5.13 and consists of a first-order plant and a ARMA(2,1) disturbance.The time delay of the plant is four samples and the plant is under Dahlin control. Along with theprocess, some of its properties are also shown in Figure 5.13. The transfer function to be estimatedis shown along with the coefficients of the F(q-1 ) polynomial. There are four coefficients in thispolynomial which corresponds to the delay and are used with the variance of the noise to determinethe minimum achievable output variance shown. Comparing the output variance to the minimumcre2 = 0.81cr2Y = 1.2116a2 = 0.9274myY(t)^1 -q-1 +0.31q-2 -0.03q-3 -0.05g-4 +0.25q-5 -0.03g-6e(t) 1 - 1.24g-1 + 0.37q-2F(g-1 ) = 1 + 0.24q-1 + 0.2376q-2 + 0.1758g-3Figure 5.13: Process used in the example41Chapter 5: Output VarianceLaguerre Gains EstimationFigure 5.14: The Laguerre gains estimated during minimum achievable output variance estimationachievable output variance, it is clear that the Dahlin controller keeps the output quite close to theminimum value.The results of one simulation using one thousand data points and using a Laguerre network withsix filters and a time scale of 0.3, are shown in Figures 5.14-5.16. This simulation used a constantreference signal, which will be the situation in a regulation problem, and in this case the referenceFigure 5.15: The (a) estimate of the F(q-1) polynomial and (b) the estimated noise autocorrelation42Chapter 5: Output VarianceFigure 5.16: The (a) FMVRE time delay estimation and the (b) estimatedminimum achievable output variance and process output varianceused was zero. The estimation of the Laguerre gains appears in Figure 5.14 and it is apparent that thegains level out after roughly four hundred points. However, these gains do not give any indicationof how well the transfer function has been modelled by the Laguerre network. To see this, equation5.63 is used to calculate the F(q-1 ) coefficients, which appear in Figure 5.15(a). The close agreementsuggests the Laguerre network was successful in modelling the transfer function. The autocorrelationof the estimated noise sequence appears to approximate white noise, as shown in Figure 5.15(b) asit has an impulse shape. This is what is desired as the sequence being estimated is white noisewhich has an autocorrelation function that has no correlation at all nonzero lags. The true delaywas correctly estimated after roughly fifty samples and the estimated minimum achievable outputvariance is compared to the present output variance in Figure 5.16(b). As the closed loop functionhas an infinite impulse response, it should be realized that the calculated output variance is only anapproximation of the true output variance. The estimated value, approaches o-,2„v , which was tobe expected based on the estimation of the noise sequence and F(q-1 ) coefficients. As mentioned,this simulation was performed using six Laguerre filters with a time scale of 0.3, this information isalways shown on the plot of the estimated F(q -1 ) coefficients, Figure 5.15(a). One simulation maynot, however, give a complete sense of how well the estimation is performing and conclusions drawnbased on a single simulation could be erroneous.43Chapter 5: Output VarianceFigure 5.17: Distribution of the results of five hundred simulationsTo arrive at a better sense of how well the estimation performs, this simulation was repeated fivehundred times using different noise sequences to obtain a distribution of the results. These resultsare presented in Figure 5.17 as a histogram of the estimated minimum achievable output variance toits calculated value. From this histogram it can seen that roughly twelve percent of the 500 testsresulted in an estimated minimum achievable output variance that was between 0.99 and 1.01 of thecalculated value. Also, the time delay percentage is listed in the upper right corner of the plot. Thevalue listed corresponds to the percentage of the five hundred tests which had the delay correctlyestimated. In this case, the delay was correctly estimated in all the simulations. The shape of thecurve gives an indication of the distribution or variance of the values estimated. As the mean valuewas 0.9599, this indicates that the results obtained in Figure 5.14-5.16 can be considered typical.Also, the mean value of 0.9599 gives a measure or indicator of how well this choice of filter numberand time scale performed.Now that all of the monitoring tools have been discussed, a comprehensive simulation studywill be presented which will examine the operation of the various estimators. The study is used toexamine the performance of the estimators and to determine any limitations of the methods.44Chapter 6: Simulation StudyChapter 6Simulation StudyThe following chapter presents various simulations to determine the performance of the mon-itoring tools. The simulations shown are by no means inclusive of all the simulations performed,but they do highlight the important trends encountered. Because the results of simulations dependon the process used, the results obtained are not conclusive. No recommendations should be applieduniversally to all processes as the conclusions are based on the specific process examined. Theconclusions derived from simulations provide a good starting approach to follow when applying themonitoring to new processes.The first study investigates the possibility of biased delay estimation occurring. Next, the effectof feedback and noise on the estimation of the static input-output relation by the ANM is examined.This is followed by estimation of two typical nonlinearities, an actuator saturation and a deadzone.The operation of the method used to determine the minimum achievable output variance is thenstudied. Various process properties are modified in order to see how they manifest themselves inthe estimation. Finally, the choice of Laguerre network parameters is examined including how theyshould be modified based on the characteristics of the closed loop noise transfer function. The resultsof these simulations can be used when the methods are applied to data from industrial processes.6.1 Time Delay EstimationThe FMVRE generally provides an estimate of the true delay in a process. However, there aresituations where the estimate of the delay is not correct. This section presents the three differentsituations found in which a biased delay estimate occurred. This knowledge will enable the user ofthe FMVRE to determine the validity of a delay estimate. Once examples of different cases havebeen introduced, each will be examined individually. The biased estimation encountered usually fallsinto one of the following three examples, which are an incorrect delay estimate of one, an incorrectdelay estimate but not at one, or an incorrect delay estimate due to insufficient excitation.45Chapter 6: Simulation StudyFigure 6.18: Process which results in correct delay estimate6.1.1 Unbiased EstimationThis example of unbiased estimation is given in order to use as a comparison in the futurediscussion of biased estimation. By comparing aspects of both the biased and unbiased estimation,it is often possible to obtain insight into the causes of the bias. The process appears in Figure 6.18and has a structure similar to the ones to be used which result in a biased time delay estimation. Ithas a time delay of seven and is controlled by a deadbeat controller. The resulting E1 function isshown in Figure 6.19 and the true delay of seven was estimated. Now, processes will be presentedthat result in biased delay estimation and this process will be referred to during the examination ofthe causes of biased estimation.Figure 6.19: E1 function from FMVRE estimation resulting in correct delay of seven46e(t)N(0,1)1 — 0.5q-1 + 0.6q-2 ^1 — 0.8q-1Ye(t)y(t)Chapter 6: Simulation StudyFigure 6.20: Process which results in a biased delay estimate6.1.2 Estimated Delay of One Sampling IntervalThe first situation to be examined results in a biased delay estimate of one sampling interval.The process under consideration appears in Figure 6.20 and it consists of the same first-order plantwith a large time delay of seven. However, now the process uses a Dahlin controller which movesthe closed loop pole to 0.4 and a second-order, ARMA(2,1) noise disturbance. The results of theFMVRE estimation appear in Figure 6.21 and show that the delay was incorrectly estimated. Thisexample illustrates that it is possible for an incorrect time delay estimation to occur.6.1.3 Arbitrary BiasThe next situation considered uses the process in Figure 6.22, which also has a first-order plantmodel, but has a time delay of two. It also uses a Dahlin controller, identical to the one usedFigure 6.21: E 1 function from FMVRE estimation resulting in a bias at one47e(t)N(0,1)1 — 0.5q-1 +0.6g-21 — 0.8g-1ye(t)y(t)Chapter 6: Simulation StudyFigure 6.22: Process which results in a biased delay estimatepreviously. The results of the FMVRE delay estimation appear in Figure 6.23. On this occasion,the time delay is also incorrectly estimated, but as four. This process shows that it is possible tohave an arbitrary bias of the time delay.6.1.4 Insufficient ExcitationThe final process used appears in Figure 6.24 and its major difference from the previous processesconsidered is its use of a PID controller. The time delay of the plant is three and the noise filteris unchanged from the last example. The results of the FMVRE estimation, shown in Figure 6.25,indicate that the delay is one. This is not the correct delay and the peculiar shape of the E1 functionindicates that the cause of this bias is different than in the case of the previous examples.Now that it has been demonstrated that it is possible to have a biased time delay estimate, thereason for its occurrence and possible solutions will be examined. It is hoped that the cause of theFigure 6.23: Ei function from FMVRE estimation resulting in a bias at four48Chapter 6: Simulation StudyFigure 6.24: PID controlled process which results in a biased delay estimatebias can be eliminated, but failing that, it should be possible to be able to indicate if the estimateddelay is biased.6.2 Causes and Solutions for Incorrect Delay Estimation6.2.1 Biased DelayThe first two examples of biased delay share a similar cause, however they affect estimation ina different manner. According to Elnaggar [12], a sufficient condition exists so that delay estimationis unbiased. This condition states that in order to ensure correct delay estimation in a closed loopprocess, the autocorrelation function of the output of the disturbance filter, rye (r) must be zerofor all lags over the interval [kmin+kft,kmax ], where kft, is the delay in the controller. If r y.(7),Figure 6.25: E 1 function from FMVRE estimation resulting in a bias at one49Chapter 6: Simulation StudyFigure 6.26: The disturbance output autocorrelation from the process (a) resulting ina delay bias, Figure 6.20, and (b) resulting in a correct delay estimate, Figure 6.18Figure 6.26(a), is examined for the process in Figure 6.20, it is evident that the correlation at lagsover the interval [0,10] is not zero. However, if ry.(r), Figure 6.26(b), is examined for the processin Figure 6.18, it can be seen that it also has correlation at lags over the interval [0,10]. Thisemphasizes that the condition on ry,(7) is only a sufficient condition and does not fully determineif a bias occurs. Therefore, something else plays a role in creating the bias. This other factorcan be determined by examining the cross-correlation between the input and the disturbance filteroutput, ru3, e (r), which is obtained by the FMVRE when the signals u(t) and y e(t) are used. TheFigure 6.27: E 1 functions for (a) the input and plant output and(b) the input and disturbance output resulting in a bias at one50654320El Function for Plant OutputE1 Function for Disturbance OutputEl (u, yo)___J-1-I_E (u ye)1^,1.50.5Chapter 6: Simulation Study2^4^6^8^10^12^ 2^4^6^ 10^12Delay Delay(a) (b)Figure 6.28: E1 functions for (a) the input and plant output and(b) the input and disturbance output resulting in a correct estimatecalculated Ei(u,y e) function, Figure 6.27(a), shows a large bias at one, when compared to the peakin the Ei(u,y0) function corresponding to the correct delay, Figure 6.27(b), shows the value at one islarger than the correct value. When the two functions are added together, a bias occurs at one andan incorrect delay is suggested. Returning to the process in Figure 6.18, its Ei(u,ye) and Ei(u,y0)functions appear in Figure 6.28 and show that the peak at seven of the E1(u,y 0) function is largerthan the peak at one of the E1(u,ye) function and so the resulting E1 (u,y) function did not have abias. So the bias in the delay occurs when the E1(u,ye) function, created by the plant, is corruptedby the noise Ei(u,ye) function.From the point of view of monitoring, it is not reasonable to expect to make changes to theprocess in order to improve conditions for better identification. However, with a little consideration,it may be possible to remove the effect of this bias from estimation. Upon examining Figure 6.21,the largest value after the bias at one, is the true delay at seven. Would it be possible to change theinterval [kff,in ,kn,.] to avoid the bias at one? The answer is yes if there is a time delay in the processbeyond that resulting from the sampling process. In this case, though the delay may vary, it wouldnot be possible for it to disappear, causing the true delay to correspond to one. Thus by using alower bound on the possible delay interval above one, this biased delay problem can be avoided.51Chapter 6: Simulation StudyFigure 6.29: E1 functions for (a) the input and plant output and(b) the input and disturbance output resulting in a bias at fourNext, the bias that resulted in the process shown in Figure 6.22, is considered, and it should benoted that it only differs from the process used previously in that the plant time delay was changed totwo. The FMVRE estimator used a range of [2,10] as the possible delay and yet the delay was stillbiased. The reason for this can be understood by examining the r uy.(7) cross-correlation function,as obtained by FMVRE estimation. Figure 6.29(a) shows the Ei(u,y e) function obtained for thisprocess. The value of the function at two is negative, so when this function is added to Ei(u,Y0),Figure 6.29(b), the correct delay value is cancelled. The maximum at four, as found in Figure 6.23,results from this cancelling effect.It is possible to correct this bias, using an approach similar to Clarke [10], if the input and outputsignals are modified using a filter Gf(ci l ), so that the resulting E1(uf,yef) function does not affect theEi(uf,yof) function and cause a bias. For processes considered here, the filtering is performed as:GT 1 (q--1) y(t ) Gp (q-1) GT1 (q-1) u(t) Gf-t (q--1)Gd (cil e(t)(6.105)yf(t) = Gp (crl ) uf(t) yef(t)which shows how the filtered signals are formed. The optimum filter to use is the disturbance filterof the process, which in this case would be:Gf (cri)^1 — 0.5cf- 1 0.6c 21 - 0.8q 1 (6.106)and shall be used here for illustration purposes. The two functions, E1(uf,y ef) and E1(uf,y0f) werecalculated and appear in Figure 6.30. They show that the negative value at two has been removed in520.3E, Function for Disturbance Output2.5E1 Function for Plant Output0.E1 (uf, yef) El 01 ft Yoôj0. 1 .5-0-0.402^3^4^5^6^7Delay(a)8^9^10^11 4^5^6^7Delay(b)8^9^10^11Chapter 6: Simulation StudyFigure 6.30: El functions for (a) the filtered input and filtered plant output and (b)the filtered input and filtered disturbance output resulting in correct delay estimateEi(uf,y,f), while El (uf,yof) still has a maximum at two, the correct value. When the two functions areadded together, the resulting E1(uf,yf) function, Figure 6.31, shows correct delay estimation. Whileit may not be possible to find the disturbance transfer function, Elnaggar [13] showed that the exacttransfer function is not required to remove the biased delay.6.2.2 Insufficient ExcitationThe final incorrect delay estimation occurred to the process in Figure 6.24 and the resultingFMVRE estimation suggested incorrectly that the plant delay was one. This is definitely wrong andthe cause of this error will be explained. It will be shown that there is nothing that can be done toFigure 6.31: Et function from FMVRE estimation using filtered signals, resulting in a correct delay estimate532 10^124^6Delay(b)2 10^124^6Delay(a)-0.05 •-0.1-0.150.2Et Function for Plant Output0.3E1 (u, yo) •0.250.20.150.10.050El Function for Disturbance Output0.3E (u ye )1^,0.25 •00.15o.0.05Chapter 6: Simulation StudyFigure 6.32: E t functions for (a) the input and plant output and(b) the input and disturbance output resulting in a bias at onecorrect the problem, however it will be shown that it should be possible to detect that the problem ispresent, thus being able to warn that the biased delay estimation has occurred.The result which suggests that the cause of this bias differs from the previous cases is shownin Figure 6.32. The FMVRE estimation using the plant output and disturbance output show that thecontribution from the E1 (u,ye) function is not the cause of the bias, as the FMVRE obtain erroneousresults from using the plant output, E1 (u,y 0). This had not occurred in the previous cases and suggestsa different cause of the bias. This led to the investigation of the excitation requirements of the inputsignal used by the FMVRE in order to obtain correct delay estimation. These requirements werefound by Elnaggar [12] and it was shown that the input must contain frequencies in the range:71- V2T,Tini11 > 7 > Ts k• Kmax kim11 )^ (6.107)to obtain correct delay estimation, where T5 is the sampling period and 7„,i, is the fastest timeconstant of the process. It should be noted that the lower bound on the frequency was determined forminimum phase stable systems with a pole excess of one, which is the case for first-order plants, thetype under consideration. The upper bound on the frequency ensures that there is only one maximumpeak in the E1 function over the range [k niin ,kmax ]. It is the lower bound however, which is of moreinterest, as this is what causes the problems experienced. The lower bound is set to ensure that thephase difference between the input and output is large enough to ensure correct delay estimation.54Chapter 6: Simulation StudyFigure 6.33: Power spectral density function for the input signal from the PID controllerReturning to the process under consideration, the lower bound on the input frequency contentis (assuming Ts=1):fe = (7r OTs r„ii,i ) --1 = (71- V2(1)(4.4823)) -1 = 0.1063 Hz (6.108)Thus for unbiased delay estimation to occur, the input must contain frequency components above thislevel. The power spectral density function of the input in this process appears in Figure 6.33. Notsurprisingly, the frequency content of the input all but disappears above 0.1 Hz. This confirms thatthe cause of the bias was indeed an input signal with too low a frequency content. For comparisonpurposes, the power spectral density function of the input for the process in Figure 6.18, where correctdelay estimation occurred, appears in Figure 6.34. It shows that the input has frequency content overa wide range, which allows unbiased estimation.Figure 6.34: Power spectral density function for the input signal from the Deadbeat controller554I I I I I I I I-I ^ 4- I- I -4---- 4--- --1- - ---I--I I I I I I I I-I^ 4. I- I -I 4. ... I-I I I I I I I II I I I I I I II,I,I,I,I,IiIIII, ,,^,,,i^,, ,III1IChapter 6: Simulation StudyWhat caused the removal of the higher frequencies from the input signal used in the process ofFigure 6.33? As the process is being driven by white noise only, the input should contain componentsat all frequencies. The answer is found upon examination of the controller used. The magnituderesponse of the controller being used, Figure 6.35 (again assuming T5=1), shows attenuation of thehigher frequencies in the input signal. This is common when a PID or PI controller is used on aprocess with large delay, as stability can only be maintained by choosing controller parameters thatignore fast changes and makes only slow adjustments.From the point of view of monitoring, little can be done to cure this problem. The controller iscausing the problem and changes to it can not be justified only to improve monitoring. However, theproblem itself suggests that consideration be given to the fact that the controller be changed. Thoughthe problem cannot be fixed, its presence can be detected. It becomes apparent after examination ofthe shape of the El function, because if it has a sinusoidal shape rather than a definite maximum,then this suggests that the input signal does not contain sufficient excitation. Also, an empirical andtherefore rough limit on the lower bound can be formulated. Assuming that the sampling occurs ata reasonable rate, that being around 4 or 5 times smaller than the process time constant, and that thesampling time is known, the lower bound is found to be:;:se (7r V.TTO 1fi^ (6.109)Thus for any Ts , a lower bound on the frequency content can be quickly determined to indicate thefrequency requirements on the process input.PID Controller Magnitude Response Function161412106420.05^0.1^0.15 0.2^0.25 0.3^0.35^0.4^0.45^05Frequency (Hz)Figure 6.35: Magnitude response of the PID controller56Chapter 6: Simulation StudyIt has been shown that the FMVRE delay estimator is not infallible. In general, when dealingwith this estimator, it is not prudent to take the suggested result without examining the Et function.In cases where this function has a large and definite maximum value, then the delay can be treatedwith some confidence. However, when the El function has two or more peaks close together, as isthe case in the above examples, then the estimate cannot be used with complete confidence. So likeall estimation, the validity of the result should be considered before it is used with abandon.6.3 ANM Estimation of Stochastic Feedback ProcessesThe ability of the ANM to work when the signals are provided by a process under feedbackcontrol and/or if the signals are noisy must be investigated. It will then be possible to conclude ifthe ANM is suitable to be used for process loop monitoring. The purpose of these simulations is totest the ability of the ANM to correctly estimate the static relationship under these conditions and totest some procedures which attempt to lessen the adverse effect of these conditions.The process to be simulated, Figure 6.36, consists of a first-order plant with a gain of five anda delay of two as well as a disturbance filter which colors the white noise sequence that drives it.A Dahlin controller is used for control, however it was poorly designed using a delay of three, toget a more active control signal. In order to get a yardstick with which to compare the followingsimulation results, an input signal, Figure 6.37, was placed through a pure gain of five only, andthe output generated was used by the ANM to estimate the static gain. The results obtained fromFigure 6.36: Process used to test ANM57Chapter 6: Simulation StudyFigure 6.37: Signal used to estimate the static input-output relationshipthe estimation appear in Table 6.1 beside the row marked Test A. The results are quite close to theexpected value and give a sense of how close the following simulations can come to the correctvalue, using the given input.The next simulation used the process in Figure 6.36, but with the disturbance filter disconnected.This was done to investigate the effect of the plant dynamics, introduced by feedback, when thereference signal changes. The reference signal used was chosen so that the steady state control signalwould be the same as the signal in Figure 6.37. The input and output signals used by the ANM formodelling of the static gain appear in Figure 6.38 and the final table values estimated appear afterTest B in Table 6.1. Comparing these results with Test A, indicates that the brief periods whereplant dynamics were present did not greatly affect the static gain estimation. It was mentionedpreviously that Astrom suggested turning estimation off after a reference change occurred, to letthe plant dynamics die out. This was done, by turning estimation off for ten samples after thereference change, thus removing the bulk of the dynamics. The results appear in Table 6.1, Test C.As expected, the estimated values improved and practically returned to the results obtained in TestInput 2 1.5 1 0.5 0 -0.5 -1 -1.5 -2Output 10 7.5 5 2.5 0 -2.5 -5 -7.5 -10Test A 9.88 7.43 4.95 2.50 -0.01 -2.47 -4.88 -7.61 -9.38Test B 9.80 7.49 5.05 2.53 -0.13 -2.43 -4.82 -7.49 -9.36Test C 9.85 7.43 4.95 2.50 -0.01 -2.47 -4.86 -7.61 -9.33Table 6.1: ANM results for estimation performed under closed loop operation58Chapter 6: Simulation StudyFigure 6.38: The (a) process input and (b) process output from closed loop operationA. While turning estimation off during reference changes improves estimation, the deterioration inthe first place, even for quite an active signal as in this case, was not that large. This suggests thatturning the estimator off may not be necessary.For the remaining tests, the noise in Figure 6.36 was added, and the estimator used the resultingsignal. This addition will change the static input-output estimation in two ways, as it will add noiseto the to the output signal, so the output is comprised of both the input and the noise. Secondly, thenoise on the output will cause the output to vary from the reference, thus causing the control to beactive, introducing dynamics into the output. The signals which are obtained from the process areshown in Figure 6.39 and they were used directly for estimation of the static gain. The results ofFigure 6.39: The (a) process input and (b) process output from stochastic closed loop operation59Chapter 6: Simulation StudyFigure 6.40: The (a) filtered process input and (b) filtered process output used by the ANMthe estimation, Test D in Table 6.2, show that the addition of the noise does indeed adversely affectestimation. This means something must be done so that the estimation can be improved.The signals of Figure 6.39 were filtered, using a low pass Butterworth filter, in an effort toremove the effects of the noise on estimation. Three simulations were done and their results are inTable 6.1, Test E, F and G. The filter used for Test E had a cutoff frequency of 0.3 of the samplingfrequency while Test F had the frequency reduced to 0.1 and Test G used a cutoff of 0.05 of thesampling frequency. The signals which resulted from the filtering are typical of those in Figure 6.40which are filtered by a filter with a cutoff of 0.05. Below a cutoff frequency of 0.05, the results ofthe estimation show little improvement. This is because the noise has low frequency componentsthat cannot be filtered without destroying the useful part of the input and output signals.Input 2 1.5 1 0.5 0 -0.5 -1 -1.5 -2Output 10 7.5 5 2.5 0 -2.5 -5 -7.5 -10Test D 9.96 7.59 4.63 3.22 -0.41 -1.96 -5.26 -7.63 -8.24Test E 9.66 7.60 4.56 3.33 -0.28 -2.28 -4.79 -7.85 -8.46Test F 9.76 7.62 4.80 2.98 -0.03 -2.47 -4.74 -7.80 -8.84Test G 9.81 7.66 4.90 2.84 -0.04 -2.42 -4.78 -7.72 -9.09Table 6.2 ANM results for stochastic process using various low pass filters60Chapter 6: Simulation StudyThese simulations show that it is possible to estimate the static input-output relationship of aprocess under feedback control. They show that it is helpful to turn the estimation off during referencechanges and for stochastic processes, the use of a low pass filter is advised.6.4 Estimating NonlinearitiesAs it has been shown that it is possible to estimate the static input output relationship of aprocess under feedback control, the ability of the ANM to estimate nonlinear relationships will beinvestigated. By completing these simulations, it will be possible to ensure the ANM's ability tocorrectly estimate nonlinearities and further determine the signal conditioning required by the ANMfor improved results.6.4.1 SaturationThe first nonlinear relationship examined is a saturation curve, which is often brought about byan undersized actuator used to drive a process. A typical saturation curve appears in Figure 6.41,which shows that once the input goes above three or below minus three, the output levels off, nolonger able to follow the input. The slope of the saturation region and the location of the start of thesaturation both depend on individual actuators being used. For this simulation, the saturation curveshown will be used with the process in Figure 6.42 which is a first-order plant under Dahlin control.As previous simulations examined the effect of noise on the estimation of the static relationship,SaturationFigure 6.41: Saturation nonlinearity to be estimated61Chapter 6: Simulation StudyFigure 6.42: Process with saturation nonlinearityno noise was used here, so only the ability of the ANM to estimate these nonlinear relations underfeedback conditions could be judged. In order to get an idea how the ANM is able to estimatethe saturation curve, the input shown in Figure 6.43(a) was placed through the saturation with theresulting output shown in Figure 6.43(b). The results of the estimation appear in Table 6.3, Test A.The estimated values vary only slightly from the expected results and confirm that the ANM is ableto identify a saturation curve.Next, the process in Figure 6.42 was simulated, using a reference signal like Figure 6.43(a), andthe signals u(t) and y(t) were used to estimate the input-output static relationship of the process.These signals appear in Figure 6.44 and exhibit an interesting characteristic which results from thesaturation curve and the integrating action of the controller, which attempts to drive the output to thereference value. The input signal ramps up to very large values, attempting to erase the error betweenthe output and reference. Despite this control signal, which could be improved by introducing ananti-windup mechanism in the controller scheme, the ANM estimates the saturation curve quite well.The results, listed beside Test B in Table 6.3, show very good estimation between three and minusFigure 6.43: The (a) process input and (b) process output from the saturation only62Chapter 6: Simulation StudyInput 5 4 3 2 1 0 -1 -2 -3 -4 -5Output 3 3 3 2 1 0 -1 -2 -3 -3 -3Test A 2.97 2.99 2.97 2.00 0.98 0.00 -1.00 -1.99 -2.97 -3.02 -2.91Test B 2.48 2.30 2.89 1.92 1.21 0.05 -0.96 -1.92 -2.97 -2.77 -2.42Test C 3.51 2.14 2.90 2.00 0.99 0.01 -0.99 -2.00 -2.97 -2.77 -2.91Test D 3.24 2.84 2.96 2.00 0.99 -0.00 -0.99 -2.00 -2.97 -2.97 -2.96Table 6.3 ANM results for estimation of saturationthree. While the estimation is not exact, this is to be expected as the control signal spends only shortperiods in the range three to five and minus three to minus five.The controller was changed to a simple gain in an effort to see how the estimation performedwithout the windup problem. The gain was chosen so that the process was stable and responded fairlyquickly, and as a result was chosen to be three. The input and output signals used for estimation areshown in Figure 6.45. The table values estimated by the ANM are listed as Test C in Table 6.3 andshow a very close estimation of the saturation curve.The final test consisted of simply turning estimation off for ten samples after a reference change.This causes the dynamics in the signals resulting from the plant to be removed. The results of theestimation appears beside Test D in Table 6.3 and shows the improved estimation as a result oftemporarily turning off the estimation.Input^ Output50 4 •30-2-3-454030201 0 •0.-10 •-20 •-30 •-40 •^50 ^ -^0 200^400^600^800 1000 1200Time(a)1400 1600 1800 2000 200^400^600^800^1000 1200 1400 1600 1800 2000Time(b)Figure 6.44: The (a) process input and (b) process output from closed loop operation with Dahlin controller63Chapter 6: Simulation StudyFigure 6.45: The (a) process input and (b) process output from closed loop operation with gain controller6.4.2 DeadzoneThe next nonlinear relationship to be examined is known as a deadzone. A deadzone ischaracterized by a range of input where no change in output occurs. Figure 6.46 displays a typicaldeadzone, which can be brought about by a number of reasons and often represents a malfunction inthe operation of an actuator. For this particular deadzone, the output does not change over the interval[0.5,-0.5] and returns to a linear slope outside of this interval. The size of the deadzone and the slopeof the curve outside the deadzone will depend on the individual situation. For these simulations,the process in Figure 6.42 was used once again, however the saturation curve was replaced by thedeadzone in Figure 6.46.Figure 6.46: Deadzone nonlineality to be estimated64Chapter 6: Simulation StudyFigure 6.47: The (a) process input and (b) process output from the deadzone onlyThe simulations performed were similar to those done for the saturation nonlinearity andstart by passing the signal in Figure 6.47(a) through the deadzone which resulted in the signal inFigure 6.47(b). The ANM estimator used the signals in Figure 6.47 and attempted to estimate thedeadzone nonlinearity. The results of the estimation, which appear in Table 6.4, Test A, confirm thatthe ANM is able to estimate a deadzone nonlinearity.Next, the process in Figure 6.42 was simulated using the deadzone and the signal in Figure 6.47(a)was used as the reference signal. The resulting input and output signals, u(t) and y(t) appear inFigure 6.48 and were used by the ANM to estimate the table values of the static input output rela-tion. The results, Test B in Table 6.4, show that the ANM was able to estimate the deadzone. Exceptfor the one value at zero, the ANM performed very well.Finally, to see if the previous results could be improved upon, as in previous simulations, theANM estimation was turned off for ten samples after a reference change. The signals that were usedInput 2 1.5 1 0.5 0 -0.5 -1 -1.5 -2Output 1.5 1 0.5 0 0 0 -0.5 -1 -1.5Test A 1.48 0.99 0.49 0.01 -0.00 -0.01 -0.47 -1.02, -1.41Test B 1.48 0.98 0.52 0.02 0.35 0.04 -0.52 -0.94 -1.49Test C 1.48 0.99 0.49 0.02 0.00 0.03 -0.50 -0.97 -1.48Table 6.4 ANM results for estimation of deadzone65Chapter 6: Simulation StudyFigure 6.48: The (a) process input and (b) process output from closed loop operation with Dahlin controllerthen had most of the dynamics removed. The result of the estimation, Test C in Table 6.4, show animprovement in the table values estimated.These two simulations show the ANM ability to model typical nonlinearities. The improvedresults obtained by turning estimation off after reference changes suggests that this be done in allsituations. There are several different nonlinearities which could appear in an industrial process.One is hysteresis, which will require some modification to the ANM so that this nonlinearity can bemodelled. It will require two separate gain curves, one which is built as the input increases overtime and another which is built as the input decreases over time. This will require that the directionof the input be tracked, however these small changes will enable hysteresis to be modelled.6.5 Effect of Process Changes on Minimum Achievable Output VarianceThe estimated minimum achievable output variance will change as process properties change.The goal of the following simulations is twofold. First, they will show that the operation of theminimum achievable output variance estimator is consistent with the aspects behind its calculation.Secondly, the results will be used to show how changes in estimated minimum achievable outputvariance can be related to changes in the process. This change in estimation can be used to determineprocess changes, as well as evaluate the performance of the present controller66Chapter 6: Simulation StudyFigure 6.49: Original process and associated propertiesThe process used for these simulations appears in Figure 6.49, which has a first-order plant witha time delay of four. The disturbance filter is first-order, ARMA(1,1) and the noise sequence hasa N(0,0.36) distribution. The process uses a Dahlin controller that compensates for the delay andmoves the closed loop pole of the process to 0.3. The minimum achievable output variance and theprocess output variance for this process are listed in Figure 6.49, along with the closed loop noisetransfer function of the process, s(q_, . This transfer function is estimated so that the minimumachievable output variance can be determined. For this process, the transfer function has a second-order denominator, while the numerator polynomial is more complicated, having order five. Lastly,the F(q-t ) polynomial whose coefficients must be found to determine the minimum achievable outputvariance, is given.The results of one simulation appear in Figure 6.50. The autocorrelation function of the estimatednoise sequence is white because correlation at all lags is near zero. The estimated noise has avariance of 0.3598, which is nearly identical to the process noise variance. The FMVRE delayestimator successfully estimated the correct delay, while the Markov parameters, that are determinedfrom estimated Laguerre gains, are very close to the process values. Based on this, the estimatedminimum achievable output variance is very close to the calculated value for this process. This singlesimulation result is confirmed by repeating it five hundred times using a different N(0,0.36) noisesequence. The histogram in Figure 6.51 shows that the mean minimum achievable output variance67Estimated Noise Autocorrelation FunctionCO006 •8 •CC= 0.3598^•4 • •2 • •0^2^4^6^8Lags(a)F(cf1 ) Coefficients10 12 14^160.90.80.70.60.50.40.30.20.1Number of Filters: 12Time Scale: .3F(cf 1 ) = 1 + 0.2914 (1 1 + 0.1699 (4-2+ 0.1234 q30 2 3^4Lags(c)5 6 7Chapter 6: Simulation StudyTime Delay Estimation4.543.5321.50.5oo 100 200 300 400 500^600Time(b)700 800 900 1000Figure 6.50: The (a) estimated noise autocorrelation function, (b) the FMVRE delayestimation, (c) the estimated F(q-1 ) coefficients and (d) 6-n2iy and ay2 for the original processOriginal Process HistogramFigure 6.51: Distribution of the results of five hundred simulations for the original processC.20161841.61210.8Number of Filters: 12Time scale^: 0.3Points per test^: 1000Number of tests : 5000.85 0.9^0.95 yom2 v 1.05^1.1^1.15^1.282inv 68Time Delay Percentage: 1002 v = 0.4034eil,= 0.4131Chapter 6: Simulation Studyis just slightly above the calculated value, while the individual tests are closely distributed aroundthe calculated value.The first property changed was the delay, and it was reduced from four to two. This reductionof delay changes the F(g-1 ) polynomial which will have two, rather than four coefficients:F(q-1 ) = 1 + 0.27q-1and then the minimum achievable output variance becomes:T 2- Inv = 0.36 (1 + 0.272 ) = 0.3863The transfer function from the noise to the output becomes:y(t)^R(g-1 )^1 - 0.7q-1 - 0.58q-2 0.28q-3e(t)^S(q-1 ) -^1 - 0.97q-1 0.20q-2(6.110)(6.111)(6.112)Figure 6.52: The (a) estimated noise autocorrelation function, (b) the FMVRE delayestimation, (c) the estimated F(q -1 ) coefficients and (d) (5- ;„ and a3 for a change in delay69Chapter 6: Simulation Studywhich changes the process output variance to:o-2 = 0.5277^ (6.113)The results of a single simulation reflects these changes, and Figure 6.52(b) shows the FMVRE delayestimator correctly estimated the change and thus there was only two coefficients in F(q -1 ). Theestimation of the noise sequence variance remained closed to 0.36. Five hundred simulations wererun to examine the consistency of the estimation. The mean estimated minimum achievable outputvariance was found to be 0.3852 and its distribution was similar to Figure 6.51, which suggestsFigure 6.52 is a typical result. This simulation shows that once the FMVRE estimates a change inthe process delay, there should be a change in the output variance and the minimum achievable outputvariance becomes smaller because F(q -1 ) has fewer coefficients.Returning the delay to four, the next change was made to the plant transfer function: B (q1 )^0.25q-4 (6.114)A1(q-1 ) — 1 — 0.75q-1which results in a new closed loop noise transfer function of:R(q-1 )^1 — 1.45q-1 0.65q-2 — 0.09q-3 — 0.70q-4 0.80q-5 — 0.21q-6 (6.115)S(q-1 ) — 1 — 1.72q-1 0.93q-2 — 0.15q-3 — 0.17q-4 0.28q-5 — 0.11q-6This change makes the process output variance:o-2 = 0.5591 (6.116)The minimum achievable output variance is unchanged, however, as the F(q -1 ) polynomial remainsthe same. The results of a simulation, Figure 6.53 confirms this, as the estimated minimum achievableoutput variance remains near the calculated value of 0.4034, despite the change in the output variance.Figure 6.53(a&c) shows that the estimate of the F(q-1 ) coefficients and noise sequence remainsunchanged for the first process, while the histogram in Figure 6.54 confirms that the estimationremains virtually unchanged by the plant variation. This simulation shows that variation in only theoutput variance suggests that the plant transfer function has changed.70Changed Plant Histogram20 Number of Filters: 12IS Time scale^: 0.316 Points per test^: 1000Number of tests : 500Time Delay Percentage: 1002om v= 0.4034a`mv =0.41451412sr641.21.0.95 62 / 2^1.05^1.1^1.15^1 2v atn V00.8^0.85 0.9Chapter 6: Simulation StudyFigure 6.53: The (a) estimated noise autocorrelation function, (b) the FMVRE delayestimation, (c) the estimated F(cf 1 ) coefficients and (d) ?IL, and av2 for a change in plantFigure 6.54: Distribution of the results of five hundred simulations for the change in the process plant71Chapter 6: Simulation StudyOnce the plant was reverted back to its original value, the disturbance filter was changed:C(q-1)^1 — 0.2q-1A2(q-1 ) — 1 — 0.67q-1which changes the F(q-I ) polynomial coefficients to:F(q-1 ) = 1 + 0.47q-1 + 0.32q-2 + 0.27q-3and then the minimum achievable output variance becomes:= 0.36(1 + 0.472 + 0.322 + 0.272 ) = 0.49152Once again, the closed loop noise transfer function changes:y(t)^R(q-1 )^1 — 0.5q-1 + 0.06q-2 — 0.70q-4 + 0.14q-5e(t)^S(q 1 )^1 — 0.97q-1 + 0.20q-2and then the process output variance becomes:CrILIV(6.117)(6.118)(6.119)(6.120)= 0.7434 (6.121)From Figure 6.55, the corresponding changes in the F(q -I ) coefficients, the minimum achievableoutput variance and output variance are evident, while the estimated variance of the noise sequenceremains the same. Figure 6.56 confirms that the mean value of the minimum achievable outputvariance moved towards 0.4915, the new level. Thus a change in the disturbance filter becomesnoticeable, as both o and ol,„ change.To see the effect the noise sequence had on the estimation of the minimum achievable outputvariance, it was changed from N(0,0.36) to N(0,1.21) in the original process. This change in thenoise variance results in new process output variance and minimum achievable output variance:2i„, = 1.356(6.122)cry = 2.055One simulation was performed and it showed that the output variance and minimum achievableoutput variance level did change to the new values. This was a direct result of the change in noisevariance, as the F(q -I ) coefficients and delay estimate remained unchanged. Five hundred simulations72201618Number of Filters: 12Time scale^: 0.3Points per test^: 1000Number of tests : 500Time Delay Percentage: 100v = 0.5175= 0.491514F6F400.8^0.85^0.9^0.95^1^1.05^1.1^1.15^1.2°Th V 6111vVn1Chapter 6: Simulation StudyFigure 6.55: The (a) estimated noise autocorrelation function, (b) the FMVRE delay estimation, (c) theestimated F(q-1 ) coefficients and (d) and o for a change in disturbance filterChanged Disturbance HistogramFigure 6.56: Distribution of the results of five hundred simulations for a change in disturbance filter73Chapter 6: Simulation Studywere performed and they confirmed this finding. As a result, if both o and o ^change while F(q-1 )remains the same, then the noise variance must have changed.The final part of the process to be changed is the controller. During the derivation of the methodto find the minimum achievable output variance, it was determined that the minimum level could notbe reduced any further by a feedback controller. The controller was changed to:N(q-1 )^0.7 — 0.5040q-1D(q-1 ) — 0.28 — 0.084q-1 — 0.196q-4^(6.123)which is a Dahlin controller which moves the process pole to 0.3 for the plant:B (q-1 )^0.28q-4Ai(q-1 ) — 1 — 0.72q-1^(6.124)This change in controller causes the closed loop noise transfer function to become:y(t)^R(q-1 )^1 — 0.7q-1 0.12q-2 — 0.70q-4 0.28q-5 (6.125)e(t)^S(q-1 ) — 1 — 0.97q-1 0.20q-2 0.12q-4 — 0.12q-5and this changes the process output variance to:CI 2 = 0.6744 (6.126)but it does not change the minimum achievable output variance, because the delay, noise and F(q -1 )coefficients remain unchanged. The results of one simulation, shown in Figure 6.57 proves that thecontroller only modifies the output variance, if compared to Figure 6.50. It is promising to note thatthe histogram in Figure 6.58 has a mean estimated minimum achievable output variance level closeto the value found in Figure 6.51, while the shape of the histogram remains basically unchanged.Various properties of the process have been changed, in order to study the effect on estimation ofthe minimum achievable output variance. This knowledge can be used to determine what has changedin the process when either o and 011V change. For example, a change in delay will be estimated bythe FMVRE, while a change in noise variance will be estimated by the RELS. Other changes in theprocess output variance become evident in a more subtle manner. If only the process output variancehas changed, this indicates that the plant characteristics have been modified. Meanwhile, if just thevalue of the F(q-I ) coefficients change, then the disturbance filter has changed values. Using thesedirect indicators of change to signal operational changes in the process, it is possible to focus attentionon a specific aspect of the process. Finally, it should be remembered that changes to the controllershould only change the output variance and not affect the minimum achievable output variance.74Chapter 6: Simulation StudyFigure 6.57: The (a) estimated noise autocorrelation function, (b) the FMVRE delay estimation,(c) the estimated F(cil ) coefficients and (d) erm2 , and o-y2 for a change in controllerFigure 6.58: Distribution of the results of five hundred simulations for a change in the controller75Chapter 6: Simulation Study6.6 Choice of Laguerre Time Scale and Number of FiltersThere are two parameters of a Laguerre network that must be chosen before it is used to modela process. They are the number of filters used in the network and the time scale of these filters. Thefollowing simulations investigate the effect these parameters have on the estimation of the minimumachievable output variance. A brief discussion on what values for these parameters are usually made,is followed by the simulations. The results are discussed and following this suggestions for generalchoices are made.As mentioned, the set of Laguerre functions can completely represent the impulse responseof a linear process only if an infinite number of filters are used. As it is impossible to use aninfinite number, the number of filters used in the Laguerre network must be chosen so that a goodapproximation of the process is obtained. It has been found that four to eight Laguerre filters usedin the network provides a reasonable approximation [15][25] of a sampled-data plant. This workshowed that good approximation also depends on the choice of the Laguerre filter time scale. It isknown that the optimum choice for the time scale will depend on the process to be modelled. If aprocess exhibits a heavily damped behavior, which suggests that it has real poles, the best choice forthe Laguerre filter time scale is a value which corresponds to the dominant process pole. However,if the Laguerre network is modelling a process which has complex poles, evident as the process willexhibit a lightly damped behavior, a different approach to choosing the time scale is required. Fushowed in [15] that the Laguerre network approximation of a sampled-data system improved as thetime scale is moved toward zero.However, the previous conclusions were drawn from situations where the input to be used forestimation was known. In this application, the transfer function to be modelled is driven by anunknown white noise sequence, that also must be estimated along with the unknown parameters. Theprocess determines the closed loop noise transfer function, which is given by:R(q--i)y(t) — S(q-1)^e(t)^ (6.127)where the numerator polynomial is often unlike any polynomial usually found in a sampled-dataplant. As a result of these differences, it has been found that the best choice for the number of filters76Chapter 6: Simulation Studyand time scale for the Laguerre network, when used with RELS estimation, depends on properties ofS(q-1 ) and R(q-1 ). The choice of filter number and time scale depends on whether the closed loopnoise transfer function has real or complex closed loop poles, that is the roots of S(q-1 ), or if thenumerator polynomial, R(q-1 ), is strictly positive real (SPR). Combinations of these situations willbe simulated to see if some general choices for the time scale and filter number can be made toobtain good estimation.6.6.1 Real Closed Loop Noise Transfer Function Denominator RootsThe first two simulations examine the case where the closed loop poles are real and the R(q -1 )polynomial is strictly positive real. These examples are the least complicated and most likely theleast encountered situations, however, they are important as they show the behavior of estimationin a uncomplicated situation.The first process examined appears in Figure 6.59 and it uses a gain as the controller, so thatR(q-1 )=C(q-1 ). The properties of the process are also shown in Figure 6.59, indicating the transferfunction to be identified, which has real roots, and the calculated minimum achievable output variancevalue which must be estimated. Because the C(q-I ) polynomial is first-order, it is strictly positive real.The results of one simulation, Figure 6.60, shows that the autocorrelation function of the estimatednoise sequence is white and has a variance close to the expected value, while the delay and F(q 1 )coefficients were correctly estimated. As a result, the estimated minimum achievable output variancee(t)N(0,0.49)y(t)zeroreferenceae = 0.4902 = 0.6435= 0.6125y(t)^R(q-1 )^1 —0.2q-1e(t)^S(q -1 )^1 — 0.7q -1 + 0.1225q-2F(q-1 ) = 1 +0.5q-1Figure 6.59: Gain controlled process77Chapter 6: Simulation StudyFigure 6.60: The (a) estimated noise autocorrelation function, (b) the FMVRE delay estimation,(c) the estimated F(cf 1 ) coefficients and (d) ?f,,2,,,, and o-y2 for process with real polesFigure 6.61: Distribution of the results of five hundred simulations for process with real poles78Chapter 6: Simulation StudyEffect of Time Scale and Filter Number on Estimation0.620.6u 0.58i0.560.540.520^0.1^0.2^0.3^0.4^0.5^0.6^0.7^0.8^0.9^1Time ScaleFigure 6.62: Effect of changing the time scale and filter number on erl, for process with real polesis very close to the calculated value, Figure 6.60(d). This single test was confirmed by Figure 6.61,which shows that after five hundred simulations using different N(0,0.49) noise sequences, the meanestimated minimum achievable output variance was just under the calculated value.Now in order to obtain an idea of how the number of filters and time scale used in the Laguerrenetwork affects the estimated minimum achievable output variance, the test shown in Figure 6.61was performed with different time scales and filter numbers. The mean value obtained after each testis then an indicator of the performance of that choice of time scale and filter number. The results arein the plot of Figure 6.62, which presents variance versus the time scale for a different number offilters, listed next to the corresponding curve. The plot also shows the calculated minimum achievableoutput variance along with minus five percent of this value. Lines are used only to join results whichshare the same number of filters used. When only two filters are used, the estimated minimumachievable output variance already overfits the calculated value and this continues as the number offilters is increased. Except for when two filters are used, the choice of time scale does not have mucheffect on the minimum achievable output variance estimated. These results suggest that the Laguerrenetwork is not behaving as expected, perhaps because it is being used with RELS. Before furtherconsideration of the results, another process will be examined to confirm these findings.79e(t)y(t)zeroreferenceChapter 6: Simulation Studycr = 1.6902 = 2.1091.2arn, = 1.7576y(t) R(q-I ) _ 1- q- I - 0.19q-2 + 0.22q-3 - 0.03q-4e(t) - s(q- I) -^1- 1.2q-I + 0.35q-2F(q-I ) = 1 +0.2q-1Figure 6.63: Dahlin controlled process with real polesFigure 6.64: The (a) estimated noise autocorrelation function, (b) the FMVRE delay estimation, (c) theestimated F(q-I ) coefficients and (d) O-;2„v and o for Dahlin controlled process with real poles8081.751.71.651.6Chapter 6: Simulation StudyThe next process examined appears in Figure 6.63 and the major difference is that this processis controlled by a Dahlin controller which makes the closed loop noise transfer function morecomplicated. Despite this, the numerator polynomial is still strictly positive real, so the resultsof this simulation can be compared to the previous results. The results of one simulation appears inFigure 6.64, which uses twelve filters and a time scale of 0.3 and the minimum achievable outputvariance is quite closely estimated. Examination of Figure 6.65 shows a similar trend as found inFigure 6.62, in that the minimum achievable output variance was overfitted for almost all choicesof filter number. It is believed that due to the increased complication of the transfer function to beidentified, the number of filters required before the overfitting begins is larger.The results of single simulations are examined in order to determine the reason why the minimumachievable output variance is overfitted. The mean value of 1.69 found from Figure 6.65 can beexplained by closely examining Figure 6.64(a&c). While the estimated noise autocorrelation functionis white, the variance is less than expected, also, the F(q -1) coefficient is larger than expected. Thetwo effects combine and result in &Iv being less than the calculated value. Similar results werefound from the gain controlled process. It appears that the Laguerre network is including the noiseEffect of Time Scale and Filter Number on Estimation0^0.1^0.2^0.3^0.4^0.5^0.6^0.7^0.8^0.9Time ScaleFigure 6.65: Effect of changing the time scale and filternumber on &Iv for Dahlin controlled process with real poles81Chapter 6: Simulation Studysequence into its estimation of the closed loop noise transfer function, resulting in an incorrect noisevariance and F(q-1 ) coefficients. However, close estimation is still possible for both cases for a widerange of time scales and filter number.6.6.2 Complex Closed Loop Noise Transfer Function Denominator RootsIn the processes considered next, the polynomial R(q -1 ) has complex poles. The time scale of theLaguerre network used is real and it has been found that the number of filters used by the networkmust be increased [15] in order to arrive at satisfactory approximation of the process. By simulatinglightly damped processes, it can be determined if the number of filters need be increased to obtaingood estimation of the minimum achievable output variance.Figure 6.66: The (a) estimated noise autocorrelation function, (b) the FMVRE delay estimation, (c)the estimated F(q') coefficients and (d) b -,2„, and o for process with complex poles82Chapter 6: Simulation StudyThe first process used, Figure 6.59, controlled by a pure gain, is used once again, but this timethe gain is two so that the roots of R(q -1 ) are complex. The effects of this change make:^(q 1 )^1 - 0.2q-1Cry2 = 0.819^ (6.128)S(q - 1)^1 - 0.7q-1 0.6q-2so the controller does not change the minimum achievable output variance. The results of onesimulation, Figure 6.66, illustrates how good the estimation performs for four filters and a time scaleof 0.4. The autocorrelation of the estimated noise sequence for this simulation, Figure 6.66(a), isnot exactly an impulse, as was found in the previous simulations, and its variance is larger than theexpected value of 0.49. The oscillatory nature of the autocorrelation function suggests that it hasbeen influenced by the complex roots in S(q 1 ). However, the estimated minimum achievable outputvariance comes very close to the calculated value of 0.6125 in this simulation. From Figure 6.67, itis apparent that the mean value of er 2 v is just slightly larger than the calculated value.The plot of the mean estimated minimum achievable output variance values, Figure 6.67, forvarious filter numbers and time scales suggests that more filters are required to obtain results similarto Figure 6.62. The effect of the complex poles becomes evident when results using two filters isexamined. The noise sequence and F(q -1 ) coefficients appear in Figure 6.68, and while the F(q -1 )Effect of Time Scale and Filter Number on Estimation0.80.750.70.650.60.550^0.1^0.2^0.3^0.4^0.5^0.6^0.7^0.8^0.9^1Time ScaleFigure 6.67: Effect of changing the time scale and filter number on er ilv for process with complex poles83B(q 1 )^0.2q- 2A(g-1 ) ^1 — 0.8q- 1e(t)y(t)zeroreferenceChapter 6: Simulation StudyFigure 6.68: Effect of reducing filter number on (a) noiseautocorrelation function and (b) the estimated F(q -1 ) coefficientscoefficient is correctly estimated, Figure 6.68(a), shows that the problem lies with the estimation ofthe noise sequence. The autocorrelation function has a more oscillatory shape than previously foundand the variance is much higher than expected. The cause of the error in the estimation resultsfrom the error in the noise variance, which is caused by the complex pole pair.The next process considered, Figure 6.69, is controlled by a Dahlin controller. This controller ispoorly designed however, using a plant transfer function slightly different from the true one:(6.129)•zs = L69452r = 2.45132Gm, = 1.7576y(t) R(q-1 ) 1 - q-1 - 0.19q-2 + 0.22q-3 - 0.03q-4e(t)^S(g -1 ) =^1 - 1.2q-1 + 0.6q-2 - 0.25q-3F(q-1 ) = 1 + 0.2q-1Figure 6.69: Process controlled by Dahlin controller with complex poles84Estimated Noise Autocorrelation Function Time Delay Estimation100 200 300 400 500^600^700Time(b)1612 14102^4^6^8Lags(a)F(q-1 ) Coefficients900800 10002.58> ,1000100^20()^300^400^500^600^700^8(0^900Time(d)5 63^4Lags(c)54.543.53e 2.521.50.5Process Variances00^Figure 6.70: The (a) estimated noise autocon-elation function, (b) the FMVRE delay estimation, (c) theestimated F(q-t ) coefficients and (d)^and ay2 for Dahlin controlled process with complex polesNumber of Filters: 8Time Scale: .3F((f l ) = 1 + 0.2214 cf l10.90.80.70.60.50.40.30.20.100^20.80.60g 0.4O0.20- — - — - — - — -- — - — Calculated Output Variance^ Output VarianceCalculated MinimumAchievable Output Variance .Minimum AchievableOutput Variance0.50^(4= 1.621Chapter 6: Simulation StudyThis mismatch causes the introduction of a pair of complex poles at 0.18 ± 0.515j. Processes of thisform, with a complex pole pair resulting from a mismatch between the true plant and the model usedfor control, are common and therefore important to study. The results of one simulation appear inFigure 6.70 and show that the autocorrelation function is once again oscillatory, though in this caseii-„2 ,v is quite close to the calculated value. The result of five hundred simulations, 1.741, found fromFigure 6.71, confirms that the mean estimated (5-1211V for this choice of time scale and filter number isnearly perfect. Figure 6.71 plots the remaining mean estimated minimum achievable output variancevalues and comparing this plot to Figure 6.65, the increase in the number of filters required toobtain a similar result is evident. When four filters are used rather than eight, Figure 6.72, &I v is852.121.91.81.7Chapter 6: Simulation StudyEffect of Time Scale and Filter Number on Estimation0^0.1^0.2^0.3^0.4^03^0.6^0.7^0.8^0.9^1Time ScaleFigure 6.71: Effect of changing the time scale and filter numberon iim2 , for Dahlin controlled process with complex polesoverestimated because the autocorrelation function of the estimated noise sequence is more oscillatorythan in the previous case, resulting in an inflated noise variance.These process simulations suggest that the number of filters must be increased in order to getthe same results as for a similar process which has real poles. However, the plots of Figure 6.71and Figure 6.65 suggests that it is possible to get estimation within plus or minus five percentFigure 6.72: Effect of reducing filter number on (a) noiseautocorrelation function and (b) the estimated F(q -1) coefficients86Chapter 6: Simulation Studyof the calculated minimum achievable output variance, with the same choice of Laguerre networkparameters, for both situations.6.6.3 Closed Loop Noise Transfer Function Numerator Not SPRThe final situation, when the numerator polynomial of the closed loop noise transfer function isnot strictly positive real, is considered. The numerator polynomial is given by:R(q-1 ) = C(q 1 )Ai(q 1 )H(q-1 )^ (6.130)so it is made up of the process disturbance numerator, the plant denominator and the denominator ofthe controller. Combining these polynomials can easily result in R(q -1 ) having a high order, whichwill then most likely mean it is not SPR. The SPR property can be explained in the following manner.A transfer function will have a SPR numerator if a white noise input results in an output which hasan autocorrelation function close to an impulse. The higher the numerator order, the more noisecoloring that occurs, and so the autocorrelation function will no longer be an impulse. Then therewill be a large chance that the numerator is not SPR.It is important that the R(q -1 ) polynomial is positive real for the estimation of the minimumachievable output variance. According to a condition of RELS estimation, it is necessary for R(q -1 )to be SPR for convergence of the estimates to be assured. However, the nature of the conditionmeans that the estimates may converge even if R(q -1 ) is not SPR. The first example illustrates thise(t)N(0, 0.36)zeroreference2a = 0.36cr2Y =0 92362amv = 0.6256y(t) = R(q-1 )^1 - 0.2q-1 -q-4 +0.2q-5e(t)^s(q-1 ) 1 -0.8q-1F(q-1 ) = 1 +0.6q-1 +0.48q-2 +0.384q-3Figure 6.73: Process controlled by Dahlin controller without SPR R(q -1 )87Chapter 6: Simulation StudyFigure 6.74: The (a) estimated noise autocorrelation function, (b) the FMVRE delay estimation, (c) theestimated F(q-1 ) coefficients and (d) 6-,2„, and cr for Dahlin controlled process without SPR R(q-1 )fact and shows how the choice of filter number and time scale are affected by the realness condition.The second process illustrates limitations of the estimator caused by this condition.The process in Figure 6.73, similar to those used previously, but it uses a Deadbeat controllerwhich results in a R(q-I) polynomial that is not SPR. Results of a single simulation, Figure 6.74, showthat the error between 6-1, and o appears to be caused by an inability to obtain a correct estimateof the noise sequence and incorrect F(q-1 ) coefficients. It has been observed for various processeswith a numerator that is not SPR, that a good value can be obtained despite poor estimation ofthe noise sequence and F(q-I ) coefficients, . The histogram constructed from the & v values obtainedfrom five hundred simulations, Figure 6.75, confirms that the mean is indeed close to the calculatedvalue. However, the distribution of the histogram differs, in that it has a larger variance, which is88Chapter 6: Simulation StudyFigure 6.75: Distribution of the results of five hundredsimulations for Dahlin controlled process without SPR R(cf1)also a common occurrence for processes with a numerator that is not SPR. The unrealness of thenumerator affects the choice of the Laguerre parameters used to obtain good estimates. A plot ofthe mean erm2 determined for different choices, indicates the number of filters required to obtaingood estimation has increased a great deal over previous situations. Figure 6.76 shows that even ifFigure 6.76: Effect of changing the time scale and filter numberon '4, for Dahlin controlled process without SPR R(e)89Chapter 6: Simulation Studytwenty four filters are used, the mean value does not drop below the calculated value. Also, thechoice of the time scale becomes more important then previously, with a value in the range [0.2,0.4]a suitable choice. It appears that if R(q-1 ) is not SPR, the estimation of the minimum achievableoutput variance is noticeably affected, but it can be compensated for by the prudent choice of a timescale and filter number.To illustrate the limitations of the RELS method, when a numerator is not SPR, the process inFigure 6.59 is used, except that the noise filter numerator has been changed. This new choice ofthe C(q-1 ) polynomial is:c(q—i)^1 + 1.5q1 0.75q-2^(6.131)Figure 6.77: The (a) estimated noise autocorrelation function, (b) the FMVRE delay estimation, (c)the estimated F(q -1 ) coefficients and (d) 6-;,2i , and ay for process without SPR R(cf1)90Chapter 6: Simulation Studywhich because the controller is a pure gain equals the R(q-1 ) polynomial. The changes to the otherprocess properties is:a2 = 6.1521 1T1? =S q-1^1+1.5q-1^q-+0.75 2 1 —0.7q-= 2.2616^F(q-1 ) = 1 + 2.2q-1This polynomial is used extensively in Pseudolinear Regression literature [21] as an example forwhich RELS estimation will not converge. Attempts to estimate the minimum achievable outputvariance, as in Figure 6.77, illustrate its failure to come close to a m2 v . Both the estimated noisesequence and F(q-1) coefficients are nowhere near their calculated values. Figure 6.78 confirms theinability of the method to converge, as the mean eir m2 „ values are not close to the calculated value.This process illustrates that estimation may not be possible for all processes, however, the problemwill be evident after the autocorrelation function of the estimated noise sequence is examined. Theautocorrelation function will have large correlation at numerous lags, as shown in Figure 6.77(a),when the nonrealness of the R(q-1) polynomial is affecting estimation. Fortunately, this problemdoes not seem to occur for the majority of cases when R(q-1 ) is not SPR.The results of these simulations suggest that some trial and error may be necessary to obtainsuitable values for the time scale and number of filters used in the Laguerre network. There is oftenEffect of Time Scale and Filter Number on Estimation54.58g 43.530^0.1^0.2^0.3^0.4^0.5^0.6^0.7^0.8^0.9Time ScaleFigure 6.78: Effect of changing the time scale and filter number on &,2,R, for process without SPR R(q-1)91Chapter 6: Simulation Studya wide range of Laguerre parameters which result in good estimation, that being plus or minus fivepercent of the calculated value. Based on the previous simulations, some general recommendationwill be made.The estimated noise sequence autocorrelation function is an important tool that should be usedwhen choosing a filter number. If the autocorrelation function has correlation at lags above zero, thisindicates the number of filters should be increased. If increasing the filter number does not improvethe autocorrelation function shape, then estimation may not converge, due to the fact the R(q -1 )polynomial is not strictly positive real. In this situation, an alternative estimator, such as MaximumLikelihood should be employed. If the noise autocorrelation function has an impulse shape, thendecreasing the number of filters used should be attempted, as long as the autocorrelation function ofthe estimated noise remains white. Using the lowest number of filters which result in an impulseautocorrelation, will provide the best estimate of cri2m, Generally, the simulations obtained goodresults for a range of Laguerre network filters, regardless of whether R(q 1 ) was SPR or S(q -1 ) hasreal or complex poles. The choice of time scale does not appear as crucial as the choice of filternumber. Except when R(q-1 ) is not SPR, the 6-12nv value found varied only slightly for different timescale values for a constant filter number. Processes where R(q-1 ) positive realness is a concern, thetime scale should be used in the range of [0.2,0.5]. Based on various simulations, this time scalerange can be used with a moderate degree of confidence, in all situations.6.7 ConclusionsThere were several different simulations presented in this chapter, all attempting to examinerelevant details encountered when using the estimation methods presented for monitoring. The resultsof each section will now be summarized in order that the important details can be emphasized.The first simulations performed were concerned with determining the possible causes of biaseddelay estimation and suggesting solutions able to remove the bias. If the autocorrelation functionof the disturbance output is not zero for lags which correspond to the possible interval of delay, abiased delay estimate may occur. The bias may occur at one sample, but this also corresponds tothe delay caused by the sampling procedure and does not need to be included in the possible range92Chapter 6: Simulation Studyof delay, if the process is known to have delay. It is also possible for the bias to occur at a valueother than one. This bias can be removed by filtering the signals as explained previously. The othersource of bias arises when the input signal does not contain sufficient excitation to estimate the delay.Unfortunately, it is difficult to justify removing this bias based solely on the fact that the delay isrequired only for monitoring. An equation was determined that only requires the sampling intervaland the lower frequency required to estimate the delay can be determined. Thus by examining thepower spectral density function of the input it is possible to determine if the excitation requirementsare met. In general, the safest course to follow when estimating the delay is to examine the Elfunction to check if the maximum value definitely outweighs the other values. If this occurs, thenthe delay can be used with confidence.The next study was concerned with determining the performance of the ANM on closed loopprocesses with noise. The first simulation examined a linear stochastic process and examined the useof low pass filters and the benefit of turning off the ANM for a short period after a reference change.The results showed the benefits of these techniques and it is suggested that they are implementedif it is not too difficult. Any nonlinearity in the process will degrade its performance, so it isimportant that they can be estimated by the ANM. Two typical nonlinearities were used, a deadzoneand saturation, and both were correctly estimated by the ANM. By including filtering to reduce thenoise and turning off the ANM after reference changes, the ANM is capable of estimating the staticinput-output relation.Focus then shifted to the estimation of the minimum achievable output variance and interest wasin determining a suitable choice for the Laguerre network parameters and determining the feasibilityof using the estimated properties, required to calculate 61,, as indicators of process change. TheFMVRE will estimate changes in delay and RELS will estimate any changes to the noise variance.Only o will change if the process plant changes, while if only the Markov coefficients change,modifying "6-2„,, then the process disturbance filter is different. Using this knowledge it is possible todetermine changes in the process. It was found that the choice of Laguerre parameters depends onthe properties of the closed loop noise transfer function. The best way to determine if the Laguerrenetwork is suitably modelling the closed loop noise transfer function is by examining the estimated93Chapter 6: Simulation Studynoise autocorrelation function. As few filters as possible should be used so that it has a whitenoise autocorrelation function. In general more filters will be required when the denominator of theclosed loop noise transfer function has complex roots or if the numerator in not SPR. There is alsothe possibility that the estimates will not converge to their true values when the numerator in notSPR. This is a limitation of RELS estimation, which can be overcome by using a more complicatedestimation algorithm. Fortunately, this problem is evident as the noise autocorrelation function willnot be near an impulse and the problem itself does not seem to occur that frequently. While sometrial and error may required to determine suitable Laguerre values, a good starting point is aroundten filters with a time scale in the range of [0.2,0.5].While investigating the method used to determine the minimum achievable output variance, theadvantages of using a Laguerre network over a time-series model were displayed. Despite changingevery part of a process, which resulted in a different closed loop noise transfer function, one choiceof Laguerre parameters was able to provide a good estimate of al2114• Meanwhile, an ARMA modelwould have to be chosen with some care if a time-series representation was used and it would likelyhave to be modified as the closed loop noise transfer function changed. There also is comfort in theknowledge that increasing the size of the Laguerre network will result in improved estimation becausethe computational load will increase. With an ARMA model there is no guarantee that the estimationwill improve though the computational load also increases. While investigating a suitable choiceof Laguerre parameters, it was found that a wide range will result in good estimation, a robustnessnot possessed by an ARMA model. While the Laguerre gains do not offer the same insight into aprocess as do the roots of the time-series polynomials, in this application it is not important as itis the Markov coefficients that are of interest and must be calculated regardless of which model isused. In conclusion, the use of a Laguerre network provides good estimation of a changing closedloop noise transfer function while avoiding the problem of choosing the ARMA model order.94Chapter 7: Industrial DataChapter 7Industrial DataThe most crucial test faced by the monitoring tools is now presented, the application to industrialdata. The data were gathered from three different processes, a pressure screen, a reject refiner and aKamyr digester. Monitoring was performed in an effort to test the operation of the tools on industrialprocesses as well as demonstrate the usefulness of monitoring. For each process presented, a briefoverview of its operation is followed by the goal of monitoring in this specific case. The outputand input signals are presented and then analyzed. The monitoring is performed, outlining what wasrequired in order that the presented results could be achieved. Finally, attempts are made to verifythe results obtained and once they are considered reasonable, the results are interpreted to evaluatethe process performance.7.1 Pressure Screen Reject FlowA pressure screen is an integral component of the pulping process. They are often used afterwood chip refining to remove unwanted components from the pulp, such as oversized fibres, dirt andbark. A screen will typically have two outlets, one for the accepted pulp and the other for rejectedpulp. The reject flow is usually held constant, removing a sufficient percentage of the unwantedcomponents in the pulp. The pressure screen reject flow loop is a basic flow loop, where the outputis the flow of pulp and the valve position is the input. In this application, the flow is controlled bya PI controller and the process is sampled every two seconds.This flow loop is typical of the type of processes which can be monitored using the methodspresented. The data will be examined to determine the performance of the PI controller. Oncethe minimum achievable output variance is determined, it can be compared to the process outputvariance to determine any possible improvement.The process data appears in Figure 7.79, and shows approximately 1200 samples or about fortyminutes of the screen operation. The flow reference was six hundred litres per minute, and theoutput stays about this value, and has a mean of six hundred litres per minute. The input is also95ry=6001000^1200200^400^600^800Time (2 seconds)Output800 ^700600 5001400c4) 300 -200 ^05048464.)^4442(a)InputChapter 7: Industrial Data40 0^200^400^600^800^1000^1200Time (2 seconds)(b)Figure 7.79: The pressure screen reject flow input and outputshown, the position of the valve, and it is given as a percentage of the fully open position. Beforeapplying the monitoring tools, the input and output signals are examined. Figure 7.80(a) shows theautocorrelation function of the process output. The process variance after 1200 samples was 5128 andthere is significant correlation at most lags, even up to fifteen. This indicates that unless the delay isfifteen samples or thirty seconds, which is unlikely for a flow loop, the process is definitely not underminimum variance control. There is an indication at this point then, that the process variance is not atthe minimum level. The variation of the output signal seems high as a standard deviation of seventyis quite large when compared to the reference. Also shown on Figure 7.80(b) is the power spectraldensity of the input signal, which will be examined to determine if it will be possible to estimate thedelay with the FMVRE. Using equation 6.109, the lower limit on the frequency is determined:fL=2 (71.V8r)-1 = (nV8(2 2 ))-1 = 0.056 Hz^(7.133)Comparing this number with the power spectral density, it can be seen that the power in the signaldies out around 0.05 Hz which indicates that there may be difficulty in estimating the delay.9630II^ II I^I50 •I^I II^II I^II I III^ I1. ^I. I^I I0.05^0.1^0.15^0.2^0.25Frequency (Hz)(b)122^4^6^8Lags(a) 10 14 16200 1000400 600Time(a)800 1200E, Function2 6Delay(d)10^124^6Delay(c)I0 12Chapter 7: Industrial DataFigure 7.80: The output autoconelation function and input power spectral density function25All Samples2015102 4^6Delay(b) 10^12FMVRE Delay Estimation01098E, FunctionE1 FunctionSamples 800-1200^•••. .Figure 7.81: The FMVRE (a) estimation of the delay and E 1 functionfor (b) all samples, (c) samples 300 to 700 and (d) samples 800 to 12009700.9. 00 0C.)-o0.5-0.5-1.5-2-2.5-3-3.51o—Y =5128•44.^- •1020Output Autocorrelation^ Input Power Spectral Density Function70 ^Samples 300-700Chapter 7: Industrial DataThe estimation of cl-11, begins by first determining the time delay of the process, using theFMVRE. The results of the estimation appear in Figure 7.81 and Figure 7.81(a) indicates that thedelay of the process is two samples. Examining the El function, Figure 7.81(b), the peak is at two,though the values at one and three are not that far off. It was shown that the power in the inputsignal dies out near the lower frequency limit required by the FMVRE to estimate the delay andthis accounts for the shape of the El function. The consistency of the estimation is checked inFigure 7.81(c) and Figure 7.81(d), by using two subsets of the data. They also suggest a delay oftwo. A delay of two shall be used throughout the remaining section for determining the minimumachievable output variance.Next, in order to find the closed loop noise transfer function and noise sequence, the numberof filters and their time scale must be decided. The filter number was chosen so that the noiseautocorrelation function had an impulse shape. The results from two different choices are shown inFigure 7.82. It can be seen that the autocorrelation function has some oscillation, but is closer toan impulse when ten filters are used. It was found that increasing above ten filters did not improvethe shape of the autocorrelation function.The minimum achievable output variance estimation is shown in Figure 7.83, for the Laguerrenetwork parameters discussed above. The results are quite similar for both choices of filter numberas the variance ratios are within five percent. It was found that despite using a wide range of filternumber (4 to 20) and time scale (0.1 to 0.7), the majority of the variance ratio lay between 2.1 and2.6. Once again the Laguerre network shows its ability to provide good estimation without too muchconcern with the parameters chosen.It is time to step back a moment and consider the results obtained, to see if they are reasonablefor this process. Some experience has been gained by controlling this process, which indicatesthat the delay is around four seconds. Also, the process is known to be quite oscillatory, which isconfirmed by the estimated noise autocorrelation function. Based on the a-priori knowledge available,the results obtained seem reasonable. As a final check of the estimates, the closed loop noise transferfunction is modelled using an ARMA(5,5) model whose parameters are estimated using MaximumLikelihood Estimation. This is the method that was used by Harris [17] and as the Maximum98Chapter 7: Industrial DataFigure 7.82: The estimated noise sequence and F(q -1 ) coefficient for (a&b) ten filters and (c&d) six filtersFigure 7.83: The estimated minimum achievable output variance using (a) ten filters and (b) six filters99Chapter 7: Industrial DataLikelihood estimator used is an off-line batch estimator, rather than an on-line recursive estimator,the estimation will not be biased and can be compared to the values estimated using RELS. Theresult of this estimation found:f{(q-1 )^1 — 1.01q-1 0.55q-3 — 0.28q-4 — 0.23q-5§ (4- 1 )^1 - 2.33q-1 1.9q-- 2 — 0.94q-4 0.5q-5t(q-1 ) = 1 + 1.316q-113-2 = 738cr2A-2Y = 2.5-The variance ratio, estimated noise variance and estimated F(q -1 ) coefficient all compare very closelywith the values obtained previously.Finally, what do the results mean? They suggest that the output variance could be reduced byover two times if a minimum variance controller was used. This is a fair amount, but if the PIcontroller is desired rather than a minimum variance controller, it should be retuned to attempt tolower the ratio. If this is not possible, then perhaps a different control algorithm, such as a Dahlincontroller should be employed. Also, it should be considered if the minimum standard deviation offorty-five is acceptable. If it is not, then changes to the process will be required to reduce it further,or feedforward control must be employed. No time should be wasted retuning the controller if theminimum standard deviation is not satisfactory.7.2 Reject Refiner Motor LoadThe next process examined is the motor load loop of a reject refiner. A reject refiner is a smallerversion of a main refiner, one which appears in Figure 7.84. The reject refiner is typically founddownstream of a pressure screen, as it refines the rejected pulp, reducing the pulp size further so thatit is suitable to be used. The mechanical operation of a double-disc refiner, which has two plates,will now be explained. The plates are rotated and pushed together by large motors and a hydrauliccylinder [111 Fibres and dilution water are fed into the center of the plate gap and a fibre padforms between the plates. The pulp works its way out the periphery of the plates, in smaller piecesthan when it entered. The energy per unit mass of wood fibres is a major factor in the resultingpulp quality and depends on the chip feedrate into the refiner and its motor load. The motor loadis changed by manipulating the gap between the plates which increases as the gap decreases. An100Chapter 7: Industrial DataFigure 7.84: A double-disc refinerinteresting property of this relation is the fact that the gain will change sign, which occurs if theplates are moved too close, as the pad of fibres breaks down, allowing the plates to move togethereasier, and ultimately they clash together. This means that there is a maximum motor load setpointwhich is obtainable. The goal of controlling the motor load is to keep it as high possible, all thewhile avoiding a pad collapse. The reject refiner motor load is controlled by an Adaptive ConstrainedMinimum Variance (CMV) controller and it is sampled every second.Though this loop will typically be in regulation mode for large periods of time, this data wasgathered during a period when the reference motor load was changed often. This is because thecontrol action was being tested by moving the motor load so the refiner operated in near the padcollapse region. The reference signal changed often over the period it was collected, so it is notpossible to examine the minimum achievable output variance. However, the ANM will be used onthis data, in order that the static input-output relation can be estimated.The input and output signals appear in Figure 7.85, and the 2700 samples corresponds to aboutforty-five minutes of operation of the refiner, including seventeen reference changes made in thisperiod. It is obvious from the output that this process is quite noisy. The motor load reference andoutput is measured in megawatts (MW) and the plate gap is given as a percentage, where 100%corresponds to the completely closed position. The motor load reference is set to values which arechosen to induce a pad collapse, but the control algorithm successfully moves the refiner operation101Chapter 7: Industrial DataFigure 7.85: The input and output signals of the reject refineraway from a possible collapse. This action provides several cases where the refiner operation mayhave entered the region of slope change, so that the ANM may be able to estimate the slope change.The first step that must be performed is to filter out the noise in the signals. This was done usinga Butterworth low pass filter with a cutoff frequency at 0.05 Hz. The signals which result from thisfiltering appear in Figure 7.86. They have a similar trend as the original signals but have no highfrequency components. Next, the interval for the ANM must be decided on, so that the static input-output relationship can be estimated. First, the interval boundaries were chosen, set so that most ofthe input falls in that range. The interval of [56.6,58.4] was decided upon and the ANM estimationfor this range appears in Figure 7.87(a). The slope of this line indicates the gain of the process overthis interval and it highlights a limitation associated with using a single estimate for the gain of thisprocess. This result suggests that the gain is constant over the interval, however, the ANM will showthis is not the case. The second curve in Figure 7.87(a), using the intervals [56.6,57.2,57.8,58.4]indicates that the slope is relatively constant below 57.8, but seems to change above this value.102Chapter 7: Industrial DataFigure 7.86: The filtered input and output signals of the reject refinerBreaking the interval down into six segments, Figure 7.87(b), gives a better indication of the areawhere the slope change occurs. In Figure 7.87(c), further breaking down of the region above 57.8occurs, where the slope changes. The ANM was also run using the unfiltered data, to see the effectof the noise on the estimation. Figure 7.87(d) compares the results with the single interval result andthough the bias created by the noise is evident, the same shape is there.As this is data from a real process, it is not possible to know if these results are correct, basedonly on the data supplied. However, the results do make sense when they are compared with whatis known about the process. For this process, a change in slope is expected, so the results of ANMestimation are consistent with this knowledge. There is an advantage in knowing where the slopechange occurs, as it can be used to operate the motor load loop. The motor load reference can bechosen based on the results of the estimation, to avoid pad collapse and plate clash. Also, the staticgain curve can be used by the a controller to signal that a reference change is required to avoida plate clash.103Chapter 7: Industrial Data7.3 Kamyr Digester Chip LevelA digester continuously cooks wood chips, the first step on their road to becoming used in paperproducts. The properties of the cooked chips depend on the extent of the reaction in the digester.The most direct effect on this reaction is the residence time of the chips in the digester and the mostcommon method of controlling this time is by manipulating the chip level. This control problem wasexamined by Allison [1] and some of the data collected during this trial is examined. The chip levelwas controlled by manipulating the blow flow using an Adaptive Generalized Predictive Control(AGPC) method. The blow flow is the rate at which cooked pulp is released from the digester.Chip level control is very difficult because it is a time-varying and operating point sensitive process.This is because the chip level is difficult to measure and typically the strain gauges used result in104Chapter 7: Industrial DataFigure 7.88: A typical Kamyr digestermeasurements high in noise and of poor resolution. Also, the chip packing in the digester depends onthe size, density, species and moisture content of the chips which will change on a regular basis. Atypical Kamyr digester appears in Figure 7.88, showing the flow through the digester. Also availableis data collected when the digester was under P control and Constrained Minimum Variance control.All the data was collected every five minutes from the process.As mentioned, the data was collected during periods when different control schemes were beingused. The various schemes will be compared using the ratio of output variance to erm2 . Also, onedata set contains approximately three days of digester operation under AGPC. These data will besplit into four sections to compare the operation of the process over the different intervals.The chip level for the various schemes appears in Figure 7.89, with only the output shown, so thatit can be compared to the reference and the variance of the input signal. Thus it is possible to see howclose the output is to the reference and how active the input was to achieve the results. Figure 7.89(a)presents thirty-three hours of operation under P control where the large dip at sample 140, corresponds105Figure 7.89: The output signals of the digester for (a)P control, (b) AGPC, (c) CMV control and (d) more AGPCto a short manual period in an attempt to bring the output closer to the reference. This was typicalof the control used on the digester, as a PI controller did not provide satisfactory performance andthe operator was required to provide the integral action. The operation was then switched to AGPCcontrol shown in Figure 7.89(b), which lasted for approximately forty hours. The mean of the outputis much closer to the reference, indicating better control. Twenty-five hours of control under CMVcontrol is shown in Figure 7.89(c), which exhibits the worst control. Finally, Figure 7.89(d) presentssome more AGPC control at a later date, about eighty-three hours of operation. The autocorrelationfunction of these signals appears in Figure 7.90, confirming the above observations, that the AGPCscheme appears to provide the best control and CMV controller seems to provide the worst control.106a ii = 259.5100^200^300^400^500^600^700Time (5 minutes)(d)800^900 1000450 500150^200^250^300^350Time (5 minutes)(b)Adaptive GPC Output40050^100^150^200^250^300Time (5 minutes)(c)0P Controller Output^ Adaptive GPC Output70656055si 50-6> 45ry 408 35Y=46.10. 1 = 228.850^100^150^200^250^300^350^400Time (5 minutes)(a)Constrained Minimum Variance Outputy=38.9I^aj= 98.2y=39.470503020100^70 6560 •55 •:3 45 ^40 •C.)3530 •25200 7020Chapter 7: Industrial Data+ - y 2 (7.135)Chapter 7: Industrial DataFigure 7.90: The autoconelation functions of the output signals of the Kamyrdigester for (a) P control, (b) AGPC, (c) CMV control and (d) more AGPCIt should also be noticed that two variances are listed in Figure 7.90, q, the variance of theoutput about its mean, while o is the variance of the output about the reference, that being the meansquare value of the control error. The reason ay must be considered can be explained by examiningFigure 7.90(c). While o-2 for the CMV control is only around seventy, which is close to seventy-six,the variance of the AGPC in Figure 7.90(b), the output response shows that the AGPC is superiorto the CMV. Comparison of the cr; values for these processes, one-hundred and six to seventy-six,gives a better indication of the superior control strategy. This approach penalizes the performanceof a control which may keep the process variance down, but does not stay near the output. The cr!value is quickly calculated using:107Chapter 7: Industrial Datawhere y is the mean of the output and y s is the reference value.Also, an apparent contradiction exists if Figure 7.90(b) and Figure7.90(d) are compared,as Figure 7.90(d) has a lower variance, yet the autocorrelation function dies out quicker inFigure 7.90(b). Based on this information, an inconsistency seems to exist, but it will be explainedafter the minimum achievable output variance of the process is calculated.Unfortunately, it was not possible to estimate the delay of this process. The noise created a biasat one, but even if this was ignored, the data gave an inconsistent delay estimate, ranging betweentwo and four. While this was in the range of likely values [1], there is no guarantee that the estimatedchange was actually a change in the delay and not caused by the disturbance. As a result, the delaysuggested in [1] as occurring most often will be used to determine 1712nv that being two samples.The data were manipulated using various filter numbers and time scales and a choice of ten filtersand a time scale of 0.2 was settled on, as it produced an estimated noise autocorrelation function thatwas closest to an impulse. The results are listed in Table 7.5 and confirm that the AGPC providesthe best control, while the CMV control does not perform near as well. The ratio for the two instantsof AGPC control are very similar, which suggests that the o is lower for the second case becausethe process has changed slightly, not that the controller is doing a better job. The slightly lower ratioof 1.32 for the first instance of AGPC control could also explain the decrease in correlation in lagsfive through ten apparent in its autocorrelation function.A closer look at the data in Table 7.5 confirms the nonstationary behavior of the digester. Eventhough the AGPC and gain control methods share the same reference, variation in the wood chipsentering the digester are evident by examining O-E2 and the F(q-1 ) coefficients. The change in operatingpoint, which occurs during the CMV control, appears to have a large effect on the properties of theController olp,2„, o gym" al t(cri)coefficient cr;2/er;,2,,vP 1.78 125.5 70.8 49.7 0.652 1.26GPC(a) 1.32 76.2 58.0 43.8 0.568 1.29CMV 4.15 107 25.74 16.7 0.737 2.73GPC(b)^_ 1.4 66.2 47.4^_ 37.2^_ 0.524 1.39Table 7.5 The results of minimum achievable output variance estimation108Chapter 7: Industrial Datadigester as the estimated noise standard deviation drops to four from seven. Comparing these figuresto knowledge available during operation of the digester, such as the type of chips used, may enableimprovements in the process performance to be made.Next, the output of Figure 7.89(d) was broken into four equal segments, each representing abouttwenty-one hours of digester operation. The purpose of this is to show how monitoring a 2 can beused to indicate changes in the process. Once these changes are associated with a cause, steps canbe taken to minimize their effect on performance.The control performance, av al, was calculated using ten Laguerre filters and a time scale of0.2. The results appear in Table 7.6. The control performance ratio improves between period oneand two, stays approximately the same for period three, but worsens in period four. Comparison ofthese ratios to the variance performance ratio, 4/61,,, shows that the this ratio stays approximatelythe same, which is encouraging as it suggests that the AGPC enjoys the same success in minimizingthe variance over the three days. However, the control fails to keep the output close to the referenceall the time, and this drift results in an inflated control ratio in period four.Some system identification was performed in order to obtain a time-series model of thedigester [1]. One day of data was used to fit the model with the digester operating at a reference offorty-five percent. The following model was decided upon:Y( t) = fq6-,u(t - 2) + 1— O —1e(t)v.5—2 —0.e2 = 31 F(q-1) = 1 + 0.5q-1These numbers are similar to those estimated in Table 7.5 and Table 7.6, while the time-varyingand operating point dependency of the digester found are to be expected, based on the operationof the digester.Period 4 lerlv cd agiv 'cil f' (q-1) coefficient cr: felv1-250 1.57 77.7 49.4 36.6 0.591 1.43251-500 1.37 61.3 448 34.3 0.552 1.35501-750 1.42 54.7 38.2 28.3 0.589 1.41751-100 2.06 72.1 35.0 28.2 0.4912 1.43Table 7.6 The results of minimum achievable output variance estimation on the AGPC control(7.136)109Chapter 7: Industrial DataUsing the methods presented in earlier chapters and applying the conclusions reached in thesimulation study, it was possible to successfully monitor some industrial processes. Results ofmonitoring the power spectral density function indicated that the performance of the PI controller istwo and a half times worse than is possible, though retuning it may result in better performance. Theapplication of the ANM to the motor load control of a reject refiner resulted in the estimation of itsstatic gain. The knowledge of the slope change is very useful and could be used either in control or inchoosing the motor load reference. The performance of various control algorithms, used to regulatethe chip level in a Kamyr digester, was determined using the monitoring tools. They indicatedwhich performed the best and this was consistent with what was suggested by the autocorrelationfunction of the outputs. The time-varying characteristic of this process was also noted, as the processproperties changed quite often. These examples point out the benefits of monitoring and the extrainsight gained about the process.110Chapter 8: ConclusionChapter 8ConclusionThe important role that control loop monitoring plays in maintaining high process performancecan not be understated. The monitoring methods presented in this thesis are concerned with theregulation control problem so that the output variance can be kept as small as possible. The ideasand results will be summarized in the order that they were presented in this thesis.It is possible to estimate the delay of a stochastic process in regulation by using the Fixed ModelVariable Regressor Estimator. Successful delay identification was obtained for a wide range ofprocesses and knowledge of the time delay can be used to guarantee that the control method correctlycompensates for it. Also, changes in delay over time will often signal the start of deterioration ofthe process operation. The delay estimate can sometimes be biased, which will be caused by thedisturbance or if the input signal lacks sufficient excitation. It is possible to correct the bias caused bythe disturbance, but the bias due to insufficient excitation can not be fixed. However an approximationof the lower limit of the input signal frequency content can be easily found with only the knowledgeof the sampling interval. The insufficient excitation problem is most often encountered when PIcontrollers are used for control.The Adaptive Nonlinear Modeller was successfully applied to closed loop stochastic processesas it was able to estimate their static input-output relationship. The use of Recursive Least Squaresestimation to find the table values was not computationally intensive, even if a large number ofintervals are used. This is because all but two of the entries in the regressor are zero, whichreduces the total number of calculations required. The Adaptive Nonlinear Modeller was able toidentify various nonlinearities including a deadzone and actuator saturation, both which will reducethe process performance. Two steps can be employed to improve the estimation of the static input-output relation. The estimator can be turned off for a short period after setpoint changes, whichreduces the plant dynamics present in the signals. Also, for stochastic processes, the use of a lowpass filter is important to reduce the effect of noise on estimation.IllChapter 8: ConclusionThe final process property determined was the minimum achievable output variance. Theestimation of the minimum achievable output variance involved the successful application of aLaguerre network. The Laguerre gains were correctly estimated for the majority of processesconsidered, using Recursive Extended Least Squares estimation. This justifies the use of this estimatorover other, more complicated ones, even though in some situations estimation may not converge. Thismay occur if the closed loop noise transfer function numerator is not strictly positive real. Simulationsshowed a wide range of Laguerre parameters can be used to obtain high estimation accuracy. If theautocorrelation function of the estimated noise sequence is close to an impulse, then the correspondingchoice of Laguerre parameters can be used with a high degree of confidence. It was also foundthat examination of the estimated noise sequence and Markov coefficients, required to determine theminimum achievable output variance, can be used to indicate if the plant or disturbance of the processchange. Besides being able to use a large range of Laguerre network parameters, a single choice ofparameters were able to model numerous changes in the closed loop noise filter. It should alwaysbe remembered that changing control will not change the minimum achievable output variance. Theminimum achievable output variance will only change with the noise properties of the process.The last stage in this thesis involved applying the monitoring tools to data from industrialprocesses. This stage moved the study from simulation to the real world and the results weresuccessful. Much useful information about the processes considered was provided by monitoring.The results show how monitoring can be used to compare various control algorithms, by determiningwhich provides the best reduction in the output variance. Also, once the minimum variance level isknown, decisions can be made as to the best approach to try in order that the process output variancecan be reduced. If the minimum achievable output variance is not low enough, further reduction canbe made by either retuning or replacing the control algorithm or by modifying the process. Also,because changes in the process plant and disturbance are detected during estimation, this knowledgecan be used to improve process performance. Lastly, this industrial data shows that PI controllersand high noise in processes can cause problems in the estimation of the time delay. The applicationof monitoring to industrial processes will result in improved process performance.There are several possible avenues envisioned that would continue this work. A very useful112Chapter 8: Conclusionexercise would involve the transportation of the monitoring tools to a distributed control system toenable the monitoring of an industrial process for a long period of time, such as days or tens ofthousands of samples. This would confirm the benefits derived from monitoring. Also, the minimumachievable output variance could be determined for other control structures, including the case wherefeedforward control is used. It was found that the number of filters in the Laguerre network hadto be increased to maintain good estimation for closed loop noise transfer functions with complexpoles. It is possible to use a Laguerre network which employs a complex time scale [24], whichshould be done if the extra computation required by the network size is a concern. It would also bepossible to use the estimated process properties with a knowledge based system. Decisions could bemade based on the results obtained from monitoring the process. Finally, attention could be focusedon the possibility of monitoring more process properties. It would be useful to have an idea of theservo performance of the process and to use this in conjunction with the variance performance toobtain a controller that satisfies both criteria. This would be very helpful when applied to cascadeloops. The inner loop in this configuration generally performs a servo role while the outer loopacts as a regulator. The availability of an indicator of the overall performance of the cascade loopwould be very beneficial.113References[1] B. J. Allison, G. A. Dumont, L. H. Novak, and W. J. Cheetham. Adaptive-Predictive Controlof Kamyr Digester Chip Level. AIChE Journal, 36(7):1075-1086, July, 1990.[2] K. J. Astrom. Introduction to Stochastic Control Theory. Academic Press, Inc., New York,New York, 1970.[3] K. J. Astrom. The Adaptive Nonlinear Modeller. Technical Report TFRT-3178, Department ofAutomatic Control Lund Institute of Technology, 1985.[4] K. J. Astrom. Assessment of Achievable Performance of Simple Feedback Loops. InternationalJournal of Adaptive Control and Signal Processing, 5(1):3-19, 1991.[5] K. J. Astrom and B. Wittenmark. Computer Controlled Systems. Prentice—Hall, New Jersey,1984.[6] K. J. Astrom and B. Wittenmark. Adaptive Control. Addison-Wesley Publishing Company,Don Mills, Ontario, 1989.[7] W. L. Bialkowski. Bleach Plant Process Control. In Tappi Bleach Plant Operations ShortCourse, Hilton Head, S.C., June, 1990.[8] W. L. Bialkowski. Dreams vs. Reality: A View from Both Sides of the Gap. In Control Systems'92: Dream vs. Reality, Whistler, B.C., September, 1992.[9] G. Box and G. Jenkins. Time Series Analysis. Holden-Day, San Francisco, California, 1976.[10] D. W. Clarke, C. Mohtadi, and P. S. Tuffs. Generalized Predictive Control Parts I and II.Automatica, 23(2):137-160, 1987.[11] G. A. Dumont and K. J. Astrom. Wood Chip Refiner Control. IEEE Control Systems Magazine,pages 38-43, April, 1988.[12] A. Elnaggar. Variable Regression Estimation of Unknown System Delay. PhD thesis, Universityof British Columbia, 1990.[13] A. Elnaggar. Direct Delay Estimation in Stochastic Closed Loop Systems. Submitted to 12thWorld Congress IFAC for publication, 1992.114[14] G. Ferretti, C. Maffezzoni, and R. Scattolini. Recursive Estimation of Time Delay in SampledSystems. Automatica, 27(4):653-661, 1991.[15] Y. Fu and G. A. Dumont. Optimum Laguerre Time Scale and Its On-line Estimation. IEEETransactions on Automatic Control, publication February 1993.[16] S. Gunnarsson and B. Wahlberg. Some Asymptotic Results in Recursive Identification UsingLaguerre Models. International Journal of Adaptive Control and Signal Processing, 5:313-333, 1991.[17] T. J. Harris. Assessment of Control Loop Performance. The Canadian Journal of ChemicalEngineering, 67:857-862, October, 1988.[18] T. J. Harris and J. F. MacGregor. An Overview of Discrete Stochastic Controllers:GeneralizedPID Algorithms with Dead-Time Compensation. The Canadian Journal of ChemicalEngineering, 60:425-432, June, 1982.[19] R. Isermann. Digital Control Systems. Springer—Verlag, Berlin, Germany, 1981.[20] L. Ljung. System Identification. Prentice—Hall, New Jersey, 1987.[21] L. Ljung and T. SiiderstrOm. Theory and Practice of Recursive Identification. MIT Press,Cambridge, Massachusetts, 1983.[22] J. P. Norton. An Introduction to Identification. Academic Press, Inc., New York, New York,1986.[23] M. Perrier and A. A. Roche. Towards Mill-Wide Evaluation of Control Loop Performance. InControl Systems '92: Dream vs. Reality, Whistler, B.C., September, 1992.[24] B. Wahlberg. Identification of Resonant Systems Using Kautz Filters. In Proceedings of the30th IEEE Conference on Decision and Control, Brighton, England, December, 1991.[25] C. C. Zervos. Adaptive Control Based on Orthonormal Series Representation. PhD thesis,University of British Columbia, 1988.[26] C. C. Zervos and G. A. Dumont. Deterministic Adaptive Control based on Laguerre SeriesRepresentation. International Journal of Control, 48(6):2333-2359, 1988.115
- Library Home /
- Search Collections /
- Open Collections /
- Browse Collections /
- UBC Theses and Dissertations /
- Control loop performance
Open Collections
UBC Theses and Dissertations
Featured Collection
UBC Theses and Dissertations
Control loop performance Lynch, Christopher 1992-09-30
pdf
Page Metadata
Item Metadata
Title | Control loop performance |
Creator |
Lynch, Christopher |
Date Issued | 1992 |
Description | Control loop monitoring provides the opportunity to maintain high process performance. In this thesis, the regulation performance of a process is monitored by examining the time delay, static input-output relation, and the minimum achievable output variance of the process. The time delay estimate is obtained from the Fixed Model Variable Regressor Estimator and the static input-output relation is estimated using the Adaptive Nonlinear Modeller. The minimum achievable output variance of the process is determined by modelling the closed loop noise transfer function by a Laguerre network and finding the Laguerre gains of the network using Recursive Extended Least Squares estimation. An extensive simulation study was performed to examine the operation, usefulness and limitations of the monitoring tools. These simulations confirmed the operation and benefit of the monitoring methods, and defined their limitations. The monitoring tools were applied to data generated by industrial processes. The successful application to industrial processes highlighted the benefits obtainable by control loop monitoring. Due to the ease of applying the monitoring tools and the valuable insight they provide, these monitoring tools should be applied to any process where regulation performance is a concern. |
Extent | 5284960 bytes |
Genre |
Thesis/Dissertation |
Type |
Text |
File Format | application/pdf |
Language | eng |
Date Available | 2008-09-30 |
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.0086131 |
URI | http://hdl.handle.net/2429/2419 |
Degree |
Master of Applied Science - MASc |
Program |
Electrical Engineering |
Affiliation |
Applied Science, Faculty of Electrical and Computer Engineering, Department of |
Degree Grantor | University of British Columbia |
Graduation Date | 1993-05 |
Campus |
UBCV |
Scholarly Level | Graduate |
Aggregated Source Repository | DSpace |
Download
- Media
- 831-ubc_1993_spring_lynch_christopher.pdf [ 5.04MB ]
- Metadata
- JSON: 831-1.0086131.json
- JSON-LD: 831-1.0086131-ld.json
- RDF/XML (Pretty): 831-1.0086131-rdf.xml
- RDF/JSON: 831-1.0086131-rdf.json
- Turtle: 831-1.0086131-turtle.txt
- N-Triples: 831-1.0086131-rdf-ntriples.txt
- Original Record: 831-1.0086131-source.json
- Full Text
- 831-1.0086131-fulltext.txt
- Citation
- 831-1.0086131.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-0086131/manifest