UBC Theses and Dissertations

UBC Theses Logo

UBC Theses and Dissertations

Dual adaptive control of chip refiner motor load Allison, Bruce James 1994

Your browser doesn't seem to have a PDF viewer, please download the PDF to view this item.

Item Metadata

Download

Media
831-ubc_1994-893249.pdf [ 6.42MB ]
Metadata
JSON: 831-1.0058582.json
JSON-LD: 831-1.0058582-ld.json
RDF/XML (Pretty): 831-1.0058582-rdf.xml
RDF/JSON: 831-1.0058582-rdf.json
Turtle: 831-1.0058582-turtle.txt
N-Triples: 831-1.0058582-rdf-ntriples.txt
Original Record: 831-1.0058582-source.json
Full Text
831-1.0058582-fulltext.txt
Citation
831-1.0058582.ris

Full Text

DUAL ADAPTIVE CONTROL OF CHIP REFINER MOTOR LOAD by BRUCE JAMES ALLISON B.Eng., McMaster University, 1983 M.Eng., McMaster University, 1986 A THESIS SUBMITTED IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF DOCTOR OF PHILOSOPHY  in THE FACULTY OF GRADUATE STUDIES DEPARTMENT OF CHEMICAL ENGINEERING  We accept this thesis as conforming to the required standard  THE UNIVERSITY OF BRITISH COLUMBIA January 1994 © Bruce James Allison, 1994  In presenting this thesis in partial fulfilment of the requirements for an advanced degree at the University of British Columbia, I agree that the Library shall make it freely available for reference and study. I further agree that permission for extensive copying of the thesis for scholarly purposes may be granted by the head of my department or by his or her representatives. It is understood that copying or publication of this thesis for financial gain shall not be allowed without my written permission.  Bruce J. Allison Department of Chemical Engineering The University of British Columbia 2216 Main Mall Vancouver, B.C. Canada V6T 1Z4  Date:  V= &, . \Q V»ft\  ABSTRACT  A characteristic of wood chip refiners is that the incremental gain between the motor load and the plate gap is subject to a slow drift due to plate wear and sudden changes in sign due to pulp pad collapse. A pad collapse can be caused by a change in operating point, or may occur suddenly due to a feed rate or consistency disturbance. This poses a problem forfixed-parameterlinear controllers which may actually accelerate pad collapse and induce plate clashing as a result of getting caught in a positive feedback loop. The objective of this thesis is to develop a reliable chip refiner motor load controller and to test it out on an industrial refiner. The problem is approached from a fault detection and control viewpoint and the proposed algorithm consists of an improved parameter estimator and a controller containing "dual" features. The role of the improved estimator is to track both drifts and jumps in the parameters. The use of active learning or probing in the controller is justified by the fact that the parameter estimates are key to identifying a pad collapse, and that probing targets a portion of the input energy at continuously identifying these parameters. Since there still does not exist a general dual controller design methodology, the main challenge was to extend existing suboptimal approaches to handle realistic dynamics including deadtime and correlated disturbances. To track both slow and fast changes in the system parameters, a multi-model approach called adaptive forgetting through multiple models or AFMM is used. A method of modifying the AFMM to include information about what to expect in the event of a pad collapse is proposed. The main contribution of the thesis is the development of the active adaptive controller or AAC, which consists  ii  of a constrained certainty equivalence approach coupled with an extended output horizon to deal with nonminimum phase systems and a cost function extension to get probing. Finally, the AAC is combined with the AFMM, and the resulting combination is tested on an industrial refiner.  iii  TABLE OF CONTENTS Page Abstract  ii  List of Tables  viii  List of Figures  ix  Acknowledgement  xii  1.  Introduction  1  1.1  General Introduction  1  1.2  Background and Literature Review  2  1.2.1  Mechanical Pulping Process  2  1.2.1.1 Unit Operations  3  1.2.1.2 Pulp Characterization  4  1.2.1.3 Chip Refining  6  1.2.1.4 Theory of Refining  9  1.2.2  1.3  1.2.1.5 Motor Load Plate Gap Relationship  12  Control of Mechanical Pulping  15  1.2.2.1 Refiner Control  16  1.2.2.2 Motor Load Control  18  The Thesis  20  1.3.1  21  Objectives  iv  Page 1.3.2  Contributions  22  1.3.3  Overview  23  Definitions  25  2.1  Introduction  25  2.2  Definitions  26  Comparison of Adaptive Controllers  29  3.1  Introduction  29  3.2  Problem Definition  33  3.3  Estimators and Controllers  34  3.3.1  Passive Adaptive Controllers  35  3.3.1.1 Recursive Maximum Likelihood (RML) Estimator  35  3.3.1.2 RML-Based Constrained Certainty Equivalence (CCE) Controller  36  3.3.1.3 Extended Kalman Filter (EKF)  38  3.3.1.4 Cautious Controller  41  3.3.1.5 Innovations Dual Controller (IDC)  43  3.3.1.6 EKF-Based Constrained Certainty Equivalence (CCE) Controller 3.3.2  44  Active Adaptive Controllers  45  3.3.2.1 Square wave Perturbation (PERT)  46  v  Page 3.3.2.2 Input-Constrained One-Step Minimization (ICM)  46  3.3.2.3 Hard Covariance-Constrained One-Step Minimization (HCCM)  46  3.3.2.4 Soft Covariance-Constrained One-Step Minimization (SCCM,) 3.4  3.5  49  Results and Discussion  54  3.4.1  Passive Adaptive Controllers  57  3.4.2  Active Adaptive Controllers  63  Summary and Conclusions  70  Active Adaptive Controller  73  4.1  Introduction  73  4.2  Problem Definition  74  4.3  Estimator and Controller  74  4.3.1  Recursive Maximum Likelihood (RML) Estimator  74  4.3.2  Active Adaptive Controller (AAC)  75  4.4  Results and Discussion  81  4.5  Summary and Conclusions  90  Application to the Refiner  91  5.1  Introduction  91  5.2  Process Description  92 vi  Page 5.3  6.  Estimator and Controller  93  5.3.1  Adaptive Forgetting through Multiple Models (AFMM)  93  5.3.2  Active Adaptive Controller (AAC)  100  5.4  Results and Discussion  102  5.5  Summary and Conclusions  118  General Summary and Conclusions  122  Nomenclature  126  Bibliography  133  Appendices  140  1.  One-Step Optimal Controller  140  2.  Equivalence Between CCE Controllers  142  3.  Expression for Future Covariance  143  4.  Variance of Open-Loop Gain Estimate  145  vii  LIST OF TABLES Page Table 1.1  Pulp and handsheet and/or paper properties.  Table 3.1  Conditions corresponding to results in Figure 3.6.  61  Table 3.2  Conditions corresponding to results in Figure 3.7.  63  Table 3.3  Conditions corresponding to results in Figure 3.11.  69  Table 4.1  Schedule of parameter changes for simulated system.  81  Table 5.1  Results of preliminary model identification experiments.  viii  5  102  LIST OF FIGURES Page Figure 1.1  Typical chemi-thermomechanical pulping (CTMP) plant.  4  Figure 1.2  Simplified diagram of a chip refiner.  7  Figure 1.3  A typical first-stage refiner plate showing bar-groove pattern.  8  Figure 1.4  Idealized operating curve showing static gain between plate gap and motor load.  13  Figure 1.5  Operating curve for the industrial refiner.  13  Figure 3.1  Block diagram of the optimal stochastic regulator.  30  Figure 3.2  Minimization of cost function when |Vw*| £ |VMJ.  52  Figure 3.3  Minimization of cost function when |Vw*| < |VuJ.  52  Figure 3.4  Histogram of output variances from 1000 runs of 1000 iterations each of the cautious controller.  56  Figure 3.5  Comparison of output variances for passive adaptive controllers.  57  Figure 3.6  Effect of parameter transition function mismatch on the passive adaptive controllers. Speed of true parameter transition functions overestimated.  Figure 3.7  60  Effea of parameter transition function mismatch on the passive adaptive controllers. Speed of true parameter transition functions  Figure 3.8  underestimated.  62  Comparison of output variances for active adaptive controllers.  64  ix  Page Figure 3.9a  A realization of the CCE-SCCM combination. The dashed lines are the true parameters and the solid lines the estimates.  65  Figure 3.9b  Continuation of Figure 3.9a.  66  Figure 3.10  Comparison of probing signals for a typical realization.  67  Figure 3.11  Effect of parameter transition function mismatch on the active adaptive controllers.  68  Figure 3.12  Variances for the CCE-SCCMj combination.  71  Figure 4.1  Control without probing (A, = 0).  83  Figure 4.2  Control with probing (X = 10).  83  Figure 4.3a  Parameter estimates for probing disabled case (Figure 4.1). The dashed lines are the true parameters and the solid lines the estimates.  84  Figure 4.3b  Continuation of Figure 4.3a.  85  Figure 4.4a  Parameter estimates for probing enabled case (Figure 4.2). The dashed lines are the true parameters and the solid lines the estimates.  86  Figure 4.4b  Continuation of Figure 4.4a.  87  Figure 5.1  Schematic diagram of hydraulic plate positioning system for a Bauer refiner using LVDT position measurement.  92  Figure 5.2  Schematic diagram of AFMM with three models running.  94  Figure 5.3  First method of dealing with input multiplicity.  101  Figure 5.4  Second method of dealing with input multiplicity.  101  x  Page Figure 5.5  Time series record from open-loop identification experiment 3 conducted immediately after a plate change.  Figure 5.6a  103  Preliminary trial of motor load controller using control strategy shown in Figure 5.3.  105  Figure 5.6b  Continuation of Figure 5.6a.  106  Figure 5.7a  Continuation of preliminary trial.  108  Figure 5.7b  Continuation of Figure 5.7a.  109  Figure 5.8a  Trial of motor load controller using control strategy shown in Figure 5.4 and modified tuning parameters.  112  Figure 5.8b  Continuation of Figure 5.8a.  113  Figure 5.9a  Trial of motor load controller showing automatic setpoint adjustment or "reset" when the setpoint exceeds the maximum load.  116  Figure 5.9b  Continuation of Figure 5.9a.  117  Figure 5.10a  Initialization phase of motor load controller.  119  Figure 5.10b  Continuation of Figure 5.10a.  120  xi  ACKNOWLEDGMENT  An undertaking like this can never be brought to fruition without the financial, technical and motivational help and support of many. I would like to take the time now to acknowledge those people. To the Pulp and Paper Research Institute of Canada, and Jim Rogers in particular, thank-you for giving me this opportunity. It simply would not have been possible without your support. Also, thanks to the Natural Sciences and Engineering Research Council of Canada for providing financial support. To my research supervisors Patrick Tessier and Guy Dumont, thank-you for your guidance, encouragement and, most importantly, your friendship. Patrick, thanks for your technical and administrative support and for your sense of humour in the difficult times. Good luck in your future endeavors. Guy, thanks for nearly eight very fruitful years of collaboration and for playing a major role in shaping my career to this point. Your unwavering confidence has been a great example to me. To Joe Ciarniello and others at MacMillan Bloedel, thank-you for your support of this project and for letting me use your mill as a laboratory. Joe, thanks for being there every step of the way and for putting your trust in me. To my co-student, co-worker and good friend Robert Gooding, thank-you for the many hours of enlightening discussion, both technical and otherwise. Your going through this at the same time has made this much easier for me and I am sure that one day I will look upon these as "the best years of my life". To my immediate family, my mom, brothers, sister-in-laws, nieces and nephews, thank-you for showing interest and providing encouragement. To my dad, who faithfully upheld the virtues of education, my one regret it that you can't be here to share this with me. Finally, to my wife Beth, two year old daughter Alex and baby number two (expected February 1994), thank-you for your patience, understanding, love and support, but most of all for keeping the kid in me alive. Yes, this means my weekends are free from now on!  xii  CHAPTER 1 INTRODUCTION  1.1 General Introduction The are two main wood pulping processes - chemical and mechanical - each producing fibres with substantially different characteristics and end use applications. Chemical pulping dissolves the lignin which binds the fibres together. The resulting pulps are composed of long, flexible, strong fibres which are easily whitened in the bleaching process. Unbleached chemical pulps are used in the production of linerboard, wrapper and paper bag, etc. - products which require high strength. Because of their high stable brightness, bleached chemical pulps are used in the production of high-grade printing and writing papers, and to add strength in newsprint and magazine grade papers. The main disadvantage of the chemical pulping process is its low yield (40-45%) - one tonne of wood yields only about one-half tonne of usable fibre. Thus, the process is not very resource conserving. Mechanical pulping can provide an alternative to the chemical pulping process. Mechanical pulping fractures the lignin binding the wood fibres together, resulting in pulps composed of fibre bundles andfibrefragmentswith some whole individual fibres. More wood constituents are retained resulting in higher yields (90-95%). On the negative side the problems are that; (i) papers made from mechanical pulps do not at present have the physical characteristics necessary for high quality printing, (ii) their colour is not stable, i.e. they yellow quickly, and (iii) for a variety of reasons which have to do with a general lack of adequate instrumentation and control, it is difficult to obtain uniformly high  1  Chapter 1. Introduction  2  quality pulp at the desired target properties. These deficiencies currently restrict the use of mechanical pulps to newsprint and other low-cost papers. A positive effect on the environment would be achieved if mechanical techniques, coupled with some chemical treatment, could be improved to point where the physical and chemical properties of the pulps are acceptable for higher grade applications. This thesis partially addresses problem (iii) above, i.e. that of obtaining uniformly high quality pulp at the desired target properties. In particular, a chip refiner motor load control algorithm is developed and tested on an industrial refiner. The purpose of this chapter is to introduce the problem and to state the thesis objectives. The chapter is organized as follows. Section 1.2 gives the relevant background material. In Section 1.2.1, the mechanical pulping process is introduced. It consists of a number of unit operations, the most important being chip refining. Chip refiner operation is described in detail and the importance of operating conditions on refining action and pulp properties is discussed. Finally, the complex nonlinear relationship which exists between the motor load and the plate gap is explained. In Section 1.2.2, the literature on mechanical pulping process control is critically reviewed. Chip refiner control has received by far the most attention. The plethora of refiner control strategies is reviewed in the light of current ideas about pulp characterization and refining theory. Finally, Section 1.3 gives the thesis objectives and contributions, and an overview of the thesis.  1.2 Background and Literature Review 1.2.1 Mechanical Pulping Process Strictly speaking, mechanical pulps are produced by two different processes - grinding and refining. Grinding was the first industrially successful way of using wood as a papermaking fibre, and its history dates back to the 18th century. In 1960, virtually all mechanical pulp was produced  Chapter 1. Introduction  3  by the grinding process and, as recently as 1975, grinding still accounted for 90% of North American mechanical pulp production. However, by 1980 approximately 50% of mechanical pulp was being produced by the refiner method. This was due, in part, to a better understanding of pressurized refiner mechanical pulping (Atack, 1972) and to the development of heat recovery systems (Huusari and Syrjanen, 1980) which reduced the net energy requirement to approximately the same level as that for groundwood. This paved the way for a rapid changeover to refiner-based methods, which have the following advantages over grinding; 1. utilization of sawmill residuals as a feedstock 2. higher pulp strength 3. reduced labour cost  1.2.1.1 Unit Operations A substantial amount of information is available (Casey, 1980; Smook, 1982; Leask, 1987) describing the sequence of unit operations involved in producing refiner mechanical and chemimechanical pulps. Generally, the process takes place in three main parts; (i) chip pretreatment, (ii) chip refining, and (iii) pulp processing, as shown in Figure 1.1. Briefly, the unit operations involved in chip pretreatment consists of: screening to remove under and oversize material; washing to remove contaminants (dirt, sand, grit, etc. - basically anything that could cause undue wear of refiner plates); and steaming to soften the lignin binding the fibres together and to drive off entrained air which is detrimental to heat recovery. In chemimechanical pulping, the wood chips are impregnated with chemicals to improve pulp brightness and strength, and to reduce energy requirements. Following chip refining (discussed in Section 1.2.1.3), the mechanical pulping process  4  Chapter J. Introduction  comprises a series of unit operations aimed at enhancing and controlling pulp quality. These are: latency removal to straighten fibres; screening and reject refining to remove debris (unrefined fibre bundles) while flexibilizing the long fibres; cleaning to remove contaminants; washing to remove wood resins and metallic ions; and bleaching to increase brightness.  1.2.1.2 Pulp Characterization Pulp quality can really only be defined with respect to its intended end use. Most mechanical pulps are used to make paper and the goal in papermaking is to produce a uniform sheet with high strength and good optical properties. However, it is not a simple matter to relate paper properties to  RECHIPPER  CHEMICALS  OVER!SIZE CHIPS  ' FINES  CHIP  CHIP  SCREEN  WASHER  H  CHIP PRE-  .  CHIP IMPREGNATION  STEAMING  "—*  CHIP  }  PRETREATMENT  }  PROCESSING  ^  LATENCY REMOVAL  CHIP STEAMING  i  i  PULP SCREEN  '  PRIMARY  SECONDARY]  REFINER  REFINER  PULP CLEANER  |  1 1 f  CHIP REFINING  PULP WASHER  1 — SAND, DIRT,  1  REJECT  GRIT  REFINER  ' SCREEN  RE-JECTS ACCEPTS  Figure 1.1: Typical chemi-thermomechanical pulping (CTMP) plant.  PULP BLEACHING  PULP  PULP  5  Chapter 1. Introduction  the properties of pulp. For example, some commonly used mechanical pulp and handsheet and/or paper properties are listed below in Table 1.1. Current industrial practice is largely based on characterizing mechanical pulps by a single Canadian Standard Freeness or CSF measurement. CSF measures the drainage characteristics of pulp and, therefore, is an important indicator of how well the pulp will perform on the paper machine. Moreover, the specific surface (bonding potential) of mechanical pulp has a high inverse correlation with freeness. The adverse effect of excessive freeness variability on papermaking and printing is well known. For example, when pulp freeness is consistently below the target value, slower machine speeds and reduced paper production are required on drainage capacity limited paper machines. Web break frequency increases when the pulp freeness consistently exceeds the target value (Hill et al., 1979), and when the pulp shive content is high (Laurila et al., 1978). For newsprint grade mechanical  Mechanical Pulp Properties  Handsheet and/or Paper Properties  freeness latency average fibre length fibre length distribution shive content long fibre content specific surface area fibre fibrillation fibre coarseness fibre flexibility fibre curl fibre kink  burst strength tear strength tensile strength fracture toughness internal bond strength stiffness fold stretch brightness opacity light scattering coefficient colour roughness/smoothness formation air permeability / porosity thickness / sheet density  Table 1.1: Pulp and handsheet and/or paper properties.  Chapter 1. Introduction  6  pulps, high shive content causes liming in printing press equipment, resulting in lost production and lower print quality (Wood and Karnis, 1977). Forgacs (1963) showed that many of the physical properties of mechanical pulps which are important to their papermaking potential (see Table 1.1) can be accurately predicted from only two measurements; afibrelength factor and a shape (specific surface or bonding) factor. Recently, several groups (Strand et al., 1989; Nobleza, 1991; Guan et al., 1991) have tried to relate the properties of paper to those of mechanical pulp by analyzing large sets of mill data using multivariate statistical analysis tools such as factor analysis, canonical correlation and partial least squares (PLS). All their results show that at least two groups of paper properties can be clearly separated, of which one is related to fibre length and the other to specific surface. From this, it is probably safe to say that mechanical pulps are characterized by at least two parameters.  1.2.1.3 Chip Refining The wood chip refiner (Figure 1.2) is the heart of the mechanical pulping operation. Wood chips are introduced by a metering screw into the open eye of the refiner. As they move through the refining zone between the revolving plates (or discs), the chips are progressively broken down into smaller particles, andfinallyinto individual fibres. Dilution water is supplied to the eye of the refiner to control pulp consistency, which is defined as the mass of dry fibre divided by the total mass of pulp. The plate gap is of critical importance, and is controlled by either an electro-mechanical or hydraulic loading system. The plates are usually tapered slightly, and the peripheral plate separation is in the range 0.1 to 0.5 mm, depending on the refining stage. A plate gap corresponding to the lower end of this range is on the order of three to four uncompressed softwood fibre diameters. Refiners are currently available with up to 24 MW of installed power capacity and plate diameters  7  Chapter 1. Introduction  close to 2 m (Leask, 1993). They come in single-disc (one revolving disc and one stationary disc), double-disc (two discs revolving in opposite directions), and twin (a revolving double-sided disc between two stationary discs) configurations. Depending on the configuration, the discs rotate at either 1200 or 1800 r/mia Chip refining is usually carried out in two stages. First-stage (primary) refining of chips may occur at atmospheric pressure or above (note, this refers to the casing pressure, and not the pressure developed between the refiner plates which is usually much higher). Primary refiner plates (Figure 1.3) are divided into three distinct sections; (i) the breaker bar section, (ii) the intermediate bar section, and (iii) the fine bar section. The plates may have dams to encourage recirculation of fibres. The  WOOD CHIPS  FEED END MOTOR  PULP & STEAM  CONTROL END MOTOR  HYDRAULIC CYLINDER ~"  DILUTION WATER PLATE POSITION  Figure 1.2: Simplified diagram of a chip refiner.  Chapter 1. Introduction  8  breaker bar section has wide bars with deep grooves and is closest to the eye of the refiner. It serves to break up the chips as they enter. The intermediate and fine bar sections consist of progressively narrower bars and shallower grooves. According to observations made on a single-disc, atmospheric discharge refiner using high speed cin6 films (Atack et al., 1984), the breaker bars shred chips into coarse pulp. This pulp proceeds to the refining zone where some of it recirculates back to the breaker bars along the grooves of the stationary plate, while the rest moves radially outward. For a short fraction of its passage through the refiner, most of the fibrous material is constrained to move in the direction of rotation of the moving plate. Some of this material is stapled momentarily in a tangential orientation across the  Figure 1.3: A typicalfirst-stagerefiner plate showing bar-groove pattern.  Chapter 1. Introduction  9  bars of both sets of plates. The immobilized fibres are then subjected to refining action between the relatively moving bars before being disgorged into the adjacent grooves. Second-stage (secondary) refining also takes place at atmospheric pressure or above. However, the plates have a shorter breaker bar section and a larger refining zone. The breaker bars align and impart centrifugal force to the partly refined pulp. Generally, the energy input is lower and the plate gap is smaller than in the primary refiner. Refining is accompanied by the generation of copious amounts of steam. Some of it flows backwards against the feed of chips and is used in the steaming process. The rest flows out with the pulp and is separated in a cyclone separator. This also cools the pulp to help prevent brightness loss. A carefully designed heat recovery system can reduce energy costs by 35-40% (Leask, 1987). Atmospheric steam recovery systems employ tubular heat exchangers to preheat boiler fed water and shower water for paper machines, and air for drying or building heating. High pressure steam from pressurized refining is converted to clean steam in heat exchanges and used in paper machine dryer sections.  1.2.1.4 Theory of Refining Chip refiners play a major role in determining the final physical and papermaking properties of the pulp, as well as its uniformity. As in many comminution processes, the extent of particle size reduction is dependent upon the applied specific energy (E), defined as the net power applied (P), i.e. the difference between the total power and the backed off power, divided by the dry pulp mass throughput (F)  where the units are; E (energy/mass), P (energy/time), and F (mass/time).  10  Chapter 1. Introduction  The effect of specific energy in determining the properties of mechanical pulp, such as freeness, is well known. However, there is an overwhelming amount of empirical evidence (Breck et al., 1975; Stationwala et al., 1979; Franz6n and Sweitzer, 1982; Strand and Harder, 1985; Croteau et al., 1990, Duever et al., 1991) that specific energy alone does not determine refining conditions and pulp properties. Refining consistency also has a major effect. Fortunately, knowledge of the refining process has improved remarkably over the last 10 years. There now exists an emerging theory of refining (Atack, 1980; Miles et al., 1980; Pearson, 1990; Miles and May, 1990) which, at least partially, explains the specific energy and consistency effects. According to theory, refining action can be described by two variables; the total number (N) of refiner bar impacts on a unit of pulp during its passage through the refiner and the average specific energy or intensity (/) of each impact. Together, these two parameters account for the expenditure of specific energy, according to the following equation E = N'l  (!-2)  where the units are; N (impacts), and / (energy/mass/impact). To be useful, N and / must be related to key refining variables, which are: pulp mass throughput; net power applied; and refining consistency. N may be determined from the rotational speed of the refiner, the number of bars per unit length of arc, the plate radius, and the residence time of pulp in the refining zone (Miles and May, 1990). However, because refining is done at high consistency, there is insufficient water to form a liquid suspension of pulp in the refiner. Therefore, the pulp residence time cannot be determined by calculating the flow using volumetric methods. Instead, it depends in a complicated way on the mechanical and aerodynamic forces exerted on the pulp by the refiner and the steam generated in it. Clearly, another measurement is required to determine N and / from operating conditions.  Chapter 1. Introduction  11  Nevertheless, the theory is useful for explaining some well known, but not generally well understood refining effects. It is a well documented fact that increased refining consistency results in mechanical pulps of improved strength at reduced levels of specific energy, but at the possible expense of optical properties (Atack and Wood, 1973; Ryti and Manner, 1977; Manner et al., 1986; Miles and May, 1990). This is explained as follows. Increasing refining consistency decreases the wet mass of pulp in the refiner. This reduces the centrifugal force acting on the pulp and, thereby, increases the residence time (Miles and May, 1990). Increased residence time is accompanied by increased plate gap; there is more pulp between the refiner plates and therefore less axial thrust is required to develop the same motor load. If the motor load, and thus the specific energy, is kept constant by reducing the axial thrust, the result is a greater number of low intensity impacts, which preserves fibre length and improves pulp strength, but at the expense of optical properties. Conversely, decreased refining consistency at constant specific energy results in decreased plate gap and a smaller number of impacts but at higher intensity. Refining under these conditions produces pulps with lower strength, but with better optical properties, lower rejects, and lower freeness (Strand and Hartler, 1985; Miles et al., 1990). At extremely small plate gap, excessive fibre shortening occurs. Recent work by Miles et al. (1990) suggests that for a constant total specific energy, the combination of high intensity-low energy refining in thefirst-stage,followed by low intensity-high energy refining in the second-stage, could produce pulps with good printing properties and with no loss in strength. Furthermore, Stationwala et al. (1993) have shown that refining energy can be reduced by roughly 25% by increasing the rotational speed of the primary refiner. The apparent reason for this is that rotational speed affects bar impact frequency, whereas consistency does not.  Chapter 1. Introduction  12  Because wood is a viscoelastic material, pulp properties depend not only on the amount of energy absorbed at impact, but also on the rate at which the energy is transferred. However, impact frequency cannot be considered a manipulated variable for refiner control as the rotational speed of industrial refiners is fixed. The main point is that the correct combination of two variables (N and /) , and not just one (£), is required to determine refining action and, therefore, the pulp properties for a given raw material. This agrees with the earlier discussion of mechanical pulp characterization (see Section 1.2.1.2). Although it is difficult to determine N and / from operating conditions, the theory does help to explain the important effects of specific energy and refining consistency.  1.2.1.5 Motor Load Plate Gap Relationship At constant throughput, the net power applied is the major factor determining the applied specific energy. The net power or load drawn by the refiner motor is manipulated by changing the hydraulic pressure, which in turn changes the gap between the refiner plates. An idealized motor load plate gap curve is shown in Figure 1.4. The nonlinearity may be explained as follows. As discussed earlier, inside a refiner wood fibres are separated from fibre networks and lined up tangentially along metal bars by passing bars on the opposite refiner plate. Fibres are worked upon by these passing bars. Decreasing (closing) the gap between the bars increases the degree of work expended on the fibres. This increases the refining action. However, eventually a point is reached at which the net power to the refiner passes through a well defined maximum. The maximum is commonly thought to be the point at which a pulp pad collapses. The mechanism by which this occurs is unknown, but a possible explanation is that at this point refiner bars "cut" fibres rather than "work" on them. Another possible explanation is that, at elevated loads,  13  Chapter 1. Introduction  47.5  48  48.5  49  49.5  50  50.5  51  51.5  52  52.5  Plate Position, u (%)  Figure 1.4: Idealized operating curve showing static gain between plate gap and motor load.  47.5  48  48.5  49  49.5  50  50.5  51  51.5  Plate Position, u (%)  Figure 1.5: Operating curve for the industrial refiner.  52  52.5  Chapter 1. Introduction  14  increased steam production contributes to feeding problems associated with large quantities of "blowback" steam. Whatever the reason, once this happens, there are no longer fibres between the gap to support the load, so the motor load drops dramatically and eventually the refiner plates break through the pulp pad. This causes metal to metal contact between the plates (a plate clash) resulting in excessive plate wear and operating disruptions. During refiner operation, the control objective is to regulate the motor load while operating solely in the linear region of the operating space to the left of the dashed line in Figure 1.4, i.e. in the "normal" operating regioa To avoid plate clashes, a controller should never attempt to operate in the "pad collapse" region. The hysteresis in Figure 1.4 reflects the characteristic that, should a controller encounter a pad collapse, it will need to open the plates beyond the point of collapse to rebuild the pulp pad. Assume that a fixed-parameter linear controller were designed to regulate the motor load in the region where the process gain (slope of the operating curve) is positive. Then, operation in the pad collapse region where the gain is negative clearly poses a problem because an unstable positive feedback loop would result. This situation would certainly occur if the setpoint were made greater than the maximum load, in which case the controller would actually accelerate pad collapse and induce plate clashing as a result of the gain change associated with traversing the critical gap. Clearly, if the nonlinear relationship shown in Figure 1.4 were deterministic, then plate clashes could be avoided by simply placing an upper limit on the plate position and/or the motor load setpoint. However, at each instant in time, the total axial thrust on the refiner plates is balanced by the reactive thrust of the pulp and steam between the plates. This reactive thrust is a complex function of refining conditions, and of the characteristics of the wood chips being fed to the refiner. Since thefibresare between the refiner bars for less than one second (Miles and May, 1990), conditions between the  Chapter J. Introduction  15  plates may vary quite rapidly, causing rapid changes in the maximum load and the corresponding critical gap at which it occurs. The relationship is also affected by plate wear, but this occurs at a much slower rate since the typical life of a refiner plate is roughly 1500 hours. Figure 1.5 shows a motor load plate gap curve for the industrial refiner used in this study, constructed from data collected over a 15 minute period of operation as discussed later on in Chapter 5. Note that the curve sketched out by the points in Figure 1.5 has the same general shape as the idealized curve shown in Figure 1.4. However, the amount of scatter in the data underscores the fact that there is a fair amount of uncertainty about the maximum load and the critical gap, and that both do, in fact, change over a relatively short period of time.  1.2.2 Control of Mechanical Pulping Wood is a biological material and, as such, it has a certain inherent variability in composition and properties. Intensive research efforts to develop control systems for refiner mechanical pulping plants have been underway for nearly 20 years. Survey papers by Michaelson (1978), Battershill and Rogers (1980), and Dumont (1986) give an overview of this effort An intense period of activity took place during the late 1970's and early 1980's, coinciding with the period of rapid changeover from grinding to refining. After a period of near inactivity in the mid 1980's, refiner control has again become an area of interest. One can attribute this to three major factors; 1. more and better on-line sensors for pulp properties and refining conditions 2. improved process understanding 3. increased computing power.  Chapter 1. Introduction  16  1.2.2.1 Refiner Control The systems of the late 1970's and early 1980's are best described as freeness control systems (Horner and Korhonen, 1980; Alsip, 1980; Rogers et al., 1980). Infrequent laboratory CSF measurements were used to periodically adjust the specific energy setpoint. This, in turn, changed the feed rate setpoint which was controlled by the metering screw speed, or the motor load setpoint which was controlled by manipulating the plate gap via the hydraulic pressure. For plate clash detection, plate gap sensors (Dahlqvist and Ferrari, 1981) and plate clash detectors using vibration monitors (Metsavirta and LeppSnen, 1979; Vespa et al., 1979; Alsip, 1980) were employed. As no sensor was available, refining consistency had to be inferred from material and energy balances. Much of this early work failed to be accepted by the industry. Laboratory measurements were too slow and imprecise and, despite early efforts to apply more advanced control techniques (Dumont et al., 1982), a reliable method of controlling the motor load could not be found. This is discussed in more detail in the next section. Recently, research and development in the area of refiner control has progressed along two separate lines. The first line of research is aimed at stabilizing refining conditions. Since a chip refiner is fed by volumetric means, both the bulk density of the chips in the metering screw and their moisture content can vary. The bulk density varies with the changing density of the wood and with the chip size distribution. The moisture content changes with wood species and the age of the chips. This alters the gravimetric flow of dry fibre to the refiner and the refining consistency. MacDonald and Guthrie (1987) proposed the use of a device to measure the total wet mass flow of chips directly inside the metering screw feeding the refiner. They concluded that, to individually control the mass throughput of dryfibreand water, the sensor would have to be combined with a chip moisture sensor. This was later confirmed by Foumier at al. (1991), who unsuccessfully  Chapter 1. Introduction  17  tried to measure chip moisture in the same locatioa Cort et al. (1992) argued that, provided the chip moisture content is relatively uniform, motor load variations are an indication of dry pulp mass throughput. Thus, they employed a control strategy which compensated for these throughput changes by using the metering screw speed to regulate the motor load and, therefore, the specific energy. This strategy reportedly obtained a reduction in freeness variability of over 50%. A similar strategy is currently under development by Cluett et al. (1992). However, it is unlikely that most mills have good control of chip moisture and, therefore, to be useful such a control strategy would have to be combined with a means of measuring and controlling refining consistency. It is interesting to note that Manner et al. (1986) used the converse argument, that motor load variations are primarily caused by chip moisture changes, to justify using the dilution water feed rate to control the motor load. Savonjousi et al. (1990) recently developed a consistency sensor for the blow line immediately down stream of the refiner. This allows direct control of refining consistency. Furthermore, when combined with the refiner motor load, the authors report that direct control of dryfibreand water flow rates through the refiner are achieved. However, even if the flow rates of fibre and water are well controlled, pulp properties can still vary because of changing wood species and/or chip properties. Different species of wood require different applications of refining energy because of differences in fibre morphology; fibre length, coarseness, etc. Chip properties can change with chip age, degree of rot, etc. The second line of research is aimed at addressing these problems by direct measurement and control of pulp properties. On-line freeness or drainage analyzers are now available with built-in latency removal and compensation for consistency and temperature (Brewster and Rogers, 1985; Tamminen et al., 1987; Kaunonen and Luukkonen, 1992). They sample the pulp stream, prepare the sample, and perform the  Chapter 1. Introduction  18  analysis in a total time of about 5 minutes. Toivonen and Tamminen (1990) recently used an on-line drainage analyzer to directly manipulate the metering screw speed to control freeness. However, in light of current thinking on characterization of mechanical pulps (see Section 1.2.1.2), this attempts to control what is actually a two parameter problem with only a single parameter. Work at the Swedish Forest Products Research Laboratory (STFI) has resulted in the pulp quality monitor or PQM (Hill et al., 1979; Nordin et al., 1980; Rydefalk et al., 1981) which measures pulp drainage, fibre length distribution (3 length classifications), and shive content. Recently, Strand et al. (1991) developed a control system based on the PQM to regulate average fibre length and freeness using plate gap and dilution water flow rate. However, their system runs in supervisory mode, i.e. the system makes recommendations and the operator decides whether to implement the required changes. While a significant step forward, this dependency on the operator does not permit the system to be commensurate with the extremely fast dynamics of the refining process.  1.2.2.2 Motor Load Control As mentioned in Section 1.2.1.5, the refiner motor load is a nonlinear function of the plate gap. Thus, a fixed-parameter linear motor load controller may actually accelerate pad collapse and induce plate clashing as a result of the gain change associated with traversing the critical gap. Although plate gap sensors and plate clash detectors using vibration monitors have proven very useful as last resort safety devices, they do not provide the information necessary to prevent the control system from inducing plate clashes. Productivity and/or pulp quality considerations may dictate operating the refiner close to the maximum load where pulp pad collapse is always a possibility. A first attempt at applying adaptive control to this problem was made by Dumont (1982). He applied an adaptive controller consisting of a recursive least squares parameter estimator employing  Chapter 1. Introduction  19  a variable forgetting factor (VFF) (Fortescue et al., 1981) and a Dahlin controller. He also used a set of rules to back out the refiner plates and restore the pulp pad when the gain changed sign. Choosing the variable forgetting factor design parameters to achieve the required trade-off between tracking slow process gain changes due to plate wear and fast changes associated with pulp pad collapse proved to be extremely difficult in practice. Although some success on an industrial chip refiner was reported, the strategy was deemed unreliable for continuous, unsupervised operation. There are several references in the literature to similar schemes (Jones and Rla, 1983; Koivo et al., 1985), although not much information is provided in either case. In another study, Dumont and AstrOm (1988) investigated several alternatives to improving the reliability of this adaptive controller, including an improved estimator due to Hggglund (1983) and a "dual" adaptive control approach. HSgglund's estimator, which is a method based on monitoring the parameter estimates and modifying the estimator gain when the parameters are thought to be changing, did not prove to be a substantial improvement over the VFF method. The development of an estimator capable of identifying both slow and fast changes was identified as an area requiring further study. However, the incorporation of caution and probing via dual control did prove to have substantial benefits for this problem. In dual control, the parameters of a linear process model are treated as state variables and therefore can change as rapidly as other system states. The dual controller will drive the output to the target value, but will also decrease the feedback gain (caution) when the uncertainties of the parameter estimates are large, and introduce perturbations (probing) when the uncertainties are large and the control error is small. This is made possible by including the parameter uncertainty (via the covariance or P-matrix) in the control law. The result is a regulator which can handle very rapid parameter changes, including frequent changes in the sign of the gain. Another potentially important contribution of this work was the proposed use of a nonlinear  Chapter 1. Introduction  20  performance index which explicitly penalizes operation in the pad collapse region. While the dual controller performed well in simulation, extensions to handle time delays and more complicated dynamics were identified as areas requiring more work before the method could be applied to an industrial refiner. As is evident from the preceding discussion, there are several limitations to basing an adaptive motor load control strategy on a linear model. First, a linear model is monotonic and, therefore, cannot characterize the motor load extremum. By extrapolation, the controller "thinks" that there is no limit to the maximum achievable load. Secondly, because of the input multiplicity, the controller must be modified (i.e. by using rules, a nonlinear performance index, etc.) to keep it from attempting to regulate the load in the pad collapse region. Recently, Dumont and Fu (1992) proposed yet another approach; this one based on a Wienertype nonlinear model obtained by expanding Volterra kernels using orthonormal Laguerre functions. An advantage of the nonlinear controller is that the motor load extremum and the input multiplicity are easily dealt with. That is, when the setpoint is less than the extremum, there are two possible inputs. The one corresponding to the larger plate gap is implemented. When the setpoint is greater than the extremum, there is only one possible input. This corresponds to the critical plate gap (see Figure 1.5). Disadvantages of the nonlinear model-based approach are that one must choose a suitable structure for the nonlinearity (this is likely not a trivial task) and that nonlinear controllers are generally more complex than comparable linear model-based controllers.  1.3 The Thesis The motivating force behind this thesis is that the motor load plate gap control loop is critically important to refiner control and, therefore, to the goal of producing uniformly high quality  Chapter 1. Introduction  21  mechanical pulp at the desired target properties. Research in pulp characterization and refining theory has demonstrated that the correct combination of two variables is required to determine refining action and, therefore, the pulp properties for a given raw material. There are only three final control elements on a chip refiner; chip feed rate, dilution water feed rate and plate gap, and, other than for the purpose of specific energy regulation which is required to stabilize the refiner at a particular operating point, the chip feed rate must be held constant to meet production requirements. Therefore, the two degrees of freedom needed for refiner control must come from the dilution water feed rate and the plate gap. The dilution water feed rate controls the refining consistency, and the plate gap, via its effect on the motor load, controls the applied specific energy. Refining consistency and specific energy manipulation are required to periodically change operating points in response to changing raw material characteristics. As yet, there does not exist a reliable method of closing the motor load plate gap control loop.  1.3.1 Objectives The objective of this thesis is to develop a reliable chip refiner motor load controller and to test it out on an industrial refiner. As originally proposed by Dumont (1982), and later refined by Dumont and AstrOm (1988), the idea is to treat the problem as a fault detection and control situation, and to use a linear model-based adaptive control approach consisting of an improved estimator and a controller containing dual features. Basing the controller on a linear model has several advantages; the linear model-based adaptive control field is much more mature in terms of theory and practical application than the nonlinear model-based adaptive control field, there is a rich collection of ideas and algorithms based on linear models already existing in the literature, and it permits the use of "black-box" models which describe the input-output behaviour of the system without having to go into  Chapter 1. Introduction  22  the detailed mechanism of how the refiner operates. The role of the improved estimator is to track both drifts and jumps in the parameters. The use of probing in the controller is justified by the fact that the parameter estimates are the key to identifying a pad collapse. Probing targets a portion of the input energy directly at identifying these parameters. The work was divided into three phases. The objective in phase 1 was to compare the output performance of a constrained certainty equivalence controller with a class of control algorithms known collectively as suboptimal dual controllers. The problem was a system with fast time-varying parameters, nonstationary disturbances, and frequent changes in the sign of the process gain. This is representative of the type of problem used to motivate the dual control approach. The outcome of phase 1 was a recommendation that the controller be based on a constrained certainty equivalence approach with a cost function extension to get active learning or probing. In phase 2, the aim was to extend the results of phase 1 to general dynamic-stochastic systems with deadtime. As a result, the active adaptive controller or AAC was developed. The goal in phase 3 was to combine the AAC with an improved estimator capable of dealing with sudden changes indicative of system faults and to test the resulting combination on an industrial refiner. To estimate the system parameters, a multi-model approach called adaptive forgetting through multiple models or AFMM was used. Another objective in phase 3 was to investigate the possibility of modifying the AFMM to include information about what to expect in the event of a pad collapse.  1.3.2 Contributions It is important to note that there have not been many attempts to compare adaptive controllers reported in the literature. Phase 1 of this thesis constitutes the first time that a constrained certainty equivalence controller has been fairly compared to a suboptimal dual controller on a problem of the  Chapter 1. Introduction  23  type used to motivate the dual control approach. Furthermore, the results show for the first time that there is no loss in output performance as a result of using the simpler constrained certainty equivalence approach. Finally, the benefits of actively probing the process are shown not to be as great as anticipated from results previously reported in the literature. However, it is shown that probing adds a degree of robustness against process-model mismatch in the parameter transition functions. Another contribution of this thesis is the AAC developed in phase 2. What differentiates the AAC from conventional adaptive controllers is that it contains dual features. The AAC exercises caution by simply penalizing large control actions. A simple, yet novel, method of estimating the effect of the current control action on the future parameter error covariance is used to build in probing. This was previously only possible for systems without deadtime. The AAC may have other applications than refiner motor load control. For example, it may find application in other fault detection and control problems. Furthermore, probing may be seen as an alternative to other methods of protecting adaptive controllers against problems caused by periods of low excitatioa Finally, probing is a potentially useful way of starting up an adaptive controller. Phase 3 constitutes the first known application of dual control principles in the process industries. Furthermore, this is the first known process industries application of the AFMM. Finally, the AAC-AFMM combination shows great potential for solving the chip refiner motor load control problem.  1.3.3 Overview The thesis is organized as follows. In Chapter 2, the linear CARIMA (controlled autoregressive integrated moving-average) model used to develop the controller is introduced and some  Chapter 1. Introduction  24  definitions are given to help clarify the relationships between some of the variables in the model. This is followed by the main body of the thesis consisting of phases 1, 2 and 3 mentioned above. These appear in Chapters 3,4 and 5, respectively. Finally, the overall summary and conclusions of the thesis are given in Chapter 6.  CHAPTER 2 DEFINITIONS  2.1 Introduction As mentioned in Chapter 1, one of the advantages of approaching the motor load control problem from the point of view of modelling the refiner as a linear dynamic system is that it permits the use of linear "black-box" models which describe the input-output behaviour of the system without going into the detailed mechanism of how the refiner operates. This is not to say that one can or should neglect the job of getting to know the process. On the contrary, many of the design choices which must be made in order to apply adaptive control systems based on black-box models (i.e. choices of presampling filter, sampling and control interval, model order, desired closed-loop response, etc.) are intimately related to prior knowledge of the process being controlled and the control objectives (Wittenmark and Astrom, 1984). That being said, Ljung and SOderstrOm (1983) give an excellent summary of the various linear models which may be used as well as methods to estimate their parameters. The estimation and control algorithms developed in this thesis are based on what has become popularly known as the CARIMA or controlled auto-regressive integrated moving-average model. The CARIMA model can characterize many types of input-output behaviour (time delays, nonminimum phase and open-loop unstable response, etc.) as well as the effect of realistic external disturbances and noise as they appear in the process output. It is used, for example, in generalized predictive control or GPC (Clarke et al.,  25  26  Chapter 2. Definitions  1987). In the following section, the CARIMA model is presented and a model to characterize how the CARIMA model parameters change with time is proposed. This model will be used to develop parameter estimators and controllers in subsequent chapters. Finally, some definitions are given to help clarify the relationships between some of the variables in the model.  2.2 Definitions Consider the controlled auto-regressive integrated moving-average (CARIMA) model A(z ~l)y(t) = B(z -1) u(t-k) + C(z -1) v(r)/V  (2- D  where y(t) is the process output or controlled variable; u(f) is the input or manipulated variable (y and u are deviation variables); v(t) is an independent normally distributed zero-mean random variable or white noise; A(z~l), B(z~l), and CCz'1) are polynomials in the backward shift operator z"1 (z'ly{t) = y(t-l)); V = 1 - z 1 is the differencing operator; and k is the deadtime. In general, the A{zl), 5(z_1), and C(z~l) polynomials are defined as A{z'1) = 1 + a.z"1 + ... + a z""° N  '  na  1  1  1  fi(z- ) = 60 + 61z-  +  ... + 6ntz-"i  (2.2)  C{z-l)= 1 + qz" 1 + ... + cmz-M If one defines the parameter vector as 8(0 - [«!»  ajt),  b0(t),...,bjt),  cx(t)  cJt)]T  (2.3)  where 0 has n = na+nb+nc+l elements, then the time dependency of the parameters may be characterized by the following Gauss-Markov type transition function 8(r+l) = 4> 0(0 + A + w(0  <2-4)  27  Chapter 2. Definitions  where <]> is the parameter transition matrix (2.5)  * = ^ [ i , K- K •»-• •«, •««] A is a gain vector used to adjust the expected mean value of the parameters A = a ,...,a , a,  a. , a  (2.6)  a  and w(f) is a vector of white noise variables w(0 = U ( 0  w (0, wb(t)  wb (0, wc(/)  wc (f)Y  (2.7)  Here, v(f) (Equation 2.1) and w<f) (Equation 2.4) are assumed independent. Formally, one can define the covariance matrices E{w(f)w\s)} - /?w(05o (2.8)  Zs{v(r)v(s)} = tfr(06„ E{wfl)v(s)} = 0 where 0 <. s <. t and b„ is the Kronecker delta  5 =  0  if t* s  A  if t = s  (2.9)  In Equation (2.8), Rw(i) is a matrix of parameter noise covariances where the parameters are assumed independent so the off-diagonal terms are zero R(t) = diag o2w(t)  a2w (0, o2 it)  o2 it), ol(i)  o2w it)  (2.10)  and Rv(i) = afa) is the variance of the output noise v(t). The system is completely characterized by  28  Chapter 2. Definitions  the input-output model (Equation 2.1) and the parameter transition function (Equation 2.4). In order to apply this model to any particular system the model orders {na, nb, nc), the deadtime (k), the parameter transition function (4> and A), and the covariances Rw(t) and Rv(f) must be specified. The decision to specify the covariances Rw(f) and Rv(i) as functions of time and the coefficients of the parameter transition function (4> and A) as constants is somewhat arbitrary. In fact, in the field of adaptive Kalman filtering they would all be considered unknown and possibly time-varying (see survey in Isaksson, 1988). While one might argue that such assumptions best reflect the complexity of the real world, the resulting parameter estimation algorithms are quite complex and may not be practical for application to a real system.  On the other hand, as mentioned in Chapter 1, a  characteristic of the refiner is that the static gain changes both slowly with plate wear and suddenly with pulp pad collapse. The parameter transition function has enough flexibility to model this situation by making the covariances Rw(t) and /?v(f) functions of time. An estimator capable of handling this situation which is not based on adaptive Kalman filtering, i.e. the AFMM, is discussed in detail in Chapter 5.  CHAPTER 3 COMPARISON OF ADAPTIVE CONTROLLERS  3.1 Introduction Adaptive control techniques provide a systematic approach to the design of high-performance control systems when the parameters of the process are unknown or changing with time. However, the problem of synthesizing adaptive controllers from a unified theoretical framework remains a major issue yet to be resolved. In a landmark series of papers, Feldbaum (1960) presented his "dual control" theory, which demonstrated that the optimal stochastic controller possess some interesting properties. Specifically, as mentioned in Chapter 1, the dual controller will drive the output to the target value, but will also decrease the feedback gain (caution) when the uncertainties of the parameter estimates are large, and introduce perturbations (probing) when the uncertainties are large and the control error is small. This is made possible by describing both the process and the parameters using a stochastic model. In essence, the parameter uncertainty (i.e. the covariance or P-matrix) appears explicitly in the control law.  The result is a regulator which can handle very rapid parameter changes.  Furthermore, it has been shown that there are at least two situations in which dual control may be especially warranted (AstrOm and Wittenmark, 1989); one is in the initialization phase of an adaptive controller when the uncertainty is large, and the other is when the parameters are varying rapidly and the static gain changes sign. Figure 3.1 shows the structure of the optimal regulator, which is composed of a parameter  29  30  Chapter 3. Comparison of Adaptive Controllers  estimator and a control law. The estimator generates both the parameter estimates (£) and their error covariances (P), i.e. the conditional probability distribution of the states given the measurements, sometimes called the hyperstate. The control law is a nonlinear function that maps the hyperstate into the space of manipulated variables. Although dual control provides a more rigorous theoretical framework, it is generally difficult to derive the controller. Conditions governing its existence are unknown. Under the condition that a solution exists, it is possible to derive a functional equation (called the Bellman equation) using dynamic programming. However, this can only be solved numerically, thereby making the optimal dual control law completely intractable, even for the simplest control problems (Stemby, 1976; AstrOm  Hyperstate -{x, P}  PARAMETER ESTIMATOR V }'  CONTROLLED PROCESS  CONTROL LAW u  Figure 3.1: Block diagram of the optimal stochastic regulator.  Chapter 3. Comparison of Adaptive Controllers  31  and Helmersson, 1982; Bernhardsson, 1989). This has led to the development of suboptimal approaches (for an excellent summary of these see Astrom and Wittenmark, 1989) that attempt to retain the dual property of the optimal controller. However, even with these suboptimal methods, there still does not exist a general design approach. In contrast, most "conventional" adaptive control schemes are based almost entirely on heuristic arguments. For example, they have been derived under the assumption that the parameters are constant but unknown. The unknown parameters of the process are estimated, and the estimates are used in place of the true values in the control law (i.e. the probability distribution in Figure 3.1 is replaced by a distribution with all mass at the conditional mean value Jf). The uncertainties are not taken into account in the design. This is called the certainty equivalence principle. Heuristics such as input weighting, estimator modifications and external perturbation signals are then added to obtain workable algorithms. Input weighting may be considered an analogue of caution since it is commonly used to achieve a degree of robustness against process-model mismatch, as well as to ensure that constraints on the variance of the manipulated variable are respected. It has been used extensively, for example in self-tuning control or STC (Clarke and Gawthrop, 1975) which, until recently, was one of the most widely applied adaptive controllers (Seborg et al., 1986). Estimator modifications such as deadzones, projection and leakage algorithms, variable forgetting factors, etc. are required to protect against problems associated with estimating parameters during periods of poor excitation (Wittenmark and Astrom, 1984; Astrom, 1987; Ydstie, 1987; Ortega and Tang, 1989). External perturbation signals, which may be added at the input or the setpoint, serve a similar purpose and may be considered an analogue of probing. However, it is important to note that there are very few, if any, guidelines on how and when to add external perturbations. Recent advances include the use of longrange prediction of the process output in determining the control actions (De Keyser et al., 1988;  Chapter 3. Comparison of Adaptive Controllers  32  Fisher, 1991), such as in generalized predictive control or GPC (Clarke et al., 1987). Some similarities obviously exist between the problem of chip refiner motor load control as described in Chapter 1 and the type of problem used to motivate the dual control approach, i.e. fast time-varying parameters, large uncertainties, and frequent changes in the sign of the static gain. Furthermore, as mentioned above, there is a duality which exists between the properties of caution and probing in dual control on the one hand, and input constraining and external perturbations in conventional adaptive control on the other. Therefore, it is of interest to compare conventional and dual adaptive control methods on an example exhibiting these features. In this chapter, the output performance of input constrained certainty equivalence and suboptimal dual adaptive controllers are compared on a process with fast time-varying parameters and under situations where the static gain frequently changes sign. Such a comparison has not appeared in the literature to date, although numerous comparisons of dual and unconstrained certainty equivalence controllers have been made previously; a survey by Wittenmark (1975b) covers much of this early work. In Section 3.2, the adaptive control problem is defined. The process under study is a first-order system with unit delay which is open-loop stable but subject to nonstationary stochastic disturbances. The parameters of the system are assumed to be time-varying according to GaussMarkov type transition functions. The various parameter estimators and control laws are presented in Section 3.3. Two passive (without probing) control laws are considered; a suboptimal LQ controller referred to here as a constrained certainty equivalence (CCE) controller and a suboptimal dual controller known in the literature as the innovations dual controller or IDC (Milito et al., 1982). Two parameter estimators are also considered. The CCE uses a recursive maximum likelihood (RML) estimator. The IDC, because of its structure, requires an extended Kalman filter (EKF). In addition, four methods of adding probing to make the passive controllers active are presented. In Section 3.4,  Chapter 3. Comparison of Adaptive Controllers  33  the various passive and active controllers are compared in Monte Carlo simulation. Sensitivity to plant-model mismatch in the parameter transition functions is investigated in both cases. Conclusions from the comparison are given in Section 3.5.  3.2 Problem Definition Jacobs and Saratchandran (1980) conjectured that simulation studies of a simple mathematical example could lead to conclusions having some general validity to higher-order systems, systems with deadtime and possibly even multivariable systems. They used a first-order system with time-varying parameters to compare the performance of several suboptimal adaptive control laws and their results are critiqued later on in this chapter. Here, a similar example is used to evaluate the various controllers mentioned above. Consider the following first-order process (na = 1, nb = 0, nc = 1) with unit delay (k = 1) (1 + a(t)z -J)v(0 = b(t)u(t-l) + (1 + c{t)z "Xfl/V  O- 1)  The parameters a(t), bit) and c(f) are time-varying according to the Gauss-Markov parameter transition function discussed in Chapter 2 and the subscripts on the parameters have been dropped for clarity. The parameter transition function (<I> and A) of Equation (2.4) and the covariance matrices (Rw and Rv) defined in Equation (2.8) are assumed known and time-invariant (these are discussed in more detail later on). One of the main differences between this example and the one used by Jacobs and Saratchandran (1980) is the nonstationary disturbance which is present due to the differencing operator V in the model. This has the effect of forcing one of the disturbance model poles to lie on the unit circle at z = 1. The net effect is that the output will tend to drift away from the target value unless  34  Chapter 3. Comparison of Adaptive Controllers  some control action is taken. This was not always the situation in the study by Jacobs and Saratchandran (1980) and this may have affected some of their conclusions. The inclusion of a nonstationary disturbance in Equation (3.1) gives the problem "process control" characteristics. Despite its simplicity, the system described by Equation (3.1) and the parameter transition function Equation (2.4) can exhibit a diverse range of behaviour. For example, the open-loop dynamics can be stable (|a| < 1) or unstable (|a| > 1). For -1 < a < 0, the open-loop response is an exponential rise to a new level which approaches a pure gain as a -» 0 and a pure integrator as a -» -1. Although the situation where 0 < a < 1 cannot arise from sampling a continuous first-order system, it is possible when the underlying continuous system is second-order and under-damped, and the sampling interval corresponds to half the period of oscillation. When the parameters are changing with time, frequent changes in the sign of b are particularly difficult to handle. This is similar to what happens with an operating refiner during pulp pad collapse, thus the interest in evaluating control strategies using this simple example.  3.3 Estimators and Controllers When the primary objective is to minimize output deviations from target, one might propose as a general cost function  /(A0=i^b2(^0M  (3-2)  where N is the cost horizon and the setpoint is assumed zero for simplicity. For large N, the minimization of Equation (3.2) with respect to Vw(f), Vw(?+1),..., Vu(t+N-1) leads to the optimal dual controller, which cannot be obtained in practice, even for simple systems. Therefore, only suboptimal  Chapter 3. Comparison of Adaptive Controllers  35  control laws and estimators which may be applied to real systems are considered. Furthermore, only "indirect" (controller parameters obtained indirectly via a design procedure) adaptive controllers are considered as they are more industrially relevant than "direct" (controller parameters updated directly) methods; indirect methods permit process information to be built in through the model and are the method of choice in most long-range predictive control (LRPC) schemes. Finally, use of the term "dual" is avoided when categorizing the controllers, except when referring to them by names used previously in the literature. This is done to avoid any confusion about what the controller actually does. Instead, the term "passive adaptive" is used to refer to controllers which are more or less cautious, but which do not actively probe the process, and the term "active adaptive" to refer to controllers which actively probe the system for information.  3.3.1 Passive Adaptive Controllers 3.3.1.1 Recursive Maximum Likelihood (RML) Estimator Several different methods could be used to estimate the parameters in Equation (3.1), depending whether c(t) is to be estimated or not. When c(t) is to be estimated, either recursive extended least squares (RELS) or recursive maximum likelihood (RML) estimators could be used. However, because of its superior convergence properties (Sflderstrom et al., 1978; Solo, 1979) the RML method is widely considered the most accurate recursive estimation method. The steps involved in implementing the RML for the first-order process under consideration are given in Ljung and SOderstrOm (1983). Define the parameter vector as 6(0 = [a(t) bit) dt)]T then the RML is implemented in the following steps;  V-V  36  Chapter 3. Comparison of Adaptive Controllers  Step 1: Measurement update Kit) = P(/|r-l)v r (0/[v(0/ , a|r-DV r (0+i? v ] e(0 = V;y(0-<p(06(f|f-l) 8(/|0 = d(f|f-l) + A:(Oe(0 ,  (3.4)  p(t\t) = /»a|r-i)-A:(i)v(o/ ak-i) r(0 = VXO-cp(OQ(r|r) Step 2: Parameter forecast (time update) 0(f+l|f) = 4>0(r|O + A P(H-1|0  (3.5)  = *P(r|0<I>r+^w  Step 3: Calculate control action Vu(t) (discussed in next section) Step 4: Update repressor and gradient vectors for next sampling instant <p(r+l) = [-VyO) VM(0 KO]  (3.6)  +  \j/(f l) = <p(t+l)-c(t+l\i)\\r(i) Step 5: Go to step 1 In the measurement update Equation (3.4), P is the parameter error covariance matrix, K is the Kalman gain, e is the prediction error and r is the residual error. In the parameter forecast (sometimes referred to as the time update) Equation (3.5), $ and A come from the parameter transition function defined in Chapter 2. In Equation 3.6, <p is a regressor vector and \j/ is a gradient vector (\j/ = -de/dQ).  3.3.1.2 RML-Based Constrained Certainty Equivalence (CCE) Controller The CCE or one-step optimal controller (Clarke and Hasting-James, 1971) is obtained by setting N = 1 in Equation (3.2), adding a constraint on the variance of the input, and employing  Chapter 3. Comparison of Adaptive Controllers  37  certainty equivalence to the minimization of the following cost function /(l) = E{[y\t+l) + pVu\t)] \t]  (3-?)  where p is an input weighting factor. The solution to Equation (3.7) involves the solution of a Diophantine equation, but is simplified by virtue of the fact that the process is invertible. It is, therefore, straightforward to show that the CCE is given by (see Appendix 1)  V  M  (0=-i[^±^l£^XO 6 2 + P (i + az-1)  (3.8)  where the time argument (r+1 \t) has been removed from the parameter estimates for clarity. Another way to derive the controller is to base it on the same model as the estimator. From the definition of tp(H-l) in Equation (3.6), it is immediately obvious that the one-step ahead forecast of the process output may be written as y(t+l \t) = $(f+1)0 + bVuif) + y(t)  <3-9)  where $0+1) = cp(H-l) | VH(/)=0. Substituting this expression for y(t+l) in Equation (3.7) and minimizing with respect to Vu(t) gives  VM(,)  = - *[«'+l)fl»X0] b +p  (3.10)  Here, the residual error r(t) (Equation 3.4) is used to reconstruct the innovation v(t), instead of using the current measurement and the previous control action as in Equation (3.8) (the equivalence between Equations 3.8 and 3.10 is shown in Appendix 2). The combination of the RML (Equations 3.4 to 3.6) and the CCE (Equation 3.10) shall be denoted RML-CCE.  Chapter 3. Comparison of Adaptive Controllers  38  3.3.1.3 Extended Kalman Filter (EKF) Two passive adaptive control laws which have come from suboptimal minimization of Equation (3.2), but without invoking certainty equivalence, namely the "cautious" controller and the innovations dual controller (IDC), are discussed below in Sections 3.3.1.4 and 3.3.1.5, respectively. These are based on a state space representation of the system and use an extended Kalman filter (EKF) to estimate the parameters. A rearrangement of Equation (3.1) yields the following one-step ahead forecast of the process output y(t+l) = -a(t+l)Vy(t) + b(t+l)Vu(t) + v(t+l) + c(t+lMi)+y(t)  C3-11)  To express the output in terms of past measurements of y and u and the a-, b- and c-parameters, the approach used by Jacobs and Saratchandran (1980) is followed. Define a new variable z(t) = y(t) v(0, i.e. a smoothed output After making the appropriate substitutions and simplifying, Equation (3.11) becomes z(r+l) - -a(t+l)Vy(t) + b(t+l)Vu(t) + c(t+l)[y(t)-z(t)] + y(t)  <3-12)  Equation (3.12) is used to derive the EKF and the control laws. The smoothed output (z) cannot be measured directly and, therefore, must be estimated along with the parameters. However, the one-step ahead forecast of z requires that one-step ahead forecasts of a, b and c have already been computed (see Equation 3.12). Therefore, the parameter and output forecasts must be made separately in two steps. Define the following state vectors  Chapter 3. Comparison of Adaptive Controllers  X (0  39  = [«(0 b(t) c(t) z(t-l)]T  (3.13)  T  x(t) = [a(i) b{t) c(t) z(t)]  Here, the state vector % is essentially an intermediate variable in the following two-part forecast, summarized below in Equations (3.14) and (3.15) X(f+D=¥*(0 + r+<B(f) *  o"  ¥ = 0 1  r=  A  (3.14)  m  co(0 =  0  0  where <t>, A and w<0 were defined previously, 0 is a vector with all elements equal to zero, and x(r + l)=Fx(f+l) + G / F=  0  -VX0 VM(0 [y(0-z(0] 0  0  G=  (3.15)  y(t)  where / is the identity matrix. The static linear relationship between the states and the measured output y is given by X'+l) = Hx(t+l) + v(t+l)  (3.16)  // = [0 1] The complete model (Equations 3.14 to 3.16) is composed of two coupled state transition functions, the first involving the parameters (a, b and c) and the second involving the smoothed output (z), and a static relationship between the states and the measured output (y). In the sequel, the two dynamic parts of the model are differentiated by referring to Equation (3.14) as the parameter transition function and Equation (3.15) as the output transition function. The state-space model presented above is nonlinear in the states (see Equation 3.12).  40  Chapter 3. Comparison of Adaptive Controllers  Therefore, an extended Kalmanfilter(EKF) is used to estimate the states from measurements of y and u. The idea behind the EKF (Jazwinski, 1970) is to extend linear Kalman filter theory to nonlinear systems by linearizing the system about the current state estimate, and then applying the Kalman filter to the resulting time-varying linear process. The EKF is summarized below Step 1: Measurement update Kit) = P(t\t-l)HTl\HP(t\t-l)HT  + R^  e(t) = y(f)-2(t\t-l) (3.17)  x(t\t) = x(t\t-l)+K(t)e(f) P(t\t) =  P(f\t-Y)-K(f)HP(f\t-\) r(f)=y(f)-i(f\t)  Step 2: Parameter forecast (time update)  n(M-i|0  (3.18)  'vypxt\ty9r+Rw  Step 3: Calculate control action Wu(t) (discussed in next section) Step 4: Output forecast  #r+l|fl-Ffc(r+l|0+G P(f+l\t) = LU(t+l\t)LT 0  F=  -vxo V«(0 KO o  G=  (3.19) /  L=  0  -VX0 VH(/) r(t) -Xs+\\i)  Step 5: Go to step 1 Equations (3.17) and (3.18) follow from the linear Kalman filter, with ¥ , r, and H as defined in Equations (3.14) and (3.16). The output forecast in Equation (3.19) makes use of the full nonlinear  41  Chapter 3. Comparison of Adaptive Controllers  model. The output transition matrices F and G, however, have been modified slightly from the definitions in Equation (3.15). This results from taking the expectation E{z(t+l)\t} in the output transition Equation (3.12), i.e. E{z(t+l)\t} = £(t+l\t)  (3.20)  = -a{t+\ \i)Vy(t) + b(t+l \t)Vu(t) + £(f-l \t)[y(t)-Z(t\t)} +y(t)-nc!(t+l \t)  The residual error r(t) = y(t) - z(t\t) (see Equation 3.17) is used in F to approximate the innovation reconstruction term y(t) - z(t) and the covariance TtC2(tt-l |0 has been added to G to account for the nonlinearity which this introduces. The error P(t+l \i) associated with this forecast is calculated by linearizing the output transition function about the current state estimates and then applying linear Kalman filter theory as suggested by Jazwinski. Thus, the bottom row in L is obtained by evaluating the following partial derivatives of the output transition Equation (3.12) at the current state estimates dz  -vxo,  dz lb  -  VM(O,  "5?  Kt),  dz  Tz  = -*#+! |0  (3.21)  3.3.1.4 Cautious Controller One approach to developing suboptimal dual control laws which has received a lot of attention in the literature is to minimize the one-step (N = 1) cost function (Equation 3.2). The main advantage of this approach is that an analytical expression for the control law is possible. Substituting for y(t+l) using Equations (3.15) and (3.16) and minimizing with respect to Vw(f) gives the control law  42  Chapter 3. Comparison of Adaptive Controllers  Efo[(l -a + c + az-l)y®-cz]]\tl  (3 22)  E{M\t] where the time arguments on the states, (t+l) for the a-, b- and c-parameters and (f) for the smoothed output z (see definition of the intermediate state vector % in Equation 3.13), have again been dropped for clarity. The conditional expectations in Equation (3.22) could be evaluated if the conditional probability distribution p(x(f+l) |0 were normal. The EKF (Equations 3.17 to 3.19) provides estimates of the mean %(t+l \t) and error covariance Il(r+1 \t) of p(x('+l)|0, based on the assumption that it is normal. Using standard results about moments of normal distributions (Bar and Dittrich, 1971), the control law is given by  =  _[{B-aB-nab  + Be + n^)^{ab  + nJz-']y{t)-{Mi 2  + b\a +  £nb^z^J  (3 23)  c  b +71,  where the time argument on the state and covariance is (t+l \t). From Equation (3.23) it is apparent that this controller uses all available information about the states, not just the estimates. The n-terms arise from the error involved in forecasting the states xO+1). and are a result of not employing the certainty equivalence principle. This has the effect of causing the controller to be cautious when the states are not well known, i.e. when %yj < n. The terms "cautious" or "neutral" are usually used to identify controllers of this form, i.e. resulting from one-step cost functions minimized without invoking certainty equivalence. It is easily shown that Equation (3.23) may be expressed in the following compact form  b +K4  Chapter 3. Comparison of Adaptive Controllers  43  where / = [0 1 0 0] and, as before, overstrike ~ denotes evaluation for V«(r) = 0. Although strategies such as this are appealing, they may exhibit undesirable properties such as turn-off (Wieslander and Wittenmark, 1971). This is when the overall gain of the controller decreases as the variance of b increases, i.e as P < nb. This will give less information about b at the next sampling instant, and the variance will increase further. The controller is then caught in a vicious circle, and the magnitude of the control actions become quite small. This tendency to be overly cautious has resulted in efforts to modify the cautious controller to maintain some of its desirable properties (i.e. robustness to modelling errors) while eliminating turn-off. This includes modifications to the objective function, an example of which is discussed below, as well as the addition of probing signals, which are discussed later on in this chapter.  3.3.1.5 Innovations Dual Controller (IDC) To avoid the turn-off phenomenon, Milito et al. (1982) suggested using the following extended cost function J{\) = E{[V \t+l) - KV 2(f+l)] 11]  (3-25)  which they call the innovations dual controller (IDC). They call the term 0 £ K <, 1 the learning weight. Their idea is to increase the variance of the prediction error at the expense of an increase in the variance of the output deviation from target. This is similar to an ideafirstproposed by Goodwin and Payne (1977). The IDC is given by  44  Chapter 3. Comparison of Adaptive Controllers  V««  H[b[FJL + G] + lLIir] b* + inb  0.26)  G = [0 (y(t)-xncz)]T where t = 1 - K. This controller can be considered a generalization of sorts. Consider the limiting cases. When i = 1, the IDC is equivalent to the cautious controller (Equation 3.24). When i = 0, the IDC is the certainty equivalence (CE) controller (see Equation 3.27 below). When 0 < i < 1, the controller strikes a compromise between these two extremes. In fact, contradictory to its name, the IDC does not achieve active learning in the way originally proposed by Feldbaum. Rather, it simply provides a means (through the adjustable parameter i) to adjust the degree of cautiousness and, therefore, it is included here along with the discussion of other passive controllers. The combination of the EKF (Equations 3.17 to 3.19) and the IDC (Equation 3.26) is denoted EKF-IDC.  3.3.1.6 EKF-Based Constrained Certainty Equivalence (CCE) Controller It is also possible to get a CCE controller using the state space structure. It is easily shown that the minimization of Equation (3.7) with certainty equivalence gives  b2 + p  (3-27)  G = [0 y(t)]T A clue to the fact that this controller may perform comparably to the IDC is the way the constraining factor p enters the control law. It is well known that b and its uncertainty (%) are most critical for control purposes. Comparing Equations (3.26) and (3.27), it may be seen that p and nb enter the control law in the same way. If the accuracy of b is of overriding importance, then one should be  Chapter 3. Comparison of Adaptive Controllers  45  able to adjust p so that CCE and IDC give similar results. This is illustrated later on in the results and discussion section. The combination of the EKF (Equations 3.17 to 3.19) and the CCE (Equation 3.27) is denoted EKF-CCE.  3-3.2 Active Adaptive Controllers There are many ways to design suboptimal controllers which actively probe the system to improve the parameter estimates, the idea being that the resulting reduced uncertainty about the parameters will lead to an improvement in the overall control performance. Astrtta and Wittenmark (1989) provide an excellent summary of the literature, and have arranged the methods into the following four general classes; 1. addition of a perturbation signal 2. constrained one-step minimization 3. extension to the one-step loss function 4. serial expansion of the loss function. The classes are listed in order of increasing computational complexity. The first simply involves the addition of an externally generated perturbation signal such as a square wave or a PRBS (pseudorandom binary sequence). In the second class, the controller is obtained by a constrained one-step minimization, either by limiting the minimum value of the control signal or the maximum value of some function of the covariance matrix. The approach in the third class is to extend the one-step loss function with a function of the covariance matrix. This leads to a loss function with two local minima, but it is possible to use an approximation to obtain an analytical expression for the input. In the fourth class, the methods involve an expansion of the loss function in the Bellman equation about the certainty equivalence or cautious controllers. The computations involved in this class are  46  Chapter 3. Comparison of Adaptive Controllers  quite complex and are not considered in this thesis. In this chapter, four methods of probing representative of the first three classes listed above are compared. Note that there are some subtle differences between what is outlined below and the references which is due to the controllers containing integral action as a direct result of the nonstationary disturbance in the system model Equation (3.1).  3.3.2.1 Square Wave Perturbation (PERT) Wieslander and Wittenmark (1971) suggested this be calculated according to Vw = Vuc + V5(-1)'  (3-28)  where 5 is the amplitude and the time argument on Vu(r) is omitted for convenience. V«e is the passive input, which in general may be CCE, IDC or any other non-probing control law.  3.3.2.2 Input-Constrained One-Step Minimization (ICM) Hughes and Jacobs (1974) suggested the input given by  Vuc Vu =  VUlsign(yu)  if |Vue| 2> V«,  (3.29)  if\Vuc\<Vut  where Vw, is a lower limit on the magnitude of input changes.  3.3.2.3 Hard Covariance-Constrained One-Step Minimization (HCCM) This type of constrained optimization was first proposed by Alster and Belanger (1974). They considered constraining the trace of the information matrix (P1) to be greater than some predefined value. However, since this involves constraining an aggregate of the parameter error variances, the  47  Chapter 3. Comparison of Adaptive Controllers  solution can become ill-conditioned when the trace of P becomes dominated by one parameter. This difficulty was obviated, but at the expense of having to choose three application dependent tuning parameters. Mosca et al. (1978) improved on this idea by constraining pb (i.e. the variance of b) only. Their method uses a time-varying threshold which depends on two user-specified tuning parameters. Here, an even simpler method is considered, i.e. one which requires only that pb <, pu and, therefore, involves only a single tuning parameter. Since this method explicitly uses information contained in the covariance matrix, the calculations are dependent on the type of estimator used. Here, the formulation for the RML is presented, for reasons that will become apparent later on. In this case, the input is calculated according to  Vuc V« =  if pb(t+2\t+l) £ pB  Vu0 + VuAsign<yuc - VMQ)  (3.30)  if pb(t+2 |f+1) > pu  where V«o and VwA are defined below. pb(t+2 \t+l) is used becausepb(t+l \t) is not affected by Vw(0 (the notation (t+2\t+l) is discussed in more detail in the next section). Substituting P(r\i) (Equation 3.4) into Equation (3.5) and advancing one sampling period in time gives  P(r+2|f+l) - 3> P{t+110-  P(t+l\t)V(t+l)y(t+l)P(t+l\t) ®T+R  (3.31)  Define the vector / = [0 1 0]. Then pb(t+2\t+l) = lP(t+2\t+l)lT and  (3.32)  48  Chapter 3. Comparison of Adaptive Controllers  \\r(t+l) = \j/(/+l)+ *V"(0  (3 33)  vj/(r+l) = #f+l) - £(t+l \i) v(0  <3-34)  "  where  and $(H-1) was defined in Equation (3.9). For convenience, the time arguments will be omitted for <p(t+l), Vw(f), 6(f+l|f) and P(t+l\f) in the sequel. Substituting Equations (3.31) and (3.33) into Equation (3.32) and simplifying (see Appendix 3) gives  plt+2\t+l) = 1 + aw * pbVu2 + 2dyu + d2 '  (3.35)  where the variables d^ are defined as  dl = IPtf d2 = d  ypyT+Rv =  s  d  (3.36)  d  Pb 2- i  d  4 = $U3  A plot of the function in Equation (3.35) is shown in Figure 3.2 as the curve marked "g" (Figure 3.2 is discussed in more detail in the next section). Note that pb decreases as \Vu \ increases. The "zerolearning" input V«o corresponds to the maximum of the curve and is found by solving dpjdu = 0, i.e. Wu0--dJPb  (3.37)  (first and second derivatives of pb are shown in Equations (3.48) and (3.49) in the next section). When the constraint in Equation (3.30) is violated (i.e. pb(t+2\t+l) > pu), the input required to meet the constraint can be calculated by replacing pb(t+2 \t+l) on the left hand side of Equation (3.35) with  49  Chapter 3. Comparison of Adaptive Controllers  pu and solving for Vu, i.e. V« = Vw0± V«A  (3.38)  where V«A is a deviation from V«o given by  (  V«A =  (3.39) 2  \  (/>-<)  and the sign of VHA corresponds to the location of VMQ relative to Vuc (see Equation 3.30). To summarize, when pb(t+2\t+\) evaluated at the passive input Vuc is less thanpu, the passive input is implemented. If not, the input is modified so that pb(t+ 2 \t+l) = pu. Vu,, is the input which will add absolutely no new information on b. V«A is the deviation from VWQ (in the direction of Vuc) which adds enough new information to meet the constraint.  3.3.2.4 Soft Covariance-Constrained One-Step Minimization (SCCM) This method is due to Wittenmark (1975a) and involves an extension of the loss function. While, in principle, any of the passive adaptive control loss functions considered in Section 3.3.1 could be used as a basis for the expansion, here the RML-CCE formulation (see Section 3.3.1.2) is used, again for reasons which will become apparent later on. Consider the following cost function  /(l) = E\[yXt+l) + pVu2(0 + W ) ] 11  (3.40)  where p is an input weighting factor, k is the probing weight and f{P) = pb(t+2) Equation (3.40), with f(P) as defined in Equation (3.41), is equivalent to  (3.41)  50  Chapter 3. Comparison of Adaptive Controllers  /(l) = E{[y2(t+l) + pVu\i)] \t} + \pfi+2\i)  (3-42)  Here, pb(t+2) can be removed from the expectation because it is a deterministic function of quantities known at time t (see previous section). However, it is not clear whether the notation should be pb(t+2 \t), as in Equation (3.42), or pb(t+2\t+l) as was used in the previous section. In fact, both are correct The former has the advantage of indicating when the calculations are performed, i.e. at the current sampling time t The later is consistent with the notation used in the estimator (Equations 3.4 to 3.6). pb(t+2\t+l) is the notation used throughout this chapter. / is minimized according to the procedure in Wittenmark and Elevitch (1985). pb(t+2\t+l) is expressed as function of Vw (see Equation 3.35). This results in an expression for / which is fourth-order in Vw and, therefore, has two local minima. / is expanded in a Taylor series up to second-order terms about a nominal input V«„ and minimized by Newton's method  Vu' =  Vu-.dJ,du2 VJIdu Vu  (3.43)  V«B is a "good" in initial estimate of Vw* and, as shown later on, is calculatedfrominformation related to the geometry of the problem. The solution is implemented in two parts. First, the cost function is minimized for X = 0. Roughly following the derivation in Section 3.3.1.2 and employing certainty equivalence to the first (quadratic) part of / in Equation (3.42) (denoted here by h) gives h = [$(f+1)9" + bVu(i)+y(t)f + pVw \t)  (3l44)  Here Equation (3.9) has been used to substitute for y(t+l \i) in Equation (3.42) (recall that $(r+l) = qp(r+l) |v«(o=o)- Minimizing h with respect to Vu(t) gives the passive input denoted by Vu£i)  51  Chapter 3. Comparison of Adaptive Controllers  Vu(t)  =  _ &[$(r+l)6+;K;)]  (3.45)  B +p This is the RML-CCE controller discussed in Section 3.3.1.2. In the second part of the solution, Vuc is used to guide the minimization of / . Denote the probing extension part of / by g (3.46)  g = pb(t+2\t+l)  An expression for pb(t+2\t+l) as a function of Vu was given in Equation (3.35) in the previous section. Ignoring terms independent of Vu gives d. 8 =  (3.47)  pbVu2 + 2d1Vu + d2  Now the derivatives in Equation (3.43) may be evaluated. Recalling that / = h + Xg, then  dJ  "5«"  = 2J 8[$(t+l)b + y(t)] + (b2 + p)V« - k  ^P" U + d^ pbVu + 2d 2  (3.48)  F" + tf  and  a2/ 2  du  - w <3 [ft v.^r-t] p Vu + 2dy»+4  = 2J (6 +p)+x_  (3.49)  2  b  To determine the nominal input VM„, refer to Figure 3.2 which shows the quadratic cost h, the probing cost g and the total cost / as functions of Vu. Identify the input minimizing h by the dashed line marked "c" (for the passive input Vuc). Vuc may not always be a satisfactory point to be  Chapter 3. Comparison of Adaptive Controllers  0.045  •0.2  0 0.2 Control Signal, Vu Figure 3.2: Minimization of cost function when |V«*| £ |VuB|.  08  0.09 0.08-  -0.3  -0.2  0.2 0.3 0 0.1 Control Signal, Vu Figure 3.3: Minimization of cost function when |Vu*| < {Vwfl|.  -0.1  0.4  Chapter 3. Comparison of Adaptive Controllers  53  expanded about because it can fall in the region where / has a local maximum. Identify the input maximizing g by the dashed line marked "0" (for the zero learning input VUQ). The global minimum is always on the same side of V«o as Vuc. Furthermore, it lies in the region where |Vw| £ |VMJ, where Vua is the inflection point of g, i.e the point at which d^g/du1 = 0, given by Vua = Vu0 + sign(Vuc-VU()fiJ3~/pb  <3-50)  and identified by the dashed line marked "a". Expansion about Vua will always give an approximate solution of the global minimum. However, an even better point may be found. First consider the situation when |Vu* | £ \Vua |, i.e. the situation shown in Figure 3.2. Identify the input Vub by the dashed line marked "b". It is reasonable to assume that J(Vub) » h(Vub) since g is small and can be neglected at Vub. Vu is between V«a and Vub, thus a suitable point to be expanded about is the nominal input Vu„ given by Vwn = (yua + VUb)/2  (3.51)  and identified by the dashed line marked "n". Vub is found by solving the second-order equation h(VUb) -J(Vua) = 0  (3.52)  When \Vu | < \Vua |, the method outlined above will still work because this most often occurs when g is small compared to h and / is nearly quadratic, as illustrated in Figure 3.3. Once the nominal input has been found, Equation (3.43) is used to solve for the optimal input W , identified in Figures 3.2 and 3.3 by the solid line marked "*". To summarize, the following steps are involved; 1. Calculate the passive input Vuc (Equation 3.45). 2. Determine vjf(r+l) (Equation 3.34). 3. Calculate d1A (Equation 3.36) and the "zero learning" input VHQ (Equation 3.37).  54  Chapter 3. Comparison of Adaptive Controllers  4. Determine the nominal input Vun (Equation 3.51) from Vwa (Equation 3.50) and Vub (Equation 3.52). 5. Use Equation (3.43) to calculate Vu, with the derivatives given by Equations (3.48) and (3.49). Note that both the HCCM and SCCM methods most closely approximate the dual effect, in the sense of Feldbaum, because they explicitly include the effect of Vu(t) on future uncertainties in the performance index.  3.4 Results and Discussion In this section, the efficacy of the previous control laws is evaluated using the time-varying process described by Equations (3.1) and (2.4). In all cases, <!>, A, R^, and Rv are assumed known. This is obviously not the situation in practice, but each controller is given the same prior information and, therefore, a comparison between the different approaches is valid. Furthermore, in a few cases the sensitivity to mismatch in <I> and Rw is investigated. Both the RML and the EKF have a number of parameters which must be specified a priori, namely; the parameter transition matrices <!> and A, the parameter covariance matrix Rw, and the output noise variance Rv. In addition, the parameter and parameter error covariance matrices must be initialized. This was done as follows. For the RML max  e, 6 ( 0 ) - 0.01  ec and for the EKF  Pa  PRML(0) =  0 0  0  0  max  0  Pb  0  max  Pc  (3.53)  55  Chapter 3. Comparison of Adaptive Controllers  6(0)  #0)-  EKF BKFffW =  P  o  PRML(0) 0  (0)  0  (3.54)  R  where 0 and Pmax are calculated from the coefficients contained in 3>, A, and Rw, and the following analytical expressions for the means and variances of the parameter transition functions  9  a.  '= T ^  Pr  =  (3.55) i-«  where i = a,b,c. Here, A = 0 and /""* = 0.1/. This value for the covariance corresponds to the required ± 3a to ensure that the parameters stay in the range -1 <. 6,- <. 1 and, therefore, that the A and C-polynomials are stable (\a\ < 1, \c\ < 1), i.e. their roots lie inside the unit circle in the z-plane. A and R„ were then used to adjust the "speed" of the transition functions, while maintaining constant variance. The output noise was set to RY = 0.1 in all cases. The first issue which needed to be resolved was how to conduct the simulations to get reproducible results. As a first attempt, 1000 runs of 1000 iterations each where conducted for a total of 1 million iterations using the cautious controller, Equation (3.24). The process described by <£ = 0.866/ and Rw = 0.025/ was used in this case (see Figure 3.9 for a typical realization but for a different control law). The mean output variance V(y) was calculated for each run and a histogram plotted. Figure 3.4 shows that the distribution is not normal; it is asymmetric with a long tail going out toward large values of V(y). Similar results were obtained when initial effects were disregarded, indicting that the occasional poor result was due to some rather complex behaviour under steady state conditions. Jacobs and Saratchandran (1980) observed similar behaviour in their simulations and conjectured that it might be due to the non-existence of a solution to the optimal controller (i.e. to Equation 3.2), or to occasional "escapes" (Hughes and Jacobs, 1974). Another explanation is that the  Chapter 3. Comparison of Adaptive Controllers  56  occasional poor result is due to the nonlinear nature of the estimation and control problem (due to correlated disturbances modelled by the C-polynomial), and is a result of approximations (linearization, etc.) which are necessary to obtain tractable estimation and control algorithms. At any rate, the nonnormality of the distribution in Figure 3.4 implies that the mean (of the individual run means) is not a representative statistic. Therefore, the median was used instead. Thus, each result shown in the sequel was obtained by conducting 1000 runs of 1000 iterations, calculating the median of the run output variances, repeating this three times, and averaging the run medians to get the final result. In effect, each result is the result of 3 million iterations of the simulated system.  200 180 160 140  ^120 "100 (D 3  £ 80 60 40 20  8  n  I I hlTL^ — 1.5  4.5 2.5 3 3.5 4 Output Variance, V(y) Figure 3.4: Histogram of output variances from 1000 runs of 1000 iterations each of the cautious controller.  2  57  Chapter 3. Comparison of Adaptive Controllers  3.4.1 Passive Adaptive Controllers Figure 3.5 compares the output performance of three estimator-controller combinations; (i) EKF-IDC, (ii) EKF-CCE, and (iii) RML-CCE. Note that the symbols in Figure 3.5 roughly correspond to the 95% confidence intervals. In general, however, error bars are shown explicitly in the results to follow. First, by comparing EKF-IDC to EKF-CCE, it may be seen that the CCE controller performs better than the cautious controller (i = 1) for a range of input weighting factors  10'  >  8"  iio° Q.  o  •-*  : EKF-IDC  ---*--e 10"  EKF-CCE j RML-CCE i  0.2 .014  0.4  i  .028  P  0.6 .042  Figure 3.5: Comparison of output variances for passive adaptive controllers.  0.8 .056  1 .070  Chapter 3. Comparison of Adaptive Controllers  58  0.007 < p < 0.07. This is because the CCE is less susceptible to turn-off. This result appears to conflict with results obtained by Jacobs and Saratchandran (1980), where a CCE controller was shown to be inferior to a cautious controller. However, their results were for the STC (Clarke and Gawthrop, 1975), which is the direct adaptive control counterpart to the indirect adaptive controller considered here. Furthermore, their results were correctly attributed to fixing the |30-parameter in the STC (which is commonly done to prevent parameter drift) even though the 6-parameter of the simulated process was time-varying and even changed sign. Continuing with the comparison of EKF-IDC and EKFCCE, it may be seen that CCE performs as well or better than IDC for this system. Furthermore, the minimum output variances were obtained for i = 0.5 and p = 0.021. The mean error covariance for the ^-parameter in the case of IDC was nb « 0.05 and, therefore, the product inb (» 0.025) was approximately equal to p in this case. Note that this is also the 6-parameter noise variance, i.e. alb = 0.025. This confirms the hypothesis stated in Section 3.3.1.6 that it should be possible to adjust p in CCE to equal the performance of IDC because of the way p enters the control law. This also implies that the covariances of the a- and c-parameters do not play an important role in the IDC. The above result is interesting because of the implications for extending dual control ideas to more general systems, particularly with deadtime where the output, the parameters and their error covariances must be forecast past time /+1. For the large general class of adaptive controllers derived by employing the certainty equivalence principle (such as CCE), there is no problem making the extension for cost horizons ofN> 1. This is because the solution essentially involves forecasting the output from a linear time series model, without concern for how uncertainty in the model parameters affects the forecast error. On the other hand, dual control methods (such as IDC) which attempt to use modelling errors (via the P-matrix) to constrain the controller are more difficult to extend for N > 1. The problem is due to the nonlinear relationship between the output forecast and the parameters  Chapter 3. Comparison of Adaptive Controllers  59  which becomes more complex as N increases (Wieslander and Wittenmark, 1971; Wittenmark, 1975a,b; Jacobs and Hughes, 1975). This difficulty may be obviated by using a direct adaptive control approach (Nahorski and Vidal, 1974; Chan and Zarrop, 1985), because the model may be reparameterized to express the predicted output as a linear function of past inputs and outputs. Drawbacks of this approach, however, are that it generally involves more parameters than the indirect approach, and it is more difficult to relate the parameter estimates to the physical system being controlled. A second observation from Figure 3.5 is that there is essentially no difference between EKFCCE and RML-CCE, i.e. the results are independent of the estimator. This also has implications for extension to systems with deadtime; it is more difficult to extend the EKF because additional delaystates are required (Jacobs and Hughes, 1975). The dimensionality of RML, on the other hand, is independent of the process delay, and only depends on the number of estimated parameters. One possible disadvantage of using CCE is that p is application dependent. The choice of t in EDC, on the other hand, does not appear to be to critical, provided it is not chosen too small. However, this is assuming the parameter transition functions (<I> and Rw) are accurately known. Figures 3.6 and 3.7 show what can happen if they are not. In Figure 3.6, the best (lowest output variance) controllers from Figure 3.5 (IDC, i = 0.5; CCE, p = 0.021) were applied under the conditions shown in Table 3.1. The "matched (f)" case denotes the situation of matched transition functions (<!> = <I>,fiw= Rw) and fast changing parameters (<I> = 0.866/, R„ = 0.025/). The results correspond to the minimums in Figure 3.5, i.e. the best CCE result is slightly better than the best IDC. In the "mismatched" case, the true parameter transition functions were altered to correspond to slow changing parameters (<E> = 0.995/, Rw = 0.001/) but the controller still assumed they were fast changing (i.e. the case of process-model mismatch in the  60  Chapter 3. Comparison of Adaptive Controllers  parameter transition functions). In this case, the output variances actually decreased, and there was little difference between IDC and CCE. In the "retuned" case, each controller was retuned (i in IDC and p in CCE) to minimize the output variance under the mismatched condition. This was achieved for i = 0.1 in IDC and p = 0.004 in CCE. Again, little difference between the two controllers was observed. Finally, the "matched (s)" case denotes the situation of matched transition functions and slow changing parameters and, therefore, a lower bound on output variance for the "slow" process. This was achieved for i = 0.35 in IDC and p = 0.001 in CCE. Conclusions which may be drawn from  0.7  0.6-  1. EKF-IDC 2. RML-CCE  0.5>  g"0.4 c .2 >  I  0.3  0.2  0.1  Matched (f) Mismatched Retuned Matched (s) Figure 3.6: Effect of parameter transition function mismatch on the passive adaptive controllers. Speed of true parameter transition functions overestimated.  Chapter 3. Comparison of Adaptive Controllers  Case  Process  Model  61  EKF-IDC l  RML-CCE P  <*>  Rw  4>  Rw  Matched (f)  0.866/  0.025/  0.866/  0.025/  0.50  0.021  Mismatched  0.995/  0.001/  0.866/  0.025/  0.50  0.021  Retuned  0.995/  0.001/  0.866/  0.025/  0.10  0.004  Matched (s)  0.995/  0.001/  0.995/  0.001/  0.35  0.001  Table 3.1: Conditions corresponding to results in Figure 3.6.  these results are that there is only a small cost associated with overestimating the speed of the parameter transition functions, and that both IDC and CCE are relatively robust to this type of process-model mismatch. In Figure 3.7, the converse situation was simulated; the controllers were tuned to the slow time-varying process (matched (s) case from Figure 3.6), and then the process was made fast timevarying according to conditions shown in Table 3.2. Error bars are not shown on the figure because they were essentially indistinguishable in this case. In the mismatched case, the output variances for both IDC and CCE were over 500,000% (note the log scale) greater than in the ideal (matched (f)) case. This large output variance was due to a lack of caution on the part of the controllers; in IDC the n-matrix (and therefore the product iIT) was undervalued due to the transition function mismatch, and p in CCE was simply too small. After retuning (i = 3.0 in IDC and p = 0.021 in CCE), the output variance was reduced substantially in both cases. However, IDC was still 330% greater than the ideal while CCE was only 120% greater, i.e. the output variance of CCE was about half that of IDC under the mismatch conditions. Furthermore, the value i = 3.0 required to achieve the best IDC result corresponds to K = -2, which brings the meaning of the IDC performance index Equation (3.25)  62  Chapter 3. Comparison of Adaptive Controllers  into question and makes IDC even more cautious than the cautious controller (the output variance for 1 = 1 was still 8700% greater than the ideal). Conclusions which may be drawn from these results are that it is more costly to underestimate than overestimate the speed of the parameter transition functions, and that CCE is more robust than IDC to mismatch between the assumed and true parameter transition functions.  10 4 F 1. EKF-IDC 2. RML-CCE  10°  >102 a>  g  •c cd  > =>  i  S-10 3  o 10u  -1  10  1 Matched (s)  Mismatched  Retuned  Matched (f)  Figure 3.7: Effect of parameter transition function mismatch on the passive adaptive controllers. Speed of true parameter transition functions underestimated.  Chapter 3. Comparison of Adaptive Controllers  Case  Process  Model  63  EKF-IDC l  RML-CCE P  3>  Rw  <f>  Rw  Matched (s)  0.995/  0.001/  0.995/  0.001/  0.35  0.001  Mismatched  0.866/  0.025/  0.995/  0.001/  0.35  0.001  Retuned  0.866/  0.025/  0.995/  0.001/  3.0  0.021  Matched (f)  0.866/  0.025/  0.866/  0.025/  0.50  0.021  Table 3.2: Conditions corresponding to results in Figure 3.7.  3.4.2 Active Adaptive Controllers The four methods of probing discussed in Section 3.3.2 were compared under the matched (f) conditions. The best RML-CCE controller (p = 0.021) was used to calculate the passive input signal. The results are shown with their 95% confidence intervals in Figure 3.8. Clearly, the SCCM method (A. = 2.0) gave the most improvement (about a 5.8% decrease) in output variance compared to no probing. Note that this improvement is over and above the best passive controller. Evaluations of active controllers are almost always based on a comparison with the cautious controller (IDC, i = 1.0), and when this comparison is made the best CCE-SCCM combination gave a 12% improvement. Wittenmark and Elevitch (1985) reported a 50% reduction in output variance when comparing their SCCM method to cautious control on an integrator with unknown gain. However, in this study essentially the same result was obtained on their example with a passive CCE controller and input weighting factor chosen equal to the gain parameter noise variance. Only a modest (further 6%) improvement was obtained by probing. This suggests that the benefits of probing may not be as great as previously thought, and that one must exercise a degree of caution when interpreting the literature in this area.  Chapter 3. Comparison of Adaptive Controllers  64  Returning to Figure 3.8, the square wave (PERT) did not improve output performance at all and similar results were observed when adding a PRBS at the input. Small improvements were observed for ICM (2.4%) and HCCM (1.3%). There are three reasons why the SCCM method gives the best performance; (i) the probing signal is continuous, (ii) it is calculated based on pb(t+2\t+l) which explicitly accounts for the effect of V«(0 on the future covariance, and (iii) it takes into account the direction of the passive control action, i.e. it tries to work with rather than against the controller. Figure 3.9 shows a realization of the CCE-SCCM combination. Note that the 6-parameter is tracked  0.75  Figure 3.8: Comparison of output variances for active adaptive controllers.  65  Chapter 3. Comparison of Adaptive Controllers  50  100  150  200  250 Time, t  300  350  400  450  500  50  100  150  200  250 Time, 1  300  350  400  450  500  50  100  150  200  250 300 350 400 450 500 Time, t Figure 3.9a: A realization of the CCE-SCCM combination. The dashed lines are the true parameters and the solid lines the estimates.  66  Chapter 3. Comparison of Adaptive Controllers  0.1 | a  0.05  ....  w  w  .... ,  r  V ^  ^  vj^  1  0  50  100  150  200  250 Time, t  300  350  400  450  500  50  100  150  200  250 Time, t  300  350  400  450  500  0.1 U  § 0.05-  I O  Figure 3.9b: Continuation of Figure 3.9a.  3  Probing Input, Vu p  Probing Input, Vu p  Probing Input, Vu p  Probing Input, Vu p  i  o  -*  en o  o o  01  o o o  fi  TJ  m  3  8 o o  o  Si o o  II  o b  68  Chapter 3. Comparison of Adaptive Controllers  more closely than the a-parameter is, and that the estimator does little to reduce uncertainty in the cparameter. This may be part of the reason why the best achievable output variance of V(y) = 0.58 is 480% greater than the theoretical minimum variance of V(y) = /?„ = 0.1, which could be achieved with perfect knowledge of the parameters and no input constraining. Figure 3.10 shows the four probing signals (Vup = Vw - Vwc) for a typical realizatioa The main difference is that the ICM and HCCM methods result in discontinuous probing signals; probing is used only when a constraint is violated. In contrast, the square wave and SCCM methods are continuous. Again the question arises as to the robustness of these probing methods with respect to  1.5  1. No Probing 2. ICM  I—i-  3. SCCM, >:  1  4. SCCM 2  CD  U  c .22 « > Q.  **  o 0.5  Matched (f) Mismatched Retimed Figure 3.11: Effect of parameter transition function mismatch on the active adaptive controllers.  Chapter 3. Comparison of Adaptive Controllers  Case  Process  Model  69  ICM  SCCMi  SCCMj  <D  Rw  4>  Rw  Matched (f)  0.866/  0.025/  0.866/  0.025/  0.80  2.0  1.0  Mismatched  0.866/  0.025/  0.995/  0.001/  0.80  2.0  1.0  Retuned  0.866/  0.025/  0.995/  0.001/  0.50  40.  1.5  Table 3.3: Conditions corresponding to results in Figure 3.11.  process-model mismatch. Figure 3.11 shows a set of results where the two best active controllers from Figure 3.8 (ICM with V«j = 0.8 and SCCM with A. = 2.0) were compared under the conditions shown in Table 3.3. A passive controller (no probing) and a modified SCCM controller (the original is denoted by subscript 1 and the modified by subscript 2 in the sequel) where also included in the comparison. In SCCMj, dA in Equation (3.36) is normalized by RJpb (Wittenmark and Elevitch, 1985). Thus, probing increases when the measurement noise increases. Also, because the covariance may vary over several orders of magnitude, probing is a function of the ratio of covariances pb(t+2 \t+\)lpb(t+\ \t), instead of pb(t+2 \t+l) alone. The overall effect is to make the choice of A, less critical (AstrOm and Wittenmark, 1989). In the matched (f) case, the performance of SCCMj with kx = 2.0 is equalled by SCCMj with 7^ = 1.0. Note that the difference between A^ and Aj is accounted for by the normalization factor, which is roughly 0.1/0.05 = 2.0 (see pb in Figure 3.9). In the mismatched case, the controller wrongly assumes slow parameter transition functions. The output variance increases for all controllers, but least so for SCCMj. Furthermore, when the various probing coefficients are retuned (see Table 3.3),"Kxmust be increased by a factor of 20, roughly the ratio of RJR\, to reach a minimum which is still slightly worse than the SCCM2 without retuning. The retuned SCCMj controller is 13% better than the retuned passive controller. One can conclude from  Chapter 3. Comparison of Adaptive Controllers  70  these results that probing is still beneficial under process-model mismatch conditions, and that the SCCMj modification adds a degree of robustness, particularly with respect to the choice of A. Figure 3.12 shows the input and output variances and the prediction error variances (output and three parameters) for the CCE-SCCMj combination, and for a range of p and A. in the matched (f) case (<!> = 0.866/, Rw = 0.0257). The objective here was to get a better feel for the interplay between input weighting and probing with respect to input-output performance and the quality of the parameter estimates. The minimum output variance occurs for p = 0.025 and A = 1.5, i.e. for values quite close to the ones found above by optimizing the performance with respect to p and A independently. The input manipulation increases with the amount of probing and decreases with the amount of input weighting. The prediction error variance follows the trend in input variance, and is related to the amount of information which is available to adapt. Estimation error variances associated with the a- and ^-parameters tend to decrease with increasing probing and increase with increasing input weighting, opposite to the trend in input variance, as one might expect. However, the estimate of the c-parameter gets better with less manipulation.  3.5 Summary and Conclusions The aim of this chapter was to compare the output performance of a constrained certainty equivalence controller with a class of control algorithms known collectively as suboptimal dual controllers. The problem was a system with fasttime-varyingparameters, nonstationary disturbances, and frequent changes in the sign of the process gain. This is typical of the type of problem used to motivate the dual control approach, which consists of two types of interactions; caution and probing. When considering "passive" control methods (i.e. without probing), the results indicate that the constrained certainty equivalence or CCE controller is as effective at minimizing output variance  Chapter 3. Comparison of Adaptive Controllers  71  2*0.06  0 0  0 0  S 0.095  Figure 3.12: Variances for the CCE-SCCM2 combination.  Chapter 3. Comparison of Adaptive Controllers  72  as a general suboptimal dual controller known as the innovations dual controller or IDC. In addition, the CCE controller is more robust to uncertainty in the parameter transition functions than the IDC. Four methods of adding probing to the passive controller to make it "active" were also evaluated. An approach which involves an extension to the one-step loss function with a function of the covariance matrix was the most effective at reducing uncertainty and improving overall performance. In terms of output variance reduction, the benefits of actively probing the process were not as great as anticipated from results previously reported in the literature. However, in terms of robustness to process-model mismatch, the results indicate that this type of probing is beneficial, and that a normalization modification to the method adds a degree of robustness, particularly with respect to the choice of the probing weight. Thus, the results presented in this chapter suggest that CCE is an effective approach to control problems where there is a large amount of uncertainty due to rapidly varying parameters and where the gain can change sign. The CCE controller exercises caution by simply penalizing control actions. Furthermore, the performance of CCE can be improved by the addition of a probing signal. This suggests an approach whereby the passive input is calculated based on constrained certainty equivalence, and probing is added by extending the loss function. In the next chapter, these ideas are extended to general dynamic-stochastic systems with deadtime.  CHAPTER 4 ACTIVE ADAPTIVE CONTROLLER  4.1 Introduction A result of the work described in Chapter 3 was a recommendation that a controller containing dual features be based on constrained certainty equivalence with a loss function extension to get probing. In this chapter, these ideas are extended to general dynamic-stochastic systems with deadtime. The main contribution of the chapter is a fairly simple and general method of adding probing to conventional adaptive controllers. To be consistent with the terminology of Chapter 3, the term "active adaptive controller" or AAC is used. The chapter is organized as follows. The problem is defined in Section 4.2. The general CARIMA model of arbitrary order is considered. The unknown parameters are assumed to vary as a random walk. In Section 4.3, the parameter estimator and controller are described. A RML (recursive maximum likelihood) estimator is used and is a straightforward extension of the RML for the first-order system of Chapter 3. The AAC consists of a constrained certainty equivalence approach coupled with a cost function extension to get active learning or probing, and is an extension of the RML-CCE-SCCMj combination discussed in Chapter 3. An extended output horizon is added to deal with nonminimum phase systems. A simulation example is used in Section 4.4 to illustrate some of the salient features of the algorithm. Conclusions are given in Section 4.5.  73  74  Chapter 4. Active Adaptive Controller  4.2 Problem Definition The general CARIMA system to be controlled was defined in Chapter 2. The parameter transition function (Equation 2.4) is taken to be 8(r+l)-8(r) + w(0  (4-J)  i.e. 4> = / and A = 0. Here, as in most applications of adaptive control, the assumption is that the parameters vary as a random walk. This is the most general assumption which can be made with respect to the unknown parameters (Ljung and Gunnarsson, 1990) and is implicit to many popular recursive estimators, including RLS (recursive least squares) with forgetting factor. Furthermore, for simulation purposes it is assumed that Rv is a known constant. Rw is adjusted in the parameter estimator (see Section 4.3.1) to achieve the usual trade-off between tracking speed and noise sensitivity.  4.3 Estimator and Controller 4.3.1 Recursive Maximum Likelihood (RML) Estimator As indicated in Chapter 3, several approaches can be taken to estimate the parameters in the general CARIMA model of Equation (2.1). However, because of some desirable convergence properties and because it is easily extended to higher-order systems with deadtime, the recursive maximum likelihood (RML) estimator is used. The RML for thefirst-ordersystem of Chapter 3 is given by Equations (3.4) to (3.6) in Section 3.3.1.1. It should be clear from those equations how the estimator for a model of arbitrary order is constructed. After combining the equations and dropping the conditional for clarity (i.e. the estimates are conditional upon information available at time t), the RML is summarized by  75  Chapter 4. Active Adaptive Controller  <P(0 «[-V;y(M)  -V^-na), Vu(t-k)  \y(t) =  Vu(t-k-nb), r(t-l)  r{t-nc))  <p(t)-ei(t)xV(t-l)-...-enc(t)w(t-nc)  e (0  = Vy(0-<P(09(0  (4 2)  '  6(r+l) = &b)+K(f)e(t) P(t+1) = P(t)-K(tMt)P(t) + Rw  K0 = W)-<pa)8('+l) where Q(f+1) and /^f+l) are the mean and covariance, respectively, of the conditional distribution p(0(r+l) \t). The other terms are as follows; <p is a regressor vector, \y is a gradient vector, K is the Kalman gain vector, e is the prediction error and r is the residual error. Rv is the output noise variance and Rw is the parameter noise covariance which must be specified to achieve the desired trade-off between tracking speed and insensitivity to measurement noise. Note that the diagonal terms in Rw (see Equation 2.10) may be used to independently adjust the rate of adaptation for each parameter.  4.3.2 Active Adaptive Controller (AAC) Consider the following cost function f RpM+k+l)) J = El [y(t+d) - yrf * pV« \i) +X " "° \  (4.3)  where E is the expectation operator, yr is the setpoint, d (d ^ k) is the output horizon, k is the deadtime, p (p ^ 0) is a control action weighting factor, and X (A, £ 0) is a probing weight. The first term in Equation (4.3) penalizes output deviations from target. By extending the output horizon  Chapter 4. Active Adaptive Controller  76  through the choice of d (Ydstie et al., 1985), the controller is able to stabilize nonminimum phase systems. The second term penalizes control moves. The input constraining factor p is a constant which makes the controller more or less cautious (see discussion of CCE in Sections 3.3.1.2, 3.3.1.6 and 3.4.1). The third term is a weight on the covariance of b0 used to introduce probing. This is a direct extension of the SCCMj method discussed in Chapter 3 (see Sections 3.3.2.4 and 3.4.2) for k > 1. Pbo(t+k+\) is used because it is directly affected by the current input Vu(/) (see Equation 4.2). RJPb0(t+k) is a normalization factor. The control objective is to minimize / with respect to V«(r). The output horizon d, the control action weighting factor p, and the probing weight X must be specified by the user. Wittenmark and Elevitch (1985) considered probing in the case of unit delay (k = 1 in Equation 4.3) only. As shown in Chapter 3, in that case an exact expression exists for Pba (t+2) as a function of Vu(i). When k > 1, no such expression exists because future covariances Pbo(t+i+\) for i ^ 2 depend on future unknown outputs and residuals. Here, a solution to this problem based on applying the certainty equivalence principle is proposed. In essence, the future unknowns in Equation (4.3) are replaced by estimates, obtained by recursively iterating through the RML equations replacing future unknowns by their minimum mean square error predictions at each step. The uncertainties in the estimates are not considered. / is then minimized according to Newton's method as originally proposed by Wittenmark and Elevitch (1985) and as discussed previously in Section 3.3.2.4. In order to avoid excessive cross referencing to Chapter 3, but with some duplication of effort, the method is shown in its entirety below. As before, the cost function is first minimized for X = 0. Again, employing certainty equivalence to the expectation in Equation (4.3) and denoting the quadratic part of / by h gives  77  Chapter 4. Active Adaptive Controller  h = [y\t+d) - yrJ + pV« \t)  <4-4)  Express the future output %t+d) as a sum of two components  where ${t+d) is the "free" response of the system (i.e. for Vw(0 = 0) and gdVu(i) accounts for the effect of the current input §d is the open-loop step response coefficient at lag d. Substituting this expression into the expression for h (Equation 4.4) gives h = [sf(t+d) + gdVu(t) - y j + p V« 2(0  (4-6)  Minimizing h with respect to VM(0 gives the passive input  (4.7)  Vu(t) §>p  The free response is obtained by recursively calculating the output increments Vj>/r+0 = cp(M-0©(M-l)  (4.8)  where / = l,2,...,d, using past estimates in the place of future unknown output increments and zero in the place of present and future unknown control actions and future residuals in the regressor veaor <p(H-0, i.e. Vy(t+i) = Vj>/f+z'), V«(r+/-1) = 0 and r(t+i) = 0 for / ^ 1 (see definition of <p in Equation 4.2). Note that the current parameter estimates Q(f+1) are used to predict the future output increments (also see thefilteringof future regressor vectors in Equations 4.17 and 4.19). This is consistent with the parameter transition function (Equation 4.1) which characterizes the parameters changes as a random walk for which the best future forecast is the current estimate. The free response is then  78  Chapter 4. Active Adaptive Controller  calculated as d  (4 9)  ?/*•**)-rtO + EVji/f+i)  -  i-i  §d is calculated recursively from coefficients of the A- and B-polynomials as follows $, = {b0 + V -1 +...) j J H k - (4, + d2z -1 +... )gl_1 where i = 1,2  (4-10)  d and s,.k is the unit step function defined as  s  i-t  =  f0  if i < k  1  if i ^ k  (4.11)  and g, = 0 for i £ 0. The second part follows Sections 3.3.2.3 and 3.3.2.4 closely. Vue(t) is used to guide the minimization of / . As before, denote the probing extension part of / by g (i.e. J -h + Xg) R ph(t+k+l) g =^ L 1 Pb(t+k)  (4.12)  where this time the normalization factor has been included in the definition of g. Note thatfromthe RML Equation (4.2), P(t+k+l) may be expressed as p(t+k+i) = p(t+k) - nt+kW(t+kMt+k)P(t+k)  + R  \j/(r+)fc)P(f+jfc)M'r(f+jfc)+/?v Define the vector / = [0  0, 1,0  0, 0  0]  (4.14)  i.e. a vector with n elements, all of which are zero except for a one in the location corresponding to b0 in 6. Then  79  Chapter 4. Active Adaptive Controller  p. (f+k+l) = lP(t+k+l)lT  (4.15)  y(t +k) = ty(t+k) + /Vu(f)  <4-16>  \p(t+k) = $(t+k) - ^(r+l)\|/(f+)t-l) -... - ^(r+l)y(t+k-nc)  ( 4 - 17 )  $(f+*) = <p(f+*)|Va(rH)  <4-18)  and  where  and  (M)-step ahead predictions of the covariance P(t+k) and the gradient \\r(t+k-l), required in Equations (4.13) and (4.17) respectively, are calculated recursively from the RML equations \\r(t+i) =  <p(t+i)-ei(t+l)\\r(t+i-l)-...-c'tK(t+l)y(t+i-nc)  K(t+i) = P(t+i)yrT(t+i)/[v(t+i)P(t+i)\\rr(t+i) + Rv]  (4 19)  -  P(t+i+l) = P(t+i)-K(t+i)\\r(t+i)P(t+i)+Rw where i = 1,2,...,M, again using past estimates in the place of future unknown output increments and zero in the place of future residuals in cp(?+0 (note, cp(t+k-l) involves only past inputs).  For  convenience, the time arguments will be omitted for ty(t+k), Vu(f), 6(f+l) and P(t+k) in the sequel. Substituting Equations (4.13) and (4.16) into Equation (4.15) and simplifying (refer to Section 3.3.2.3 and Appendix 3) gives  Pb (t*k+l) =  _ J +< <* phVu2 + 2d,Vu + d7 r  where the variables d ^ are defined as  b„  1  2  (4.20)  80  Chapter 4. Active Adaptive Controller  dx = ipy (4.21)  d2 = VJr/^ + fl, d =  3 P„A~dl  As before, Pbo(t+k+l) (Equation 4.20) is minimized for large values of VM, Also, Pbo(t+k+l) is maximized for VM„ = -dJPbt, i.e this value of VM adds absolutely no new information on b0. Substituting Equation (4.20) into Equation (4.12) and ignoring terms independent of VM gives  *  (4.22)  phVu2 + 2d,Vu + d.  where dA = R^dJPbo. Now the derivatives needed for minimization by Newton's method (see Equation 3.43) may be evaluated. Again, recalling that J = h + Xg, then  r i •> 8/ = 2) 6 [S d f<t+d) -y] * (tf * p) VM - X  d. p.Vu + d,  ±± !i pftVM + 2d1VM + d 2 ] 2  (4.23)  2  and ^ V M + rfl]2-d3] 82/ = 2J (A? + p) + A 2 3 du2 PbVu + 2di7u + d2]  (4.24)  The nominal input VM„ is determined exactly as shown in Section 3.3.2.4. To summarize, the following steps are involved; 1. Calculate the passive input VMC (Equations 4.7 to 4.11). 2. Determine (£-l)-step ahead predictions of the gradient \y(t+k-l) and the covariance  81  Chapter 4. Active Adaptive Controller  P(t+k) (Equation 4.19). Calculate \J/(r+£) (Equations 4.17 and 4.18). 3. Calculate dXA (Equation 4.21) and the "zero learning" input VMQ. 4. Determine the nominal input Vu„ (Equation 3.51) from Vua (Equation 3.50) and Vub (Equation 3.52). 5. Use Equation (3.43) to calculate Vu, with the derivatives given by Equations (4.23) and (4.24).  4.4 Results and Discussion In order to illustrate some of the salient features of the AAC, consider the following system y(t) = (b0 + bxz -'MtS) + v(r)  (4.25)  where the AAC is based on the model (1 + axz -l)y(t) - (b0 + bxz -1)a(r-3) + (1 + txz "XO/V Parameters b0 and bx are changed according to the schedule shown below in Table 4.1.  t  K  b,  0-150  -0.5  1.0  151-250  0.5  -1.0  251-350  -0.5  1.0  351-450  0.5  -1.0  451-500  -0.5  1.0  Table 4.1: Schedule of parameter changes for simulated system.  (4-26)  82  Chapter 4. Active Adaptive Controller  For these values, the system is nonminimum phase and the static gain changes sign every time b0 and bx are changed. Another feature of this system is that the disturbance is white noise. Therefore, unless a setpoint change is made, the best control in terms of minimizing output deviationsfromtarget is achieved for V«(r) = 0 (i.e. no control). This presents a problem for a conventional (passive) adaptive controller because, unless the setpoint changes coincide with the parameter jumps, there will not be enough information to adapt When a setpoint change actually is made, chances are that the controller will have the wrong sign and initially send the process output in the wrong direction. This is illustrated in Figure 4.1 for the AAC with probing disabled (K = 0). Here, the output horizon was made equal to the combined effects of the system deadtime and the inverse response (i.e. d = 4). The input constraining factor was somewhat arbitrarily set to p = 1. The parameter and regressor vectors were 6(0 = [0,(0. &„('). &,(0, c,(f)lr 1  i  J  ( 4 2 7 )  <p(0 = [-Vy(M). Vu(f-3), VtKt-4), r(t-l)] respectively. The parameter noise covariance matrix was Rw = diag[0.0OOOl, 0.01, 0.01, 0.00001 ]  (4-28)  and Rv = o*v =0.01. The estimator was initialized by setting 6(0) = 0, 60(0) = 0.5 and f(0) = /. Note that the setpoint changes are made in between the parameter jumps and, therefore, the manipulated variable initially goes in the wrong (positive) direction and the output starts to move away from the setpoint  Of course, as shown in Figure 4.3, the output eventually recovers as a result of the  information generated by the setpoint change, which allows the controller to adapt. Figure 4.2 shows the AAC with probing enabled (X = 10). Note that, each time the setpoint is changed, the manipulated variable goes in the correct (negative) direction and the output responds  83  Chapter 4. Active Adaptive Controller  ft/V^I—'iw  |0 ^^^•fiitt^wi^n^ O 0  1°  50  100  150  250 Time.t  300  350  400  r  •4h#-  -1 0  200  50  100  150  200  450  500  «»»n i. w y W "  250  300  350  400  450  500  250 Time, t  300  350  400  450  500  250  300  350  400  450  500  Time, t Figure 4.1: Control without probing (A, = 0).  | 0 •tyMSM^ftj*^O -2  0  50  100  150  200  -1 0  50  100  150  200  Time, t Figure 4.2: Control with probing (A. = 10).  84  Chapter 4. Active Adaptive Controller  <D  E is Q_  250 Time, t  500  to  u c  CD  o  o  50  100  150  200  250 Time, t  300  200  250 Time, t  300  250  300  350  400  450  500  a>  E  cd  iS  500  <D  o a) > o O  150  200  350  400  450  500  Time, t  Figure 4.3a: Parameter estimates for probing disabled case (Figure 4.1). The dashed lines are the true parameters and the solid lines the estimates.  85  Chapter 4. Active Adaptive Controller  0.5  I 01 GL  -0.5  8  £ §  50  100  150  200  10" Urni::  !!!!!!!!!!!!!!  10'nr\  Mi!!!!!!!!!!! !!!!!!!!! i^iiiiiiiiiiiiiii  250 Time, t  10'  0  10u  111 i 11  350  400  450  500  !!!!!!!!!!!!!!!!! !!!!!! Hi I l i i 5 1 i ! i i l l l i i S i I I i l ! I i i i i l l l i i i I l M III M Mill  :::::::::::::::;:::::::::;::::::::::::::::::::]  Him iii !!!!!!!!!!!!!!!!!!!!!!!!!?!!!!!!!!!!!!!!!!!!![!  ... r v i ^  10'2 n  300  TjTnTrnTrri  •M4444JJJ.  ::::;::::::::::•:::::::::;::::::::::;:::::::::]  liHirMM  Mim iii !!!!(!!!!!!!!!!ii!!i!!!!!N!!j!!!i!!i!!!ii!!!n  •  '  •  i  i  •  i  i  i  50  100  150  200  250 Time, t  300  350  400  450  500  50  100  150  200  250 Time, t  300  350  400  450  500  H  !!!!!!!!J!!!!!!!!!!!3!!!J!!J!!!!!!H  m ! MM IMMMMm!«!!!!! M n n ! ! ! ! m ! M H ! ! ! M ! ! m ! ! ! ! m ! m ! ! ! ! ! ! ! m ! ! ! ! ! n ! H H ! ! MM'  150 Figure 4.3b: Continuation of Figure 4.3a,  200  250 Time, t  300  350  400  450  500  Chapter 4. Active Adaptive Controller  86  l  \  E CD  is Q.  -2  j -  50  100  j -  150  200  250  300  350  400  450  500  Time, t  200  250 Time, t  300  500  200  250 Time, t  300  500  200  250 Time, t  300  CD  §  * 10 «  10  o 10 2 0  50  100  150  350  400  450  500  Figure 4.4a: Parameter estimates for probing enabled case (Figure 4.2). The dashed lines are the true parameters and the solid lines the estimates.  87  Chapter 4. Active Adaptive Controller  10°  ![!!!!!!!!!!!!!!!!!!! H!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !!!!!!  8  10"1h  i  io"2tn ! ! ! ! ! 1 1 i I.I ! ! 1TTT**tt^4iJjilM!!!!!Nn!!!!!ii!!^!!!M!i!NMM! !!!iii 10"  « \ l  •  1 1 lil 1 •  • • K I I I  Mji !:!!:! i l l ! ! ! : ! ! ! ! ! ! ! ! ! ! * ! ! ! ! ! ! ! ! i i j i i i i M i ! i U i H H M . . . . i .  1  0  I!!!!!!?!!!!!!!!!!!!!!!!!!  50  1  100  1  1  150  1  200  1  250 Time, t  1  1  1  1  1  —  400  Q.  1  450 500  •  | -O  is  _>  -1  :  : 1 ! 1 ! ! ! L! ! i 1 11 ! 1 ! ! j ! 1 ! ! 1 1 ! n  350  ••  ,^  i i « i i  f : : : : : : : : : : ; : : : : : :  1  300  "1  : : : : : : :  ,  _.Jt. 1  1  250 Time, t  300  350  400  450 500  300  350  400  450 500  1  -1.5  50  l!!  100  150  200  1  !!!!!!J!H!!!J!!!!!!!Jj!M!!!!!!!!!!!!H!!!:!^  h Ml!!!!!!:!!!!!!!!!!!!!!!!! !!!m!!!!!!!!i!!!!!!!M  0  50  100  150  Figure 4.4b: Continuation of Figure 4.4a.  200  250 Time, t  Chapter 4. Active Adaptive Controller  88  accordingly. The effect of probing is clearly seen in the manipulated variable (compare to Figure 4.1). As shown in Figure 4.4, the estimates track the true parameters much better than in the probing disabled case (compare to Figure 4.3) and have completely identified the jumps prior to the setpoint changes. Also, the covariance of both B0 and b1 have been substantially reduced, even though the covariance of b0 only is penalized in the cost function (Equation 4.3). A negative consequence of probing can be seen to be a slight increase in the output variance when the setpoint is constant. This, of course, is the information needed to adapt in the absence of a "persistently exciting" feedback signal. Lack of persistent excitation (PE) and what to do about it has been, and continues to be, an active area of research in adaptive control (Wittenmark and AstrOm, 1984; Seborg et al., 1986; Astrom, 1987; Ydstie, 1987; Ortega and Tang, 1989). While a formal quantitative definition of PE exists (see AstrOm and Wittenmark, 1989), here the term is used strictly in a qualitative way to refer to an input signal's suitability, or lack thereof, for parameter estimation. The problem stems from the fact that, as illustrated in the previous example, there is no guarantee the feedback signal will be PE. In fact, under some circumstances unique parameter estimates cannot be found. The lack of PE may cause the /"-matrix to become large resulting in parameter drift and, in severe cases, "bursting" in the parameters and the process output. In general, there is a feeling that external perturbations, which may be added to the setpoint or, as in the case of the AAC, to the control signal, will introduce unwanted disturbances and, therefore, are not acceptable. Thus, most research has been aimed at protecting the estimator from the effects of low excitation. Deadzones, projection and leakage algorithms, variable forgetting factors, etc. have all been used with some success. However, the problem with switching off the estimator is that this counters the symptom (parameter drift) and not the cause (lack of excitation). Furthermore, one could argue that this is not an effective solution in fault detection  Chapter 4. Active Adaptive Controller  89  situations, such as refiner motor load control, where the parameter estimates are used to diagnose a fault and, therefore, switching off the estimator is simply not possible. Obviously, the ideas contained in the AAC are not unique to a particular estimator or controller structure. For example, an ad hoc forgetting factor approach could be used in the estimator in place of Rv and Rw. Furthermore, virtually any "one-step" controller could be used to calculate the passive input V«c. Usually, a Diophantine equation is used to execute predictive control laws (see Section 3.3.1.2 for an example of this). However, one of the advantages of the implementation described in this chapter is that the same equations are used to generate the current prediction error and parameter error covariance in the estimator as to predict the future output and covariance in the controller, i.e. the recursions used to predict the output increments (Equation 4.8) and the future covariance (Equation 4.19) come directlyfromthe RML (Equation 4.2). This allows for an intuitively appealing implementation and helps to underline the close relationship between estimation and control in the AAC. A potential criticism of the AAC is its use of predicted output increments and residual errors in computing future parameter error covariances. No theoretical analysis has been carried out to examine the precision of the future covariance estimates or the propagation of error through the predictive equations. However, one should bear in mind that in all model-based controllers the control actions are based on the predicted response of the real system, without explicitly accounting for prediction errors (i.e. the certainty equivalence principle). Using the same predictions to calculate future parameter error covariances was simply the easiest and most obvious approach to developing the AAC. As shown above and in the next chapter, simulation and practical experience with the AAC indicate that it works and is simple enough to be applied to real systems.  Chapter 4. Active Adaptive Controller  90  4.5 Summary and Conclusions The purpose of this chapter was to build on the idea expressed in Chapter 3 that a controller containing dual features be based on constrained certainty equivalence with a cost function extension to get active learning or probing. As a result, the active adaptive controller or AAC was developed. The AAC represents a fairly simple and general method of adding probing to conventional adaptive controllers. It is an extension of the RML-CCE-SCCNlj combination discussed in Chapter 3 and is based on the popular CARIMA model structure. An extended output horizon was added to deal with nonminimum phase systems. A simulation example was used to illustrate a difficult yet common situation in adaptive control where the feedback signal does not provide enough "excitation" to adapt, and probing is required to continuously identify the process. While switching off the estimator is another way of dealing with periods of low excitation, this may not be possible in fault detection situations, such as refiner motor load control, where the parameter estimates are used to diagnose a fault.  CHAPTER 5 APPLICATION TO THE REFINER  5.1 Introduction As mentioned in Chapter 1, it was decided at the onset that the desired motor load control algorithm should consist of an improved estimator and a controller containing dual features. So far, the thesis has dealt with the controller. In this Chapter, the active adaptive controller or AAC, described in Chapter 4, is combined with an estimator capable of dealing with sudden changes indicative of system faults, and the resulting combination is tested on an industrial refiner. The chapter is organized as follows. In Section 5.2, the refining process and computer hardware and software setup are briefly described. The parameter estimator and controller are discussed in Section 5.3. The parameter estimator consists of a multi-model approach called adaptive forgetting through multiple models (AFMM) due to Andersson (1985). As discussed later on, it can track both drifts and jumps in the parameters. Also, the AFMM was modified for this application to include information about what to expect in the event of a pad collapse. The AAC, as described previously in Chapter 4, consists of a constrained certainty equivalence approach coupled with an extended output horizon to deal with nonminimum phase systems and a cost function extension to get active learning or probing. Two methods of dealing with the input multiplicity (see Chapter 1) are proposed. In Section 5.4, an industrial application of the AAC to a reject refiner is presented. Conclusions are given in Section 5.5.  91  92  Chapter 5. Application to the Refiner  5.2 Process Description Although the ultimate goal is to develop a motor load controller for chip refiners, due to production and product quality concerns, the initial trials reported in this thesis were carried out on a reject refiner. In general, however, reject refiners behave much like chip refiners in as much as both processes operate at high consistency (^ 20%) and identical machines are often used in both applications. The reject refiner used in this work is a Bauer model 489-5 double-disc refiner. The term "double-disc" refers to the fact that both refiner discs are driven by motors, as shown in Figure 1.2. The motors are rated at 3000 hp (roughly 2.25 MW) each and run at a fixed speed of 1200 rpm. Each disc is approximately 52 inches (130 cm) in diameter. The discs are positioned hydraulically  BEARING. HOUSING  OPENING PORT  CLOSING PORT HYDRAULIC CYLINDER  COMPENSATION LVDT  CONTROL LVDT  Figure 5.1: Schematic diagram of hydraulic plate positioning system for a Bauer refiner using LVDT position measurement.  Chapter 5. Application to the Refiner  93  and are controlled to a position setpoint. Plate position is indicated by two LVDTs (linear variable differential transformer) located on the movable shaft at the control end. However, because of plate wear and thermal expansion, this does not provide an absolute measurement of the gap. The LVDT signal is compared to a reference value and used to adjust the flow of oil into or out of the opening and closing ports of the hydraulic cylinder, Figure 5.1. This causes the position of the control end disc to change in relation to the position of the feed end disc (see Stebel and Aeby, 1980 for additional details). As discussed in more detail later on, the plate positioning system determines the refiner dynamics. Due to the lack of a high level programming environment in the distributed control system or DCS (Honeywell Inc., TDC 2000), the control algorithm discussed below was implemented in a personal computer or PC (Toshiba Model T5200, 20 MHz). The PC was connected to the DCS through a personal computer serial interface or PCSI (Honeywell Inc.) using a real time software package (The Fix™ v3.0, Intellution Inc.). The software provided the ability to incorporate user supplied C language programs, so the control algorithm was written in C. The sampling and control interval in the PC was determined by the "cycle time" of the DCS, which was one second. Again, as discussed in more detail later on, the lack of synchronization between the PC and the DCS data bases had the net effect of adding one sampling interval (1 s) of delay to the open-loop response of the process.  5.3 Estimator and Controller 5.3.1 Adaptive Forgetting through Multiple Models (AFMM) The literature on parameter estimation is extensive to say the least, and a comprehensive survey is well beyond the scope of this thesis. However, an excellent survey of recent results in this  94  Chapter 5. Application to the Refiner  field is given in Ljung and Gunnarsson (1990). This was used to narrow the search to a small number of diverse methods which, when taken as a whole, represent a wide range of available algorithms. When the parameters change quickly, the time variation of the parameter noise covariance RJf) can be used to adjust the estimator gain. Methods which estimate Rw{t) directly are called adaptive Kalman filters (see survey by Isaksson, 1988). A drawback of these methods is that they are extremely complex and difficult to implement They typically require running two estimators in a boot-strap approach; one to estimate the noise covariances and the other to estimate the parameters. Furthermore, extreme care must be taken to ensure that Rw is positive semi-definite. Therefore, methods based on adaptive Kalman filtering were not considered here.  yW.FV.Rv / 31(t).P1(l).a1(l)  MO.PsflWt)  -^81(t+1),P1(t+1),a1(t+1)  KALMAN FILTER  -810+1).Pl(t+1)+Rw,,q1 ^-§2(t+1 ),P2(t+1 ).a2(t+1)  POSTERIOR PROBABILITY p(y(t)|t) a,'(t)  Figure 5.2: Schematic diagram of AFMM with three models running.  Chapter 5. Application to the Refiner  95  Another way of dealing with abrupt system changes, that is not based on direct estimation of Rw(i), is to use the definition of w(t) as Gaussian random variables with covariance Rw(t) (see Equations 2.7 to 2.10). Here, Rw(t) is either R^ for drifting changes or R]w for jump changes, but it is not known when it switches between these two values. It is known, however, that for N data points, the true sequence RJf) is one of 2N possible combinations of R^ and R]w. It is obviously not possible to run all 2N combinations, so Andersson (1985) formulated a strategy that runs a fixed number of possibilities in parallel. The method is called AFMM (adaptive forgetting through multiple models). Figure 5.2 shows a schematic of the AFMM with three models running. Each model consists of a parameter vector (6,), a parameter error covariance matrix (P,) and a weight (a,). At the sampling instants, each model is assigned a new weight (aQ based on its posterior probability, and the parameter estimates and their error covariances are updated using a standard Kalman filter with Rw = Ri. The models are then ordered according to the weights. The least likely model (the one with the smallest weight) is terminated and the most likely model (the one with the largest weight) gives birth to a new one by allowing for a possible jump from it. This is done by increasing its covariance matrix by an amount Rl and giving it a small initial weight qx. The weights are then normalized and used to calculate a composite model. In this thesis, the AFMM is based on the RML (recursive maximum likelihood) estimator implemented in Kalman filter form. The RML was given previously in Chapter 4 (Equation 4.2). However, here /?v(0 replaces Rv (an estimator for fiv(t) is discussed in detail later on) and Ri replaces Rw, which is interpreted as the noise covariance for drifting parameters. As indicated above, a fixed number of models (A/) are created purely from the data (by running M estimators in parallel). At the sampling instants, each model is assigned a weight cx^(0 (i = \,...M) based on its posterior probability, calculated recursively using the standard Bayesian formula  96  Chapter 5. Application to the Refiner  afct) =  c,2(0  .exp  ,  oXO  (5.1)  2tp.(0/ ,.(0cpf(0 + ^ v ( 0  ^.(OP-COcpfCO+^CO  Here, the term in the square brackets is proportional to the probability density of e,<0 assuming a normal distribution, i.e. efj) e N(0,a«,) where ^ = cp,{0i°,(0<pl(0 + $/*) and a,(r) is a normalized weight for model i determined at the previous sampling interval (i.e. the sum of the cc^O's is one). Thus, when the square of the prediction error e](t) is small relative to the variance of the distribution 2  °"«,, e,(0 is n e a r  me  centre of the distribution and the probability of model i being the correct model 2  is high. On the other hand, when e\(t) is large relative to °"ej, then e,(i) is somewhere in the tail of the distribution meaning that the probability of model i is low. The importance of biasing the probabilities by the weights Oj(0 determined at the previous sampling interval will become evident later on. In essence, the family of models lives on with births and deaths, with the least likely model being killed and the most likely model giving birth to a new one by allowing for a possible jump from it at each sampling instant. This is done according to the following procedure proposed by Andersson (1985). First, the models are ordered according to the weights ij = arg  mm /,,. .  a,(0  (5.2) iM = arg  max /,,, . o (f) ;  so that i\ is the index of the worst model, i2 is the index of the second worst model, etc. Note that iM is the index of the best model. Next, the worst model is terminated and replaced with the parameter estimates and covariance of the best model  97  Chapter 5. Application to the Refiner  6.(f+l) = 6,.(r+l) (5.3)  This gives birth to a new model by increasing the covariance by RJW to allow for a jump. The sequences of past residuals and gradient vectors corresponding to the worst model are also replaced by those of the best  r,(.H) = r. (t-j)  (5.4)  ¥,('-/') = v,('-;) where j = 0,1  nc. Finally, the weights are normalized and an initial weight is given to the new  model M  a.(f+l) =  (!-<?,) aft)  (5.5)  a(r+l) = ql where qx (0 £ qx <, 1) is the prior probability of a jump. Once the worst model has been reinitialized, the following composite quantities may be calculated (see Blom and Bar-Shalom, 1988)  6(^1) = ^ a / f + l ^ + l ) i-i  P(f+1) = f)a(.(r+l)[/'.(f+l)  +  [6 ( .0+l)-9(^l)][6,a+l)-Sa+l)] r ]  i-1  (5.6)  r(0 = f>.(f+l)r.(r) i-i M  vK0 = E a .< r + 1 M« i-i  These quantities are used to execute the control law as discussed later on.  98  Chapter 5. Application to the Refiner  Ljung and Gunnarsson (1990) suggested that the AFMM may be extended (XAFMM) to include information about what changes to expect Forsman et al. (1991) used the example of a fighter jet trajectory problem where one may always want to keep a branch active that corresponds to a particularly dangerous action by the enemy, even though the posterior probability of that branch may be low. They went on to propose a general framework for combining "reasoning" and filtering using cluster analysis and the AFMM. The idea of keeping a branch active for a dangerous situation lead to the introduction of a new value of 0 into the family of models in the AFMM which is consistent with pad collapse. This was incorporated at the expense of killing the second least likely model 0'2) as follows B\(r+I) = e>+1) b.(t+\) = -b. (r+l)  where./ = 0,1  (5.7)  nb. Here, the parameters of the second worst model are replaced by the parameters  of the best model, but with the sign of the ^-parameters being reversed. The covariance of the second worst model is simply replaced by the covariance of the best model and, as in Equation (5.4), the sequences of past residuals and gradient vectors are also replaced. Note that when this branch is reinitialized, it is done after the worst model has been reinitialized and the composite estimates calculated. Furthermore, it may require a finite number of sampling intervals for the posterior probability of the pad collapse branch to build up, so it is not reinitialized at every sample (every fifth sample was used in this application). Therefore, when the second worst model is reinitialized, the weights must be re-normalized and initial weights reassigned according to  99  Chapter 5. Application to the Refiner  Y1  ( M  aff+1) = £a^(0-<(0-<(0 *—i V a.(r+l) = q1  a-q1-q2)ai(t) (5.8)  a.(r+l) = q2 where q2 (0 £ q2<, 1) is the prior probability of a pad collapse. Equations (4.2) and (5.1) to (5.8) comprise the method. The number of models M, the probability of a jump qlt the probability of collapse q2, and the expected covariance of drifts Ri and jumps RJW in the parameters are user specified design parameters for the AFMM. Initial values for a,<0), 6,(0) and Pt(P) are also required. Andersson (1985) showed that the AFMM is extremely sensitive to the accuracy of the noise variance Rv. Therefore, here Rv is estimated using the method outlined in his paper J$,(0 = t/?v(r-l) + (1 -x)[e\t) - <p(t)P(t)VT(t)]  (5-9)  where, as before, e(i) is the prediction error, calculated from the composite model. In general, either the prediction error e{t) or the residual error r(t) (see Equation 4.2) may be used here. In this work, e(f) was used because r(t) tended to underestimate the true value of RY in simulations. This occurred when the estimator gain was large causing the parameter estimates to over respond to the measurement noise. This reduced the magnitude of r(t), and thus /?„($, further increasing the estimator gain. The final result was that the estimator got caught in a vicious circle and eventually went unstable. Note that Rv(t) is calculated before the AFMM is executed. The filter time constant x (0 £ T <, 1) must be specified by the user as well as an initial value for £v(0).  Chapter 5. Application to the Refiner  100  5.3.2 Active Adaptive Controller (AAC) The AAC was discussed in detail in Chapter 4. For its application to the refiner, it was implemented exactly as described there except for the precautions required to deal with the input multiplicity (discussed next) and, as in the case of the AFMM above, $v(r) replaces Rv and R% replaces Rw. As mentioned in Chapter 1, when the adaptive control strategy is based on a linear model, the controller must be modified to keep it from attempting to regulate the motor load in the pad collapse region. Two methods of doing this were tried. In the first method, the sign of the step response weight gd (see Equation 4.10) was used to decide between controlling the motor load and retracting the plates, as shown in Figure 5.3. When §d < 0, the plates are retracted at a fixed speed of -0.1 %/s (this "retract rate" was a parameter used in the existing control system to automatically retract the refiner plates when triggered by an operator). This is essentially the strategy used by Dumont (1982). In the second method, shown in Figure 5.4, the feedback was simply made positive by changing the sign of the control error when §d < 0 and yr < $ft+d) (see Equation 4.7). This causes the plates to open at a rate proportional to the error. Both methods depend on an appropriate choice of the output horizon d. For the controller to be able to stabilize nonminimum phase systems, which are possible in this application because the 5-polynomial is first-order (see Equation 5.10 below), d must be chosen large enough so that §d has the correct sign. For a deadtime of k = 3 s and a worst case time constant of 1.7 s (Oj = -0.55, see Table 5.1), the process settling time is at most 10 s. Therefore, d = 10 is a safe choice. Note that for this value of d, gd is essentially equivalent to the open-loop gain of the process.  Chapter 5. Application to the Refiner  101  Motor Load (y)  Figure 5.3: First method of dealing with input multiplicity.  ft.P  DESIGN  -  ESTIMATOR  NEGATIVE FEEDBACK CONTROLLER  Setpoint (y,)  REFINER 4  POSITIVE FEEDBACK CONTROLLER  V  Plate Position^ 8d<0 &  (u)  yr<yt(t+d)  Figure 5.4: Second method of dealing with input multiplicity.  Motor Load (y)  102  Chapter 5. Application to the Refiner  5.4 Results and Discussion As a first step toward developing a refiner motor load controller for an operating refiner, a set of experiments were carried out to characterize the relationship between plate position and motor load. This was required to verify the choice of CARIMA model structure (Equation 2.1) and to investigate the effect of operating conditions and plate wear on model parameters. Each identification experiment involved manipulating the refiner plate position (LVDT) setpoint in the manner of a PRBS (pseudo random binary sequence). Plate position setpoint and motor load were logged once per second. Models were developed from the time series records using methods described in Box and Jenkins (1976) and MATLAB™'s identification toolbox. Table 5.1 shows operating conditions and modelling results for three representative experiments. A typical time series record corresponding to the 3rd experiment is shown in Figure 5.5. The identified models all had the same structure, given by (l+^z 1 )**) = ib0 + b1z-1Mt-3) + {l-c1z-1)v(f)/V  (5-10)  where, as before, y is the process output (motor load in MW), u is the process input (plate position in %), and v is an independent, normally distributed (white) noise. Note that this model characterizes the motor load as nonstationary (through the differencing operator V) even though the time series  Exp. No.  Plate Age (hours)  «i  b0  1  1460  -0.101  0.086  2  1460  -0.345  3  0007  -0.553  Ci  <*  0.109  -0.946  0.0113  0.179  -0.036  -0.984  0.0049  0.144  -0.031  -0.999  0.0047  *i  Table 5.1: Results of preliminary model identification experiments.  103  Chapter 5. Application to the Refiner  record shows that it is stationary. In effect, an integrator was forced into the model to get integral action into the controller. Experiments 1 and 2 were conducted at constant plate age (1460 hour old plates) just prior to a plate change. There are two main differences between the models derived from these experiments; (i) the 5-polynomial in model 1 indicates a fractional period of delay which was not present in model 2 and (ii) the residual variance in model 1 is more than double the residual variance in model 2. The varying delay may be due to a lack of synchronization between the PC and the DCS data bases. It is more difficult to assign a single cause to the difference in noise variances, although operators know that refiners are more inherently noisy at some operating points than others. Experiment 3 was conducted with new (7 hour old) plates. Model 3 is similar to model 2, indicating that the plate age does not affect the process dynamics as much as operating point. As a result of these identification experiments, the model structure was fixed atna= 1, nb = 1, nc = 1 and k = 3.  r s 4.5-  1 300 Time, t (s)  &  39  §38.5 £  600  38|-  I 37.5  400 500 600 300 Time, t (s) Figure 5.5: Time series record from open-loop identification experiment 3 conducted immediately after a plate change. 100  200  1  104  Chapter 5. Application to the Refiner  The parameter and regressor vectors then become 8(0 = K(r), b0«), bfl), L  Cl (0l J  r  (5.11)  qK0-[-Vy(M), Vu(r-3), Vu(r^), r(r-l)] respectively. Since the nonlinearity represented by the motor load plate gap curve is static, these dynamics should also apply in the pad collapse region. Figures 5.6 and 5.7 show part of a trial using the control strategy shown in Figure 5.3. Here, the design parameters for the AFMM were; M = 5 and q1 = q2 = 0.01. The parameter noise covariance matrices were Rj = diag[0.00O0l, 0.0001, 0.0001, 0.00001 ] (5.12) Ri = diag[0.0, 0.25, 0.25, 0.0] and the filter time constant in the output noise variance estimator (Equation 5.9) was T = 0.98. The estimator was initialized by setting o,(0) = 1/M, 0,(0) = 0, P,<°) = L and A>(0) = 0.005 (see o\ in Table 5.1). For the controller; d = 10, p = 0.1 and X = 0. This strategy (i.e. constant retract and no probing) roughly approximates the one used by Dumont (1982). For the period from roughly 950 s up to 1300 s in Figure 5.6, the plate gap was manually manipulated to bring the motor load up to about 3.6 MW. At this moment, the controller was switched on and setpoint changes were then made to 3.5 MW at roughly 1350 s and to 3.75 MW at roughly 1500 s. Note that, up to this point, the process gain was positive, i.e. increasing plate position (decreasing plate gap) increased motor load. At roughly 1700 s, the setpoint was changed to 3.9 MW. The controller could not achieve this new setpoint and the load began to decrease, indicating the onset of pad collapse (i.e. increasing plate position decreased motor load). The estimator quickly identified the gain change as shown in the plot of gd. The plates were then retracted enough for the pulp pad  105  Chapter 5. Application to the Refiner  &4-  u  minor pad collapse  1000  1100  1200  1300  1400  1500 1600 Time, t (s)  1700  1800  1900  2000  o  P °- 50  turrc-off  ^ 4 8 J=£ 1000  1100  1200  1300  1400  1  O)  1500 1600 Time, t (s)  1  1  1  1  1700  1800  1900  2000  ...  Eat  i  Q. 0)  CO  1  1000  1100  1200  1300  1400  1000  1100  1200  1300  1400  1500 1600 Time, t (s)  1700  1800  1900  2000  1700  1800  1900  2000  cc 0.02  1500  1600  Time, t (s) Figure 5.6a: Preliminary trial of motor load controller using control strategy shown in Figure 5.3.  Chapter 5. Application to the Refiner  -i  r  106  WU'WJui^WMJW  ^KIT  "~*Z^F  (0  <E -0.1 £  -0-2 J  I  J 1000  u 1100  J  L  1200  1300  1400  1500 1600 Time, t (s)  1700  1800  1900  2000  1000  1100  1200  1300  1400  1500 1600 Time, t (s)  1700  1800  1900  2000  1000  1100  1200  1300  1400  1500 1600 Time, t (s)  1700  1800  1900  2000  1000  1100  1200  1300  1400  1700  1800  1900  2000  -0.3  £ -0.4f I -0.5 ie -0.6 Q.  1500  1600  Time, t (s) Figure 5.6b: Continuation of Figure 5.6a.  Chapter 5. Application to the Refiner  107  to rebuild and the motor load to recover. Note that since the discrete time constant parameter a1 shown in Figure 5.6 was essentially zero over most of the run, it was possible to construct an approximate operating curve for this refiner (shown in Figure 1.5) by simply time shifting the inputoutput data to compensate for the delay. This helped to verify the nonlinearity. Another minor pad collapse occurred at roughly 1800 s. All things considered, this type of behaviour is acceptable. These minor pad collapses were not serious enough to clash the plates and the controller recovered very quickly (within 15 s). However, the more severe pad collapse which occurred at roughly 2300 s in Figure 5.7 was not acceptable. This occurred during a setpoint change to 4.0 MW. As before, the controller could not achieve this setpoint and the load began to decrease. In this case, the estimator was slow to identify the gain change (for reasons which are discussed below). This prompted the controller to decrease the plate gap even more, causing further decreases in motor load. Without the ability to adapt in this situation, a conventional fixed parameter controller would be caught in a positive feedback loop and would continue to decrease the plate gap in an attempt to increase the motor load, thereby accelerating pad collapse and eventually causing a plate clash. However, the adaptive controller eventually identified the gain change and began retracting the plates to restore the load. In the meantime, the setpoint was manually adjusted back to 3.75 MW so as to not induce another pad collapse. The main reason why the controller failed in this situation is because of a large oscillation at the Nyquist frequency (0.5 Hz) which appeared in the motor load signal at these operating conditions. This confounded the estimator because it uses motor load increments (see Equation 4.2), which are extremely sensitive to disturbances at this frequency. Although the motor load signal was being conditioned with an analog anti-aliasing filter, it had insufficient attenuation for a disturbance of this magnitude. This problem was subsequently corrected by redesigning the filter. Despite the fact that  Noise Variance, R„ o o o b o o P  !  <J> T -  00  M  ro  I.  ro  8  I a,  •i  ro  8  8  8  ro  ro  ro  ro  8  8  8  8  ro ro  ro ro  ro ro  ro ro  8  8  ro ro  ro ro  I s IS d 3 CD  •  ro  ro  • = *  I  Motor Load, y (MW)  * 8 « «  ro  1  p  o  *  Plate Position, u (%)  Step Weight, gd 9  U\  ro  ** u  S8 ro  ro ro  d 3 CD  -  ro  i  d 3  8  CD  =  f-  §8 ro  3^  -  I  8 ro CO  of  §8 ro  ro  o  o  o  K  ro  ro  ro  ro  ro in o  ro  ro  ro  ro  ro  8  8  o  o  ro oi o o  A  J  8  " * CO  to  o  m is  Q.  o  ro  8 o  i  £  109  Chapter 5. Application to the Refiner  1 1  Tp^^ry "S/— Vrt^" " ^"  O H " " P" 'Vr-"'Wii U.  2 -0.1 f "0-2  •Q  2100  2150  2200  2250  2300 Time, t (s)  2350  2400  2450  2500  2100  2150  2200  2250  2300 Time, t (s)  2350  2400  2450  2500  2100  2150  2200  2250  2300 Time, t (s)  2350  2400  2450  2500  2100  2150  2200  2250  2300 Time, t (s)  2350  2400  2450  2500  1-  OD  I-  Figure 5.7b: Continuation of Figure 5.7a.  Chapter 5. Application to the Refiner  110  the problem was not caused by the control algorithm, several lessons were learnt in this and other trials which prompted some design changes. For example, after examining the control increments during the large pad collapse at roughly 2300 s in Figure 5.7, it was observed that the controller made several large moves while the motor load was still dropping and it had not yet identified the gain change. Furthermore, once the gain change was identified, it took nearly 60 s for the controller to recover, i.e. to back out the plates. Two modifications were made to avoid large control actions. First, the control action weighting factor was increased by a factor of ten to p = 1. Second, limits were placed on the control increments, i.e. -VM""" £ Vu(f) £ ViT**, where Vif" = 0.2 %. To avoid the slow recovery, the second method of dealing with the input multiplicity was used (see Figure 5.4). This had the effect of making the retract control actions a function of the control error instead of being a constant, thereby allowing the controller to recover more quickly from a collapse. It was also observed that at times (for example during the period from roughly 1850 s to 1950 s in Figure 5.6) the control increments were quite small. This is the turn-off phenomenon mentioned in Chapter 3 (Wieslander and Wittenmark, 1971). Here the overall gain of the controller is small because §d is small relative to p. This results in small control actions and little additional information to improve the parameter estimates. The controller is then caught in a vicious circle. Since the success of the application is so critically dependent on timely identification of the gain, two changes were made to improve the quality of the estimates. First, probing was enabled by setting K = 10. Second, the parameter noise covariances were changed to R* = <W0.00001, 0.001, 0.001, 0.00001 ] (5.13) Ri = diag[0.0, 0.1, 0.1, 0.0J This put more emphasis on drifts and less on jumps in parameters b0 and bx. Thus, when the motor  Chapter 5. Application to the Refiner  111  load is moving slowly towards a pad collapse, the estimator is more likely to detect the gain change before a large load decrease occurs. This also allows the estimator to recover more quickly after a jump. To summarize, the main changes were to detune the controller, retract at a variable instead of a fixed speed, add probing, and increase the estimator gain. Figure 5.8 shows a trial with the modified controller. No pad collapses are encountered until the setpoint is raised to 5.2 MW at 700 s. Then four minor collapses are encountered, roughly 100 s apart. In each case, the controller recovers quickly from these (in the worst case the load recovers in 10 s). Note that the maximum load is roughly 1 MW higher than for the previous trial. This is mainly because the pulp mass throughput and consistency were higher in this case. Note also that the critical gap corresponds to a different plate position. This is mainly due to plate wear which causes the LVDT signal to drift. These trials were done roughly two weeks apart. However, the main difference between the two trials is in the control actions. Probing enhances the quality of the parameter estimates as seen in the gain estimate (§d), which is more consistent and recovers more quickly from the jumps. Also, there is no occurrence of turn-off. The tendency for the controller to cycle indefinitely when the motor load setpoint is higher than the maximum load is clearly undesirable. This occurs because, as discussed in Chapter 1, there are some drawbacks to basing the control strategy on a linear model. The response could be improved if the controller could identify the maximum load thefirsttime it encounters a pad collapse, and then use this information to set an upper limit on the setpoint, as originally suggested by Dumont (1982). However, since the maximum load is a function of operating conditions, provisions need to be made to periodically adjust the upper limit accordingly. One possibility is to reset the motor load setpoint to an estimate of the maximum motor load based upon "statistically significant" excursions of the gain estimate into the pad collapse (negative  112  Chapter 5. Application to the Refiner  w  |Pif V^ 600  700  800  900 1000 Time, t (s)  1100  1200  1300  1400  600  700  800  900 1000 Time, t (s)  1100  1200  1300  1400  600  700  800  900 1000 Time, t (s)  1100  1200  1300  1400  600  700  800  O)  E .2> 0 <D Q.  ,  2-2-  GC <D O  0.02  c to  0.015  >  0.01  0) (A O  0.005  1200 1400 1100 1300 900 1000 Time, t (s) Figure 5.8a: Trial of motor load controller using control strategy shown in Figure 5.4 and modified tuning parameters.  Chapter 5. Application to the Refiner  03  113  0-  (0 Q.  600  700  800  900 1000 Time, t (s)  1100  1200  1300  1400  600  700  800  900 1000 Time, t (s)  1100  1200  1300  1400  c OF a3 0_  -1  r^*?^^ps\^^k^^  foF a. -1  600  700  800  900 1000 Time, t (s)  1100  1200  1300  1400  600  700  800  900 1000 Time, t (s)  1100  1200  1300  1400  Figure 5.8b: Continuation of Figure 5.8a  114  Chapter 5. Application to the Refiner  gain) region. Define the open-loop gain as  i +< -I This notation is used to help make the point that $„ is the asymptotic value of the step response coefficient $, as i -» °°. As shown in Appendix 4, the variance of g„ is given by A  A.  (bn + byp  -2(b,+b,)(l+a,)(p  +p h) + (l+a.y(p, + 2phh+p.)  d+^)4 A statistical test to determine if §„ is less than some critical value gc is made by calculating the cumulative sum (CUSUM) denoted by s (Page, 1954) 5(0 = ^ - 1 ) + ^ - ^ - 5 / 2  (5-16)  where 5 is the size and sign of the shift to be detected. Usually, the CUSUM is intended for process monitoring and control situations, in which case $„ would represent the controlled variable and gc the target value. The shift parameter 5, then, would represent the maximum tolerable deviation from target. At the sampling instants, s would be compared to a critical boundary h to determine if a significant change had occurred. In this application, there is no "target" value; one is only interested in determining when the gain estimate has gone significantly negative. Furthermore, since a nominal value for g„ is not precisely known, it is not a simple task to determine suitable values for gc and 5. To get around this problem, one possibility is to set 5 = -2gc, in which case these two terms cancel one another out. Furthermore, it makes sense to normalize the CUSUM by the standard deviation of gn (this information is normally used in the calculation of h in the standard CUSUM procedure). After these changes, the CUSUM becomes  Chapter 5. Application to the Refiner  s(f)-s(t-l)*8jfrUJ  115  (5-17)  The CUSUM is calculated based on a user specified number of previous sampling intervals, the idea being to avoid setpoint adjustments based on old data or data from another operating point which may not represent the current state of the refiner. This also helps to reduce false alarms. A significant change is signalled when J < -h, where 3 <. h <. 5 (see Lucas, 1976). Once this happens, the motor load setpoint yr is made equal to the estimated maximum motor load /"*, provided that y™ is less than yT. Otherwise, the setpoint remains unchanged. The maximum motor load is estimated by filtering the measured motor load with a first-order Butterworth filter as follows yf(t) = O.82y'(M) + O.18()<0+y(M))/2  (5-18)  and then picking the maximum filtered value, again over a user specified number of previous sampling intervals. Figure 5.9 shows the results when the above procedure is used to adjust the setpoint The window length was set to 60 s and the critical boundary to h = 4. On two occasions, one during the setpoint change to 5.2 MW at 1600 s and the other during the setpoint change to 5.3 MW at 1800 s, the controller identified a gain change, the maximum load was identified, and the setpoint was modified. Note that the estimated maximum load is virtually identical in both cases. However, clearly the maximum load changes faster than every 60 s, as some minor pad collapses were encountered even after the setpoint had been adjusted. Some possibilities to improve on this would be to shorten the window length or to build in a margin of safety, i.e. modify the setpoint to, say, 95% of the estimated v*". A patent application has been filed on this (see Allison et al., 1993). Another potentially important application of probing which may benefit other adaptive control  116  Chapter 5. Application to the Refiner  1500  1600  1700  1800  1900  2000  2100  2200  2300  Time, t (s)  c-56ro  o  I 53  tc  1500  1600  1700  1800 1900 Time, t (s)  2000  2100  2200  2300  1500  1600  1700  1800 1900 Time, t (s)  2000  2100  2200  2300  0.02  •••yW"-i  § 0.015 >  )  0.01  § 0.005  2000 2200 2300 1800 1900 2100 Time, t (s) Figure 5.9a: Trial of motor load controller showing automatic setpoint adjustment or "reset" when the setpoint exceeds the maximum load.  1500  1600  1700  117  Chapter 5. Application to the Refiner  <D  0-  a. 1500  1600  1800 1900 Time, t (s)  1700  I  1  1  2000  , ._  2100  2200  2300  ...-!  I a. . 1  I  1500  1600  1700  1800 1900 Time, t (s)  1500  1600  1700  1800 1900 Time, t (s)  . 2000  2100  2200  2300  2000  2100  2200  2300  r  I -0. n,  is a. -1  f  I  1500  i  1600  1700  Figure 5.9b: Continuation of Figure 5.9a.  1800 1900 Time, t (s)  2000  2100  2200  2300  Chapter 5. Application to the Refiner  118  applications is during the start-up phase, which is shown in Figure 5.10. At roughly 130 s the controller is switched on and the setpoint is initialized with a bumpless transfer. Note that all the parameters start at zero. The first control action, therefore, is due solely to probing. The direction of this first move must be prespecified; the negative direction was chosen just in case the refiner was already close to a pad collapse. No other design modifications are required to initialize the controller, which has been routinely started up in this manner. Finally, some comments on a few of the "checks and balances" required during the application are in order. First, to protect the estimator, a check was made to ensure that |^ | <. 1. Otherwise, the gradient vector (see \j/(f) in Equation 4.2) may become unbounded. In the event that \tx\ > 1, the algorithm sets dr = \/dl (i.e. replaces C(z~l) by its spectral factor) in each of the M models in the AFMM Furthermore, since it is known that the true value of cx is negative (see Table 5.1), in the event that 6r > 0, the algorithm sets tx = -tx, again in each AFMM model. Constraints were also placed upon the output noise variance estimate, i.e. 0.002 <, $v(r) <. 0.02. Second, to protea the controller, the composite estimate of parameter ax is constrained to the range -1 <. ax <. 1, since it is known that the process is open-loop stable. Again, in the event that \ax \ > 1, the algorithm sets ax = l/dv  Note that by constraining the composite values only, the estimates in the AFMM are  unaffected. Simulation results indicated that limits on parameters b0 and Br (again, the composite values used in the controller only) could help prevent turn-off in the event of a large jump. Constraints, therefore, where placed on these parameters as well, i.e. -1 <. b0 and bx<.\.  5.5 Summary and Conclusions In this Chapter, the active adaptive controller or AAC, described previously in Chapter 4, was combined with an improved estimator and tested on an industrial refiner. To estimate the system  119  Chapter 5. Application to the Refiner  £6 2 >. 5  1 v^v^h* ***^*^* a 5 o 50  100  150  200  200  250 300 Time, t (s)  250  350  400  450  300  500  500  Time, t (s)  O)  Q. „  S -2 CO  50  250 300 Time, t (s)  350  400  450  500  250 300 Time, t (s) Figure 5.10a: Initialization phase of motor load controller.  350  400  450  500  C  100  150  200  0.02  8 S0.015h(0 •c  a >  0.01  $ Z 0.005  50  100  150  200  120  Chapter 5. Application to the Refiner  CO  5  J[itf?!Wrrvrtv^^^„  ca Q.  50  100  150  200  250  300  350  400  450  500  Time, t (s)  —r  ••"  iof  i  cB 0.  " •"  i  i  |  50  100  150  200  50  100  150  200  50  100  150  200  " 250 300 Time, t (s)  i  350  400  450  500  250 300 Time, t (s)  350  400  450  500  250 300 Time, t (s)  350  400  450  500  0)  ®  0  0.  Figure 5.10b: Continuation of Figure 5.10a.  121  Chapter 5. Application to the Refiner  parameters, a multi-model approach called adaptive forgetting through multiple models or AFMM was used. This provided the ability to track both drift and jump changes in the parameters. Also, the AFMM was modified to include information about what to expect in the event of a pad collapse. A set of preliminary open-loop time series identification experiments were used to chose a model structure and to investigate the effect of operating conditions and plate wear on model parameters. The results indicate that a fixed model structure is adequate and that the model parameters are more greatly affected by operating conditions than plate wear. Preliminary results of applying the AAC to the industrial refiner were generally successful. Two methods of dealing with the input multiplicity were compared. A method which involves making the feedback loop positive under a certain set of conditions proved to be effective and simple to implement The controller tended to perform better when the response to drifts (as apposed to jumps) were emphasised in the estimator and when the controller was configured to put more emphasis on probing and less on reducing output deviations from target.  This tendency toward favouring  estimation at the expense of control might have been anticipated given the multi-objective, faulttolerant nature of the control problem. The potential usefulness of monitoring the gain estimate to decide when the motor load setpoint is too high and should be lowered was demonstrated.  CHAPTER 6 GENERAL SUMMARY AND CONCLUSIONS  A characteristic of wood chip refiners is that the motor load is a nonlinear, time-varying function of the plate gap. The incremental gain between the motor load and the plate gap is subject to a slow drift due to plate wear and sudden changes in sign due to pulp pad collapse. A pad collapse can be caused by a change in operating point, or may occur suddenly due to a feed rate or consistency disturbance. This poses a problem for fixed-parameter linear controllers which may actually accelerate pad collapse and induce plate clashing as a result of the gain change associated with traversing the critical gap. Although plate gap sensors and plate clash detectors using vibration monitors have proven very useful as last resort safety devices, they do not provide the information necessary to prevent the control system from inducing plate clashes. Despite some effort to resolve this problem, there still does not exist a reliable method of controlling the motor load and, therefore, the applied specific energy. As a result, few industrial refiners are operated under closed-loop control, and wide swings in process conditions and product quality are permitted to occur. The objective of this thesis was to develop a reliable chip refiner motor load controller and to test it out on an industrial refiner. The idea was to treat the problem as a fault detection and control problem, and to use a linear model-based adaptive control approach consisting of an improved estimator and a controller containing dual features. While basing the control strategy on a nonlinear model has several advantages, mainly due to the fact that the motor load extremum and input  122  Chapter 6. General Summary and Conclusions  123  multiplicity are explicitly characterized and therefore easily dealt with, it was felt that the immaturity of the nonlinear adaptive control field was a distinct disadvantage. On the other hand, the relatively mature state of the linear adaptive control field meant that there was a larger theoretical and practical foundation upon which to base the control strategy. The idea of using a controller containing "dual" features is related to the fart that, in dual control, the parameters of a linear process model are treated as state variables and therefore can change as rapidly as other system states. The dual controller will drive the output to the target value, but will also decrease the feedback gain (caution) when the uncertainties of the parameter estimates are large, and introduce perturbations (probing) when the uncertainties are large and the control error is small. This is made possible by including the parameter uncertainty (via the covariance or Pmatrix) in the control law. The result is a regulator which can handle very rapid parameter changes, including frequent changes in the sign of the gain. The use of probing is further justified by the fact that the parameter estimates are the key to identifying a pad collapse, and that probing targets a portion of the input energy at continuously identifying these parameters. Finally, the role of the improved estimator is to improve tracking of both drifts and jumps in the parameters. Since there still does not exist a general dual controller design methodology, the main challenge was to extend existing suboptimal approaches to deal with deadtime and correlated disturbances, and to apply the resulting controller to an operating refiner. The work was divided into three phases. The objective in phase 1 was to compare the output performance of a constrained certainty equivalence controller with a class of control algorithms known collectively as suboptimal dual controllers. The problem was a system with fast time-varying parameters, nonstationary disturbances, and frequent changes in the sign of the process gain. This is representative of the type of problem used to motivate the dual control approach, and is similar to  Chapter 6. General Summary and Conclusions  124  what happens during a pad collapse. This constitutes the first time that a constrained certainty equivalence controller has been fairly compared to a suboptimal dual controller. The results show for the first time that there is no significant loss in output performance as a result of using the simpler constrained certainty equivalence approach. Furthermore, the benefits of actively probing the process are shown not to be as great as anticipated from results previously reported in the literature. However, probing does add a degree of robustness against process-model mismatch in the parameter transition functions.  The outcome of phase 1 was a recommendation that the controller be based on a  constrained certainty equivalence approach with a cost function extension to get active learning or probing. In phase 2, the aim was to extend the results of phase 1 to general dynamic-stochastic systems with deadtime. As a result, the active adaptive controller or AAC was developed. What differentiates the AAC from conventional adaptive controllers is that it contains dual features. The AAC exercises caution by simply penalizing large control actions. A simple, yet novel, method of estimating the effect of the current control action on the future parameter error covariance is used to build in probing. This was previously only possible for systems without deadtime. Finally, an extended output horizon was used to deal with nonminimum phase systems. The goal in phase 3 was to combine the AAC with an improved estimator capable of dealing with sudden changes indicative of system faults, and to test the resulting combination on an industrial refiner. To estimate the system parameters, a multi-model approach called adaptive forgetting through multiple models or AFMM was used. The AFMM was modified to include information about what to expect in the event of a pad collapse. The AAC-AFMM combination was then applied to a reject refiner. The controller was generally successful at regulating the motor load and dealing with pad collapse. Two methods of dealing with the input multiplicity were compared. A method which  Chapter 6. General Summary and Conclusions  125  involves making the feedback loop positive under a certain set of conditions proved to be effective and simple to implement The controller tended to perform better when the response to drifts (as opposed to jumps) were emphasized in the estimator and when the controller was configured to put more emphasis on probing and less on reducing output deviations from target. The potential usefulness of monitoring the gain estimate to decide when the motor load setpoint is too high and should be lowered was demonstrated. The AAC may have applications other than refiner motor load control. For example, it may find application in other fault detection and control problems. Furthermore, probing may be seen as an alternative to other methods of protecting adaptive controllers against problems caused by periods of low excitation. It is also a potentially useful way of starting up an adaptive controller. This thesis constitutes the first known application of dual control principles in the process industries, and it is the first known process industries application of the AFMM. Future work will focus on applying the AAC to a chip refiner and on incorporating it with existing refiner controls.  NOMENCLATURE  A(zl)  = polynomial in process model, Equation (2.1)  a-^a^a^,...  = coefficients of A(z~l)  B(z~l)  = polynomial in process model, Equation (2.1)  b0,bltb2,...  = coefficients of Biz'1)  C(z~l)  = polynomial in process model, Equation (2.1)  cltc2,c3,...  = coefficients of C(zl)  d  = output horizon in AAC, Equation (4.3)  rf!,^,^,,^  = intermediate quantities used in calculation of future value of pb, Equation (3.35)  e  = prediction error, Equation (3.4)  E  - applied specific energy, Equation (1.1)  E(zl)  - polynomial in Diophantine identity, Equation (A1.2)  f0  = first term in definition of gm Equation (A4.2)  /i  = second term in definition of gm Equation (A4.2)  F  = pulp mass throughput, Equation (1.1)  F  = transition matrix in EKF output transition function, Equation (3.15)  F(z~l)  - polynomial in Diophantine identity, Equation (A 1.2)  g  = probing part of /, Equation (3.46)  gc  = critical gain, Equation (5.16)  126  127  Nomenclature  gd  = open-loop step response coefficient at lag d, Equation (4.5)  g„  = open-loop gain, Equation (5.14)  G  = gain vector in EKF output transition function, Equation (3.15)  G(z_1)  = polynomial in identity, Equation (A1.5)  h  = quadratic part of/, Equation (3.44)  h  = critical boundary in CUSUM  H  = observation matrix in EKF, Equation (3.16)  ix,i2,i3,...  = index used to rank models in AFMM, Equation (5.2)  /  = average specific energy or intensity of impacts, Equation (1.2)  I  = identity matrix  /  = cost function, Equation (3.2)  k  = discrete deadtime, Equation (2.1)  K  = Kalman gain vector, Equation (3.4)  /  = vector with n elements all zero except for one in the location corresponding to b0, Equation (3.24)  L  = linearized transition matrix in EKF output transition function, Equation (3.19)  M  = number of models in AFMM  n  = the number of elements in 8  na  = order of A(zl), Equation (2.2)  nb  = order of B(zl), Equation (2.2)  nc  = order of C(zA), Equation (2.2)  N  = total number of refiner bar impacts on a unit of pulp, Equation (1.2)  N  = cost horizon, Equation (3.2)  = upper limit on covariance offe-parameterin HCCM, Equation (3.30) = net power applied, Equation (1.1) = parameter or state error covariance matrix, Equation (3.4) = elements of P = probability of a jump, Equation (5.5) = probability of a pad collapse, Equation (5.8) = residual error, Equation (3.4) = output noise variance, Equation (2.8) = diagonal matrix of parameter noise covariances, Equation (2.8) = diagonal matrix of parameter drift noise covariances in AFMM = diagonal matrix of parameter jump noise covariances in AFMM, Equation (5.3) = variable used in definition of Kronecker delta, Equation (2.8) = CUSUM, Equation (5.16) = unit step function, Equation (4.10) = discrete time, Equation (2.1) = process input or manipulated variable, Equation (2.1) = point at which g" = 0, Equation (3.50) = point at which KVub) - /(Vwa) = 0, Equation (3.51) = passive input used in active controllers, Equation (3.28) = lower limit on magnitude of input changes in ICM, Equation (3.29) = nominal (approximate) input, Equation (3.43) = probing input, VM - Vuc = zero learning input, Equation (3.30)  129  Nomenclature  VMA  = deviation input in HCCM, Equation (3.30)  Vu  = optimal input, Equation (3.43)  v  = white noise forcing output disturbance in process model, Equation (2.1)  iv  = column vector of white noises forcing parameter transition function, Equation (2.4)  wa,wb,we  = individual white noises in w  x  = state vector, Equation (3.13)  y  = process output or controlled variable, Equation (2.1)  yf  = free response of system, Equation (4.5)  yr  = setpoint, Equation (4.3)  z  = smoothed output in state space representation of process model, Equation (3.12)  0  = zero vector  Greek Symbols a,<t)  = weight in AFMM, Equation (5.1)  oc-(t)  = posterior probability in AFMM, Equation (5.1)  A  = gain vector in parameter transition function, Equation (2.4)  a a ,a f t ,a c  = elements of A  T  = gain vector in EKF parameter transition function, Equation (3.14)  Hz1)  = polynomial in identity, Equation (A 1.5)  5  = square wave amplitude in PERT, Equation (3.28)  5  = size and sign of shift to be detected in CUSUM, Equation (5.16)  ba  = Kronecker delta, Equation (2.8)  9  = parameter vector, Equation (2.3)  130  Nomenclature  t  = one minus learning weight in IDC, Equation (3.26)  K  = learning weight in IDC, Equation (3.25)  K  = probing weight, Equation (3.40)  n  = intermediate state error covariance matrix, Equation (3.18)  na,nb,nc  = elements of II  p  = input or manipulated variable weighting factor, Equation (3.7)  a  = standard deviation, Equation (2.10)  x  = time constant in output noise variance estimator, Equation (5.9)  3>  = transition matrix in parameter transition function, Equation (2.4)  $a>$b>tyc  = diagonal elements of 4>  X  = intermediate state vector, Equation (3.13)  vj/  = gradient vector in RML, Equation (3.4)  *P  = transition matrix in EKF parameter transition function, Equation (3.14)  co  = column vector of white noises forcing EKF parameter transition function, Equation (3.14)  tp  = regressor vector in RML, Equation (3.4)  Mathematical Operators C(-)  = covariance, Equation (A4.8)  E{ •}  = expectation operator, Equation (2.8)  9  = partial derivative operator, Equation (3.21)  V  = differencing operator, Equation (2.1)  V(-)  = variance, Equation (A4.8)  Nomenclature  = backward shift operator, Equation (2.1)  Superscripts /  = filtered value  max  = maximum value  T  = matrix transpose  -1  = matrix inverse  Overstrikes = estimate s vector with Vu(t) element equal to zero » mean value  Acronyms AAC  = active adaptive controller  CCE  = constrained certainty equivalence  CUSUM  = cumulative sum  EKF  = extended Kalman filter  HCCM  = hard covariance-constrained one-step minimization  ICM  = input-constrained one-step minimization  IDC  = innovations dual controller  PERT  = square wave perturbation  RML  = recursive maximum likelihood  132  Nomenclature  SCCM  = soft covariance-constrained one-step minimization  BIBLIOGRAPHY  Allison, B.J., J.E. Ciarniello, G.A. Dumont, and P.J.-C. Tessier, "Automatic Refiner Load Control", U.S. Patent Application No. 08/111,526, Filed Aug. 25, 1993. Alsip, W.P., "Digital Control of a TMP Operation", CPPA Process Control Conf., Halifax, 27-35 (June, 1980). Alster, J., and P.R. Belanger, "A Technique for Dual Adaptive Control", Automatica, 10, 627-634 (1974). Andersson, P., "Adaptive Forgetting in Recursive Identification through Multiple Models", Int. J. of Control, 42(5), 1175-1193 (1985). AstrOm, K.J., and A. Hehnersson, "Dual Control of a Low Order System", Department of Automatic Control, Lund Institute of Technology, Lund, Sweden, Report TFRT-7244, 1982. AstrOm, K.J., "Adaptive Feedback Control", Proc. IEEE, 25(2), 185-217 (1987). AstrOm, K.J., and B. Wittenmark, "Adaptive Control", Addison-Wesley, Reading, Massachusetts, 1989. Atack, D., "On the Characterization of Pressurized Refiner Mechanical Pulps", Svensk Papperstidn., 75(3), 89-94 (1972). Atack, D., and P.N. Wood, "On the High Consistency Operation of Large Double Rotating Disc Refiners", Svensk Papperstidn., 76(17), 639-644 (1973). Atack, D., "Towards a Theory of Refiner Mechanical Pulping", Appita, 34(3), 223-227 (1980). Atack, D., M.I. Stationwala, and A. Karnis, "What Happens in Refining", Pulp Paper Can., 85(12), T303-308 (1984). Bar, W., and F. Dittrich, "Useful Formula for Moment Computation of Normal Random Variables with Nonzero Means", IEEE Trans. Automatic Control, AC-16(3), 263-265 (1971). Battershill, J.W.T., and J.H. Rogers, "On-line Computer Installations: A Survey", Pulp Paper Can., 81(8), T177-182 (1980).  133  Bibliography  134  Bernhardsson, B., "Dual Control of a First-Order System with Two Possible Gains", Int. J. Adaptive Control and Signal Processing, 3(1), 15-22 (1989). Blom, H.A.P., and Y. Bar-shalom, "The Interacting Multiple Model Algorithm for Systems with Markovian Switching Coefficients", IEEE Trans. Automatic Control, AC-33(8), 780-783 (1988). Box, G.E.P., and G.M. Jenkins, "Time Series Analysis: Forecasting and Control", Holden-Day, San Francisco, 1976. Breck, D., W.D. May, and M.N. Tremblay, "Thermomechanical Pulping - A Preliminary Optimization", Int Mechanical Pulping Conf., San Francisco, 87-95 (June, 1975). Brewster, D.B., and J.H. Rogers, "Analysis of On-line Pulp Drainage Testers", Pulp Paper Can., 86(7), T186-190 (1985). Casey, J.P., "Pulp and Paper Chemistry and Chemical Technology", John Wiley & Sons, New York, 1980. Chan, S.S., and M.B. Zarrop, "A Suboptimal Dual Controller for Stochastic Systems with Unknown Parameters", Int. J. Control, 41(2), 507-524 (1985). Clarke, D.W., and R. Hasting-James, "Design of Digital Controllers for Randomly Disturbed Systems", Proc. IEE, 118, 1503-1506 (1971). Clarke, D.W., and P.J. Gawthrop, "Self-Tuning Controller", Proc. IEE, 122, 929-934 (1975). Clarke, D.W., C. Mohtadi, and P.S. Tuffs, "Generalized Predictive Control - Part I. The Basic Algorithm", Automatica, 23, 137-148, (1987). Cluett, W.R., J. Guan, and T.A. Duever, "Control and Optimization of TMP Refiners", Control Systems 92, Whistler, B.C., 159-162 (Sept., 1992). Cort, J.B., W.F. Keller, J.-P. Larose, and M. Goulet, "Automass™ Control Results in More Uniform Pulp Quality at Daishowa Inc., Quebec City", 78th CPPA Annual Meeting, Montreal, B165-171 (Jan., 1992). Croteau, A.P., G.C. Nobleza, and A.A. Roche, "Mill-Wide Analysis of the Sources of Paper Quality Variations using Time Series Analysis", Control Systems 90, Helsinki, 249-269 (Sept., 1990). Dahlqvist, G., and B. Ferrari, "Mill Operating Experience with a TMP Refiner Control System based on True Disc Clearance Measurement", Int. Mechanical Pulping Conf., Oslo, Session III: no. 6, 1-14 (June, 1981). De Keyser, R.M.C., PH.G.A. Van De Velde, and F.A.G. Dumortier, "A Comparative Study of Selfadaptive Long-range Predictive Control Methods", Automatica, 24(2), 149-163 (1988).  Bibliography  135  Duever, T.A., G.D. Mcflwraith, A.N. Hrymak, P.A. Taylor, P. Whiting, and E.W. Read, "The Application of Expert Systems to TMP Refiner Control", 77th CPPA Annual Meeting, Montreal, B259-267 (Jan., 1991). Dumont, G.A., "Self-tuning Control of a Chip Refiner Motor Load", Automatica, 18(3), 307-314 (1982). Dumont, G.A., N.D. Legault, J.S. Jack, and J.H. Rogers, "Computer Control of a TMP Plant", Pulp Paper Can., 83(8), T224-229 (1982). Dumont, G.A., "Application of Advanced Control Methods in the Pulp and Paper Industry - A Survey", Automatica, 22(2), 143-153 (1986). Dumont, G.A., and K.J. AstrOm, "Wood Chip Refiner Control", IEEE Control Systems Magazine, 8(4), 38-43 (1988). Dumont, G.A., and Y. Fu, "Nonlinear Adaptive Control via Laguerre Expansion of Volterra Kernels", 2nd Workshop on Adaptive Control: Applications to Nonlinear Systems and Robotics, Cancun, Mexico (Dec., 1992). Feldbaum, A.A., "Dual Control Theory. I-IV", Automation and Remote Control, 21, 874-880; 21, 1033-1039; 22, 1-12; 22, 109-121 (1960-61). Fisher, D.G., "Process Control: An Overview and Personal Perspective", Can. J. Chem. Eng., 69(1), 5-26 (1991). Forgacs, O.L., "The Characterization of Mechanical Pulps", Pulp Paper Magazine Can., 64(C), T89118 (1963). Foreman, K., L. Ljung, M. Millnert, and A. Skeppstedt, "Merging 'Reasoning' and Filtering in a Bayesian Framework - Some Sensitivity and Optimality Aspects", Int J. of Adaptive Control and Signal Processing, 5, 93-106 (1991). Fortescue, T.R., L.S. Kershenbaum, and B.E. Ydstie, "Implementation of Self-Tuning Regulators with Variable Forgetting Factors", Automatica, 17(6), 831-835 (1981). Fournier, M., H. Ma, P.M. ShaUhorn, and A.A. Roche, "Control of Chip Refiner Operation", Int. Mechanical Pulping Conf., Minneapolis, 91-100 (June, 1991). Franz6n, R., and R. Sweitzer, "Refining Forces in High Stock Concentration Pulping", Appita, 36(2), 116-121 (1982). Goodwin, G.C., and R.L. Payne, "Dynamic System Identification: Experiment Design and Data Analysis", Academic Press, New York, 1977. Guan, J., W.R. Cluett, and T.A. Duever, "Chip Refiner Modelling using the PLS Algorithm", 41st  Bibliography  136  Can. Chem. Eng. Conf., Vancouver, Paper 53-3, (Oct., 1991). Hagglund, T., "New Estimation Techniques for Adaptive Control", Ph.D. Dissertation, Department of Automatic Control, Lund Institute of Technology, Lund, Sweden, Report TFRT-1025, 1983. Hill, J., H. Westin, and R. BergstrOm, "Monitoring Pulp Quality for Process Control Purposes", Int. Mechanical Pulping Conf., Toronto, 111-125 (June, 1979). Horner, M.G., and S. Korhonen, "Computer Control of Thermo-Mechanical Pulping", 66th CPPA Annual Meeting, Montreal, A299-309 (Jan., 1980). Hughes, D.J., and O.L.R. Jacobs, "Turn-off, Escape, and Probing in Non-linear Stochastic Control", Preprints IFAC Symp. on Stochastic Control, Budapest, 343-352 (Sept., 1974). Huusari, E., and A. Syrjanen, "Today and Tomorrow - Jylha Tandem TMP", Pulp Paper Can., 81(3), T63-66 (1980). Isaksson, A., "On System Identification in One and Two Dimensions with Signal Processing Applications", Ph.D. Dissertation, Department of Electrical Engineering, Linkoping University, LinkOping, Sweden, Report 196, 1988. Jacobs, O.L.R., and D.J. Hughes, "Simultaneous Identification and Control Using a Neutral Control Law", 6th IFAC World Congress, Boston, Mass., Paper 30.2 (Aug., 1975). Jacobs, O.R.L., and P. Saratchandran, "Comparison of Adaptive Controllers", Automatica, 16(1), 89-97 (1980). Jazwinski, A.H., "Stochastic Processes and Filtering Theory", Academic Press, New York, 1970. Jones, R.E., and A. Pila, "Thermo-Mechanical Pulping Process Control", 69th CPPA Annual Meeting, Montreal, B105-111 (Jan., 1983). Kaunonen, A., and M. Luukkonen, "Practical Experiences using Continuous Freeness Measurement", Tappi J., 75(3), 159-164 (1992). Koivo, H.N., P. Sammalisto, M. Perkola, and J. Renfors, "A Microprocessor-Based Control System for a TMP Refiner", Mini and Microcomputers Conf., Sant Feliu, Spain, 93-96 (June, 1985). Laurila, P., G. Smook, K. Cutshall, and J. Mardon, "Shives - How They Affect Paper Machine Runnability", Pulp Paper Can., 79(9), T285-289 (1978). Leask, R.A., "Pulp and Paper Manufacture Vol. 2: Mechanical Pulping", Joint Textbook Committee of the Paper Industry, Atlanta and Montreal, 1987. Leask, R.A., "18th Annual Survey: Refiner Capacity Now Stands at 94,120 Tonnes a Day", Pulp Paper Can., 94(3), 14-16 (1993).  Bibliography  137  Ljung, L., and T. SOderstrOm, "Theory and Practice of Recursive Identification", MIT Press, Cambridge, Massachusetts, 1983. Ljung, L., and S. Gunnarsson, "Adaptation and Tracking in System Identification - A Survey", Automatica, 26(1), 7-21 (1990). Lucas, J.M., "The Design and Use of V-Mask Control Schemes", J. Quality Control, 8(1), 1-12 (1976). MacDonald, J.E., and J.J. Guthrie, "Chip Mass Flow Meter", Int. Mechanical Pulping Conf., Vancouver, 187-192 (June, 1987). Manner, H., M. Perkola, and E. Huusari, "Consistency Control - The Fundamental Factor in TMP Refining", Tappi Pulping Conf., Toronto, 223-227 (Oct., 1986). Metsavirta, A., and M. LeppSnen, "Improved Profitability of Thermomechanical Pulping with Computer Control and Effective Heat Recovery", Int. Mechanical Pulping Conf., Toronto, 373-381 (June, 1979). Michaelson, R.B., "A Summary of Thermomechanical Pulping Plant Advanced Control Applications", ISA PUPID/PMCD Symposium, Portland, 65-78 (Apr., 1978). Miles, K.B., H.R. Dana, and W.D. May, "The Flow of Steam in Chip Refiners", International Symposium on Fundamental Concepts of Refining, Appleton, Wisconsin, 30-42 (Sept., 1980). Miles, K.B., and W.D. May, "The Flow of Pulp in Chip Refiners, J. Pulp Paper Science, 16(2), J63-72 (1990). Miles, K.B., W.D. May, and A. Karnis, "Refining Intensity, Energy Consumption and Pulp Quality in Two-Stage Chip Refining", Tappi Pulping Conf., Toronto, 681-690 (Oct., 1990). Milito, R., C.S. Padilla, R.A. Padilla, and D. Cadorin, "An Innovations Approach to Dual Control", IEEE Trans. Automatic Control, AC-27(1), 132-137 (1982). Mosca, E., S. Rocchi, and G. Zappa, "A New Dual Active Control Algorithm", Preprints 17th IEEE Conf. on Decision and Control, San Diego, 509-512 (1978). Nahorski, Z., and R.V.V. Vidal, "Stabilization of a Linear Discrete-Time System using Simultaneous Identification and Control", Int. J. Control, 19(2), 353-364 (1974). Nobleza, G.C., "Multivariate Analysis of Newsprint Properties", 41st Can. Chem. Eng. Conf., Vancouver, Paper 53-4, (Oct., 1991). Nordin, J.A., J. Hill, and L. Nelvig, "The SCS Control System for the Control of the Mechanical Pulping Processes", CPPA Process Control Conf., Halifax, 121-133 (June, 1980). Ortega, R., and Y. Tang, "Robustness of Adaptive Controllers - A Survey", Automatica, 25(5), 651-  Bibliography  138  677 (1989). Page, E.S., "Continuous Inspection Schemes", Biometrika, 41, 100-115 (1954). Pearson, A.J., "A Unified Theory of Refining", Joint Textbook Committee of the Paper Industry, Atlanta and Montreal, 1990. Rogers, J.H., S.C. Vespa, N.D. Legault, D.J. Butler, and C. Mills, "Automatic Control of Chip Refining Results from Mill Trials on an Open-Discharge Refiner", Pulp Paper Can., 81(10), 89-96 (1980). Rydefalk, S., T. Pettersson, E. Jung, and I. Lundqvist, "The STFI Optical Fibre Classifier", Int. Mechanical Pulping Conf., Oslo, Session HI: no. 4, 1-16 (June, 1981). Ryti, N., and H. Manner, "An Approach to the Design of a Control Strategy for a TMP Process", Int. Mechanical Pulping Conf., Helsinki, No. 11, 1-23 (June, 1977). Savonjousi, A., R. Stenros, and K. Saarinen, "Automatic Control of the Consistency and Density of Thermomechanical Pulp (TMP)", Control Systems 90, Helsinki, 173-184 (Sept., 1990). Seborg, D.E., T.F. Edgar, and S.L. Shah, "Adaptive Control Strategies for Process Control: A Survey", AIChE J., 32(6), 881-913 (1986). Smook, G.A., "Handbook for Pulp & Paper Technologists", Joint Textbook Committee of the Paper Industry, Atlanta and Montreal, 1982. SOderstrOm, T., L. Ljung, and I. Gustavsson, "A Theoretical Analysis of Recursive Identification Methods", Automatica, 14(3), 231-244 (1978). Solo, V., "The Convergence of AML", IEEE Trans. Automatic Control, AC-24(6), 958-962 (1979). Stationwala, M.I., D. Atack, J.R. Wood, D.J. Wild, and A. Karnis, "The Effect of Control Variables on Refining Zone Conditions and Pulp Properties", Int. Mechanical Pulping Conf., Toronto, 93-109 (June, 1979). Stationwala, M.I., K.B. Miles, and A. Karnis, "The Effect of First-Stage Refining Conditions on Pulp Properties and Energy Consumption", J. Pulp Paper Science, 19(1), J12-18 (1993). Stebel, J.H., and C.F. Aeby, "Optimization Mode Control for Double-Disk Refiner Systems: Design and Experience", International Symposium on Fundamental Concepts of Refining, Appleton, Wisconsin, 30-42 (Sept., 1980). Steraby, J., "A Simple Dual Control Problem with an Analytical Solution", IEEE Trans. Automatic Control, AC-21, 840-844 (1976). Strand, W.C., and N. Harder, "Modelling and Optimization of Full Scale Chip Refining", Int.  Bibliography  139  Mechanical Pulping Conf., Stockholm, 46-54 (May, 1985). Strand, W.C., A.V. Mokvist, and M. Jackson, "Improving the Reliability of Pulp Quality Data through Factor Analysis and Data Reconciliation", Int. Mechanical Pulping Conf., Helsinki, 362-375 (June, 1989). Strand, W.C., O. Ferritsius, and A.V. Mokvist, "Use of Simulation Models in the On-Line Control and Optimization of the Refining Process", Tappi J., 74(11), 103-108 (1991). Tamminen, J., M. Perkola, H. Poutanen, and E. Viljakainen, "Mill Scale Trials with Three On-Line Freeness Testers in a TMP Plant", Int. Mechanical Pulping Conf., Vancouver, 53-60 (June, 1987). Toivonen, H.T., and J. Tamminen, "Minimax Robust LQ Control of a Thermomechanical Pulping Plant", Automatica, 26(2), 347-351 (1990). Vespa, S.C., C. Mills, D.J. Butler, N.D. Legault, and J.H. Rogers, "Automatic Control of Chip Refining", Int. Mechanical Pulping Conf., Toronto, 85-91 (June, 1979). Wieslander, J., and B. Wittenmark, "An Approach to Adaptive Control Using Real Time Identification", Automatica, 7(2), 211-217 (1971). Wittenmark, B., "An Active Suboptimal Dual Controller for Systems with Stochastic Parameters", Automatic Control Theory and Application, 3(1), 13-19 (1975a). Wittenmark, B., "Stochastic Adaptive Control Methods: A Survey", Int. J. Control, 21(5), 705-730 (1975b). Wittenmark, B., and K.J. AstrOm, "Practical Issues in the Implementation of Self-tuning Control", Automatica, 20(5), 595-605 (1984). Wittenmark, B., and C. Elevitch, "An Adaptive Control Algorithm with Dual Features", 7th IFAC/EFORS Symposium on Identification and System Parameter Estimation, York, UK, 587-592 (July, 1985). Wood, J.R., and A. Karnis, "Towards a Lint-Free Newsprint Sheet", Paperi ja Puu, 59(10), 660-674 (1977). Ydstie, B.E., L.S. Kershenbaum, and R.W.H. Sargent, "Theory and Application of an Extended Horizon Self-Tuning Controller", AIChE J., 31(11), 1771-1780 (1985). Ydstie, B.E., "Adaptive Process Control", Chem. Eng. Res. Des., 65, 480-489 (1987).  APPENDIX 1 ONE-STEP OPTIMAL CONTROLLER  Consider the process described by Equation (3.1) (1 + az ~l)y(t) = bu(t-\) + (1 + cz _1)v(0/V  (AL1)  where the parameters are assumed known and constant A one-step ahead prediction or forecast of the process output may be obtained by solving two identities. First, the polynomials E(zA) and F(zl) are obtained by solving the following Diophantine identity (1 + cz -1) - E(z "OCl + az -J)V + z -lF(z -1)  (A1-2)  By equating coefficients of z"1, it is straightforward to show that E(z -1) - 1,  F(z -1) = (1 - a + c) + az -»  (A1-3)  Multiplying Equation (Al.l) by Eiz^Vz1, substituting E(zl)(\+azl)V from Equation (A1.2) and simplifying gives $(t+l |0 = (1 - a + c + az -r)y f(t) + bVu f(t) + v(r+l \t)  (A1-4)  where superscript / denotes filtering by l/Cl+cz"1)- The final term in Equation (A1.4) is a future unknown error which may be eliminated since its expected value is zero. Minimization of the cost function Equation (3.7) is in terms of Vu(0, not Vifit). Therefore, Equation (A1.4) must be modified to separate past known filtered control actions from the present unfiltered control action yet to be determined. Consider the following identity b = G(z -J)(l + cz -1) + z "T(z -1)  ( A 1 -5)  Again, by equating coefficients of z"1 Giz-1) = b,  Hz"1) - -be  Substituting Equation (A1.5) into Equation (A1.4) gives  140  (A1-6)  Appendix 1. One-Step Optimal Controller  %t+\\i) - (l-a + c + az-^yXV + bVutf-bcVuXt-l)  141  <A1/7)  Substituting Equation (A 1.7) into the cost function Equation (3.7) and nunimizing with respect to V«(0 gives the control law v  _ b[(l - a + c + az -Qy '(t) - bcVu '(r-1)]  =  ( A 1 8)  After some algebra, this controller can be written as v  =  _ b[(\-a + c) + az-i] ft' + pd+cz" 1 )  which is the one-step optimal controller shown in Equation (3.8).  (A19)  APPENDIX 2 EQUIVALENCE BETWEEN CCE CONTROLLERS  The constrained certainty equivalence (CCE) controller was given in Equation (3.10) by Vu(t)  =  _ b[w+l)Q + y(t)]  ( A 2 .i)  2  b +p  where the parameters are assumed known and constant. Multiplying both sides of Equation (A2.1) by (b2+p) and expanding the productty(t+l)Qinto its individual terms gives (Z>2 + p)Vu(0 = -b[-aVy(t) + cr(t)+y(t)}  (A2.2)  The system Equation (3.1) may be used to reconstruct the residual error r(r) as r(0 = Vy(0 + aV)fr-l)-&VK(r-l) 1+cz"1  (A2.3)  Substituting Equation (A2.3) into Equation (A2.2) gives  (6 2 + p)Vw(0 = -b -aVy(t) + y(f) + c\Vy(i) + aVy(t-l) 1- bVu(t-l)  (A2.4)  1+cz"  Multiplying both sides of Equation (A2.4) by (1-cz1) and simplifying, after some algebra gives V«(f)  =  _ ^ [ ( l - a + tO + az-1] y(t) 2 1 b + p ( l +CZ' )  which is the controller shown in Equation (3.8).  142  (A2.5)  APPENDIX 3 EXPRESSION FOR FUTURE COVARIANCE  From Equation (3.32), the covariance of parameter b one sampling period into the future is given by (A3.1)  pb{t+2\t+l) = lP(t+2\t+l)l where / = [0 1 0] and the future covariance is given by Equation (3.31) ,r)_  P(t+2\t+l) = ^  P{t*l\t)^{t*lMt+l)P(t*l\i) &+R yvit*l)P(t*l\ty^(t*l)*Rv  (A3.2)  The gradient vj/(r+l) may be partitioned as shown in Equation (3.33) (A3.3)  vj/(r+l) = <J/(/+l) + /Vu(0  where \}>(f+l> was defined in Equation (3.34) and $(/+l) was defined in Equation (3.9). For convenience, the time arguments will be omitted for \j/(H-l), Vu(t) and P(r+110 in the sequel. Substituting Equations (A3.2) and (A3.3) into Equation (A3.1) gives P($rT + lTVu)($ + lVu)P  pb(f+2\t+l) = / $ P-  (y + lVu)P($ir + lTVu)+Rv  & +R I  (A3.4)  After several pages of algebra, this reduces to pMPV+Rj-OPyy  pb(f+2\t+l) = 4:  2  _pbVu + 2lPtyVu + tyPy +Rv Upon making the following substitutions  143  + o\.  (A3.5)  144  Appendix 3. Expression for Future Covariance  d1 = lPtyT  d2 = ty/^+fl,  (A3.6)  d  s =Pbd2~d?  d* = fi,d, Equation (A3.5) may be written as d pb(t*2\t+l) = _ ^ + ot * pbVu2 + 2d1Vu + d2 which is the expression shown in Equation(3.35).  (A3.7) b  APPENDIX 4 VARIANCE OF OPEN-LOOP GAIN ESTIMATE  Define the open-loop gain as Ooo  (A4.1) 1+a,  The parameter estimates are assumed to be normally distributed with mean value 9 and covariance matrix P. Therefore, g„ can also be approximated by a normal distribution with mean value $„ and variance V(#M). To derive expressions for g„ and V^J, let 5oo  f0 =  /o  J\  (A4.2)  b0{l+ar  f^b^\+axr Then, L = E{gJ  (A4.3)  = £{/„} + £{/,}  As shown in Equation (A4.2), f0 and fx are nonlinear functions of the parameters. Therefore, evaluation of the expected values requires and approximation. For this, one can use the following Taylor series expansion Ei&ajb)} = A&,b) 1  da2  "  dadb y"b  db2  (A4.4) "  tj>  Here, the first-order terms are all zero and the second-order terms involve the covariances. Substituting f0 for/and simplifying gives  145  146  Appendix 4. Variance of Open-Loop Gain Estimate  Eif0} = 4 ( l + ^ ) - 1 + 6 0 (l + d 1 )- 3 p a -(l + d 1 )- 2 p aA  (A4.5)  EifJ = 6 1 ( l + ^ i r 1 + 6 1 (l + <3 1 )-V ar (l+^ 1 )- 2 P aA  <A4-6)  Similarly,  The covariances are usually small in relation to the parameter estimates and, therefore, may be ignored in the calculation of $„. Thus, substituting Equations (A4.5) and (A4.6) into Equation (A4.3) and ignoring the covariances gives VV  A  I. - i ^ 1+4*1  (A4-7)  This is the expression shown in Equation (5.14). To evaluate V(£J, one can use the well known relation V(f0 +/,) - V(f0) * V{f) * 20JJ,)  (A4.8)  where C(-) denotes covariance. V(f0) is evaluated as follows. First define l^(/0) = £{(/0-£{/0})2} (A4.9) = £{/02}-£{/0P E{f0} is already known (see Equation AIV.5). Using Equation (A4.4) £{/02} = 35 0 2 (l + ^ 1 )^ r 4^ 0 (l + d 1 )- 3 p a A + (502+^Xl+<31)-2  (A4.10)  Substituting Equation (A4.5) and (A4.10) into Equation (A4.9) and simplifying gives  Wo) - V a + ^ x - ^ o a + ^ p - A + a + ^ X  (A4 U)  -  Similarly,  V(fJ = £2(1 + ^ ) X - 20,(1 + ^)X A + d + <*i>X The final term in Equation (A4.8) may be evaluated as follows. First, from Equation (A4.9)  (M 12)  -  147  Appendix 4. Variance of Open-Loop Gain Estimate  CifJ,) = Ei^-EifJW-EifJ))  (A4.13)  -EitfJ-EifjEW E{f0} and E{fx} are already known (see Equations AIV.5 and AIV.6). isi/o/i} is a function of three variables, i.e. (A4.14)  £{/„/;}-Ety&a+<!,)-*} Again, one can use a Taylor series expansion E{f(a,b,c)}=f{a,b,c) + } - VflpJbfi) „ da.+2  d%a,b,c) * dadb  +  + ah  VflflJb*) n „ VKafifi) dc: db: d2jja,b,c) *\  -\  dadc  •* Ac  (A4.15)  +dfja,b,c) ~v» -\  dbdc  * fc  *Ae  Substituting f0fx for/and simplifying gives  Eifj} = 3b0b1(i+a1rPai-2(b^afii+b1pafip+d1r3HbA+p^p^lr2  <A4-16>  Substituting Equations (A4.5), (A4.6) and (A4.16) into Equation (A4.13) and simplifying gives  cifjt = b0bl{\+a,rp-{b^.+blp.){\^xrHi-alrlpb..  (A4 17) »A  -  Finally, substituting Equations (A4.ll), (A4.12) and (A4.17) into Equation (A4.8) and simplifying gives the following expression for the variance of #„ (A4.18)  which is the expression shown in Equation (5.15).  

Cite

Citation Scheme:

        

Citations by CSL (citeproc-js)

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>
                        
                    
IIIF logo 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-0058582/manifest

Comment

Related Items