Quantitative A d a p t i v e Robust C o n t r o l Based on a Robust Frequency D o m a i n E s t i m a t o r by Hongtao Wang M.Sc. University of Science and Technology of China, P R C , 1985. A THESIS SUBMITTED IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF DOCTOR OF PHILOSOPHY in THE FACULTY OF GRADUATE STUDIES DEPARTMENT OF ELECTRICAL ENGINEERING We accept this thesis as confonning to the required standard THE UNIVERSITY OF BRITISH COLUMBIA October 1992 © Hongtao Wang, 1992 In presenting this thesis in partial fulfillment of the requirements of an advanced degree at the University of British Columbia, I agree that the library shall make it freely available for reference and study. I further agree that permission for extensive copying of this thesis for scholarly purposes may be granted by the Head of my Department or by his or her representatives. It is understood that copying or publication of this thesis for financial gain shall not be allowed without my written permissioa Department of Electrical Engineering The University of British Coliunbia 2356 Main Mall Vancouver, B C Canada V 6 T 1Z4 Date: O - / . 9 9-^ . Abstract It is widely appreciated that robust control and adaptive control focus on dealing with system uncertainty but emphasize different aspects. Adaptive control is efficient in eliminating structured uncertainty but has difficulty handling unstructured uncertainty. On the other hand, robust control guarantees the stability for the worst case but ends up with a very conservative design. This thesis will be mainly concerned with the development of a frequency domain quantitative adaptive robust control strategy, which can be viewed as a good combination of adaptive control and robust control to obtain the best closed-loop performance of a stable closed-loop system. The research will consist of two parts. Part one is to develop a new estimation method called frequency domain robust estimator which can improve the frequency function estimate and quantify system uncertainties in the frequency domain Part Two is to implement the quantitative adaptive robust controller based on the information obtained from the robust estimator. Table of Contents Abstract ii List of Figures v List of Abbreviations 1 2 viii Acknowledgments x Introduction 1 1.1 General Control Review 1 1.1.1 1 Qassical Control Viewpoint 1.1.2 Robust Control Viewpoint 2 1.1.3 5 Adaptive Control Viewpoint 1.2 Motivation and Thesis Outline 7 1.3 Contributions of This Thesis 9 Frequency Domain Robust Estimator 12 2.1 Introduction 12 2.2 Methods of Frequency Function Estimation 13 2.2.1 Empirical Transfer Function Estimation 13 2.2.2 Transfer Function Estimation via Kalman Filtering 16 Estimation of Transfer Function and Uncertainty Bounding Function 24 2.3 2.3.1 A Stochastic Approach in Describing the Modelling Uncertainty for Continuous-Time Systems 2.32 24 Frequency Function and Unstructured Uncertainty Estimation with Laguerre Functions 28 2.3.3 Estimation of Frequency Function and Uncertainty Bounding Function via Recursive Digital Fourier Transform 44 2.4 2.5 3 Least-Squares Algorithm 52 2.4.1 Frequency Domain Recursive Least-Squares Algorithm 52 2.4.2 External Disturbance Spectrum Estimation 54 Conclusions 60 Quantitative Adaptive Robust Control in Frequency Domain 62 3.1 Introduction 62 3.2 Stability of Adaptive Control 63 3.3 Stability of Quantitative Adaptive Robust Control 66 3.4 Robust Performance of Quantitative Adaptive Robust Control 71 3.4.1 Quantitative Feedback Theory 71 3.4.2 Robust Performance Criteria for Quantitative Adaptive Robust Control 75 Implementation of Quantitative Adaptive Robust Control 81 3.5.1 Off-Line Design: Open Loop Shaping 81 3.5.2 On-Line Design: Controller Design 86 3.5.3 Q A R C Design Procedures 91 3.5.4 Simulations and Discussions 92 3.5 4 External Disturbance Spectrum Estimation via Frequency Domain Recursive 3.6 Discussion about Non-Minimum-Phase Systems 102 3.7 Conclusions 106 Conclusions 108 4.1 Summary of Thesis 108 4.2 Suggestions for Further Research 109 References 111 A Estimation of Pure Time Delay 117 B Estimation of A R M A Model from Frequency Domain Estimation 122 List of Figures 2.1 Empirical Transfer Function Estimation 16 2.2 Input, Plant, Noise and Output 17 2.3 Frequency Response Estimate Using Dead-beat Observer 22 2.4 Frequency Response Estimate Using Kalman Filter 22 2.5 Frequency Response Estimate Using DFT with External Disturbance 23 2.6 Frequency Response Estimate Using Kalman Filter with External Disturbance 23 2.7 Laguerre Filter Networic 33 2.8 Frequency Response Estimate by Laguerre Filters, N=3 35 2.9 Frequency Response Estimate by Laguerre Filters, N=6 35 2.10 Estimation of Frequency Response and Bounding with Laguerre Filters, N=3; M=5 . . . . 38 2.11 Estimation of Frequency Response and Bounding with Laguerre Filters, N=6; M=5 . . . . 38 2.12 Estimation of Frequency Response and Bounding for a 2nd Order System with Laguerre Filters, N=6; M=5 39 2.13 Estimation of Frequency Response and Bounding for a 2nd Order System with Laguerre Filters, N=10; M=5 40 2.14 Estimation of Frequency Response and Bounding for Rohrs' model at t=4 sec 41 2.15 Estimation of Frequency Response and Bounding for Rohrs' model at t=20 sec 41 2.16 Estimation of Frequency Response of Rohrs' Model for L=10 (N=10; M=0) 42 2.17 Magnitude Plot of Frequency Response and Uncertainty Bounding for Rohrs' Model . . . 43 2.18 Estimation of Frequency Response and Bounding via R D F T 51 2.19 Estimation by Eq.(2.133) witii Single w Noise 57 2.20 Estimation by Eq.(2.137) with Single oj Noise 57 2.21 Noise Spectnim Estimation for Single w Noise 58 2.22 Estimation by Eq.(2.133) with White noise 59 2.23 Estimation by Eq.(2.137) witii White Noise 59 2.24 Noise Spectnim Estimation for White Noise 60 3.25 Block Diagram of a Self-Tuning Regulator 63 3.26 Error Model Drawn as a Feedback Loop 64 3.27 Inner Loop Block Diagram 67 3.28 Nyquist Plot 67 3.29 Nyquist Band 69 3.30 Block Diagram of 2-Degree-of Freedom Control 72 3.31 Bounds on Response in Tune Domain 72 3.32 Bounds on Loop Transfer Functions 74 3.33 Bound A(u) versus a> 83 3.34 Bounds b(') and T{e^'^) = l/Z(eJ'^) on tiie Nichols Chart 84 3.35 Design Loci b(») and L{e^'^) on the Nichols Chart 85 3.36 Nyquist plots 88 3.37 Magnitude Bode plots 88 3.38 Estimation of Frequency Response and Bounding with Laguerre Filter N=3; M=5 93 3.39 Quantitative Robust Property of Q A R C 94 3.40 Adaptive Property of Q A R C 95 3.41 Frequency Function Estimation in a 3-D representation 96 3.42 Uncertainty Bounding Estimation in a 3-D representation 96 3.43 Performance with a second order controller 97 3.44 Performance with a fifth order controller 97 3.45 Magnitude Plot of Frequency Response and Uncertainty Bounding for Rohrs* Model . . . 98 3.46 Qosed-Loop Step Response from Q A R C 99 3.47 Qosed-Loop Step Response from Regular A R C 99 3.48 Switch From Rohrs' model to 1st order system 100 3.49 Switch From Rohrs' model to 1st onier system with N=2 and M=8 101 3.50 Switch From 1st order system to Rohrs' model 101 3.51 Switch From 1st order system to Rohrs' model with N=2 and M=8 102 3.52 Performance with a second order controller 105 3.53 Performance with a fifth order controller 105 A.54 Time Delay Estimation 120 List of Abbreviations A R C : Adaptive Robust Control; A R M A : Auto Regressive and Moving Average; A R M A X : Auto Regressive and Moving Average Exogenous; DFT: Digital Fourier Transform; E H P C : Extended Horizon Predictive Control; ETFE: Empirical Transfer Function Estimation; FIR: Finite Impulse Response; GPC: Generalized Predictive Control; IDFT: Inverse Digital Fourier Transform; IMC: Internal Model Control; L Q G : Linear Quadratic Gaussian; L S : Least Squares; LTI: Linear Tune Invariant; LTF: Loop Transfer Function; LTR: Loop Transfer Recovery; M I M O : Multiple Input and Multiple Output; M R A C : Model Reference Adaptive Control; N M P : Non-Minimum Phase; PRBS: Pseudo Random Binary Series; Q A R C : Quantitative Adaptive Robust Control; QFT: Quantitative Feedback Theory; R A C : Robust Adaptive Control; RDFT: Recursive Digital Fourier Transform; RHP: Right Half Plane; R L S : Recursive Least Square; SISO: Single Input and Single Output; S O A C : Self-OsciUating Adaptive Control; SPR: Strictly Positive Real; STR: Self-Tuning Regulator, T F L : Target Feedback Loop. Acknowledgments I am deeply indebted to my supervisor. Dr. Guy A . Dumont, for his guidance, kind support, valuable discussions and patience throughout the long process of completing this thesis. I appreciate the advice and helpful suggestions of the members of my supervisory and examining committee. Dr. M . S. Davies, Dr. S. Salcudean, Dr. A . C. Soudack, Dr. C. C. H . Ma, Dr. E . V . Bohn, Dr. K . P. Lam, Dr. F . Sassani and Dr. S. L . Shah. M y thanks are to the members of the control group of the Pulp and Paper Center, and the staffs of the Pulp and Paper Center for their help in many aspects of this research. Among them, my thanks are especially to Dr. Y . Fu for his helpful discussion and questions. The financial support i n forms of the Research Assistantship from Dr. Guy A . Dumont, Fellowship from Pulp and Paper Research Institute of Canada and Scholarship from the Academy of Sciences of China is greatly acknowledged. Deeply in my heart, there is always another part of the world, my family and friends i n my home covmtry - China, whom I have missed for so long and so much. I dedicate this work to them for their love and encouragement Chapter 1 Introduction 1.1 General Control Review 1.1.1 Classical Control Viewpoint In control system design, the goal is to achieve certain desired requirements, such as stability and performance, in the face of imcertainty and disturbance. Although modem control characterized by M I M O (multiple input and multiple output) systems has become dominant in the recent literature, classical control characterized by SISO (single input and single output) feedback control is still perceived as a collection of dependable design tools. In fact even at this time, the great majority of control designs are accomplished by SISO feedback control theory. However, the weakness of classical feedback theory is that the design process is a trial and error process involving graphical methods. Simplicity comes from the fact that performance specifications are only approximately posed in the frequency domain. Intuitively, the classical control design process may be broken down into three steps: modeling, analysis, and synthesis; each of which may be carried out via a combination of time and frequency domain techniques. In engineering practice, however, the three steps arc loosely matched to one another. The design philosophy of classical control first in the time domain involves ordinary differential equations and their associated characteristic algebraic equations with performance specifications given in terms of simple quantities as rise time, percentage overshoot, settling time, steady-state error, etc. Then, specifications on these quantities are translated to the frequency domain as constraints on the pole-zero locations or on the frequency response of the open-loop system. This translation is only approximate and is usually performed by assuming a dominant second-order model to represent the loop transfer function. Stability of the system in the frequency domain is determined either by the location of the closed-loop poles or by the Nyquist stability criterion [Ogata, 54]. The Nyquist-Bode diagram approach to closed-loop system stability and behavior makes the frequency domain method the core of classical control theory [Nyquist, 53; Bode, 10]. It is interesting to note that the application of Nyquist's stability criterion does not depend on the availabiUty of a system model in the form of a differential equation or characteristic polynomial. Furthermore, the form of the Nyquist locus gives an immediate and vivid indication of how an unstable or poorly damped system feedback performance could be improved by modifying its open-loop gain versus frequency in an appropriate way. A very useful contribution due to Bode [Bode, 11] was his rules for the optimum shaping of the loop-gain frequency function. Also, Bode introduced the concept of maximum achievable performance and a quantitative description of the trade-offs by introducing his integral relationships. A design procedure in this spirit is also advanced by Horowitz and coworicers [Horowitz, 31; Krishnan and Cruichshanks, 41; Gera and Horowitz, 22], and will be discussed later in this thesis. 1.1.2 Robust Control Viewpoint The ability of a control system to cope with uncertainty is generally termed as the robustness of the system. Uncertainty in a control system arises either as disturbance or noise, or as a model uncertainty resulting from linearization, neglect of high frequency dynamics, plant parameter variation, plant aging etc. Hence, no nominal model should be considered complete without some assessment of its uncertainties. Various types of uncertainties arise in physical systems and so-called "imstructured uncertainties" are singled out as generic errors which are associated with all design models [Doyle and Stein, 16]. A very intuitive description of unstructured uncertainty is the kind of uncertainty lacking any phase or directionality information. The structured imcertainty, on the other hand, is the uncertainty with phase information, e.g., when it comes from parameter variations of the model. Generally speaking, a control system design based on an unstructured imcertainty gives a more conservative result Robustaess is commonly broken down into two categories. The first is robust stability and is applied when the closed-loop system remains stable for the range of model uncertainty considered. The second is robust performance and is used when acceptable performance is retained in the presence of plant variations and other disturbances. Feedback control systems design is nontrivial in the sense that a trade-off must be made when the design objectives are conflicting. For example, small tracking errors, insensitivity of the output to plant perturbation and disturbance is achieved by increasing the loop gain. On the other hand, increasing the loop gain may cause instability or saturation of the control signal due to sensor noise, etc. These trade-offs, as well as the robusmess characteristics of the control system, are perhaps best revealed in the frequency domain. One reason for this is that frequency response descriptions preserve the operations of system addition and multiplication [Macfarlane, 47]. The evolution of modem control or modem robust control started with an effort to extend classical control design tools and techniques to multivariable systems. Later some real multivariable techniques emerged when linear system theory matured. Thus emerged the standard treatment of the " L Q G " (Linear Quadratic Gaussian) optimal control problem which became the key to the state space treatment of multivariable control [IEEE Special Issue, 1]. Design specifications are posed in terms of a quadratic performance index, penalizing tracking or regulation error and control input actioa Under the mild conditions of stabilizability and detectability of the plant, minimization of this index produces a unique feedback controller which asymptotically stabilizes the system. In general, L Q G has unfortunately been shown not to possess any nice robustness properties, except in the restricted case when the state vector is directly available and thus does not have to be reconstructed from output measurement [Safonov and Athans, 64]. Because of the robustness concem for multivariable control systems, the classical sensitivity theory seems the main focus even in M I M O control systems. The use of singular value defines in a precise way the robustness properties for multivariable systems. It becomes apparent that good performance and stability properties can be achieved by shaping the singular values of appropriate transfer function matrices of the system, and this concept can be viewed as an extension of the loop shaping techniques of the classical control theory. The L Q G cost in frequency domain, or as it is commonly called Wiener-Hopf performance index [Youla, Bongiomo and Jabr, 73; Safonov and Athans, 64], is interpreted as a weighted average of the singular values of die systems sensitivity and complementary sensitivity transfer function matrices. Then loop shaping can be accomplished by judiciously choosing the weighting matrices which will require the knowledge of system and measurement disturbances. Another procedure following L Q G design is the idea of "sensitivity recovery" or "robustness recovery", a better name is perhaps "full-state loop transfer recovery" [Doyle and Stein, 16]. In this approach, the weighting matrices in the L Q G cost are chosen such that the guaranteed robustness properties of full-state quadratic optimal feedback are asymptotically recovered when the state has to be reconstructed [Doyle and Stein, 16]. The combination of "WienerHopf ' and "full-state loop transfer recovery" becomes the LQG/LTR (Linear Quadratic Gaussian / Loop Transfer Recovery) design procedure [Stein and Athans, 65; Athans, 6]. The basic design problem of L Q G / L T R is interpreted as a weighted H2 —^trade-off between transfer functions in the frequency domain [Birdwell, 7]. This problem can be solved with a modified version of the standard L Q G methodology. The burden is to design a target feedback loop (TFL) with realistic stability robustness and performance properties. This involves the selection of the design parameters to shape the singular values of die T F L in the frequency domain. Since the design assumes that any plant perturbation is possible within a spherical region of imcertainty and that all loop singular values are identical at both low and high frequencies, this makes the M I M O system have approximately the same response speed in all directions. L Q G / L T R and, loosely speaking, the use of the singular values within M I M O closed-loop system, give a conservative design. However, when uncertainty is structured, i.e., includes phase or directionality information, singular-values conditions are only sufficient and can be arbitrarily conservative. Much research effort has been devoted to diminish the conservativeness of singular value conditions. To some extent the problem is alleviated with the introduction of a canonical block diagonal perturbation model and block diagonal perturbation stability margins. These margins were designated "strucmred singular values" or —^values" [Doyle, 15]. In this approach, it is in general not yet possible to compute n —^values exactly and real parameter variations cannot be accommodated. A direct synthesis procedure for shaping ^ —^values seems to be rather distant Another branch of robust control is H^o optimal control introduced in [Zames, 74]. The subject is approached from the view point of classical sensitivity theory, with the difference that feedback will not only reduce but actually optimize sensitivity in an appropriate sense. In this approach, the infinity norm or the operator norm of a weighted sum of the system's transfer function matrix is minimized. This objective can be interpreted in the time domain as a minimization of a combination of error and control signals for the worst case out of a bounded class of reference signals [Helton and Sideris, 30]. But the main function of Hoo control is as an algorithmic synthesis tool used to appropriately shape die singular values of the system transfer function matrices. This ensures good performance and robusmess properties. The difference from LQG/LTR is that this task is accomplished in a relatively straightforward manner but with increased computation. Despite the unquestioned success of the Wiener-Hopf, L Q G / L T R and H^o approaches, the classical methods, which rely on lead-lag "compensators" to reduce sensitivity, continue to dominate many areas of design. In fact, all modem robust designs are rooted in the classical sensitivity theory. Under the assumption of linearity, frequency domain methods reveal in a transparent way the various trade-offs among conflicting design objectives. Because robustness properties of the system are transparently revealed, classical control and even further frequency domain methods are in some sense synonymous with robust control. This is certainly the case with SISO systems and it becomes clear that the same principles hold for M I M O systems. 1.1,3 Adaptive Control Viewpoint Adaptive control is a special type of non-linear feedback control in which the states of the process can be divided into two phases based on their changing rates. The fast changing states belong to an ordinary feedback loop. The slowly changing states belong to a second loop or updating loop, i.e., process parameters. In the fast changing state space, existing feedback control theory, such as state-space and classical control theory, can be applied immediately. In the parametric space, stochastic control theory and filter theory contribute very much to the establishment of estimation theory. Both control theories and signal processing theories have contributed to the development and understanding of adaptive control. There are two main branches in adaptive control theory, namely M R A C (model reference adaptive control) and STR (self tuning regulator). In the M R A C approach, the design was originally based on a continuous time and state space model. The design addresses the problem of adaptively adjusting the parameter of a controller to follow a predefined reference model; and at the same time, the stability of the closed-loop system is guaranteed by satisfying the Lyapunov function under very restrictive assumptions [Parks, 59]. The first M R A C can be traced back to the well known " M I T " rule [Whitaker, Yamrom and Kezer, 69]. Further developments of M R A C , such as the stability concept of the design procedure and the treatment of multivariable system has been also rigorously investigated [Landau, 43; Landau and Courtiol, 44; Narendra, Yin and Valavani, 52]. The STR approach has been proposed mainly as a means for tuning digital controllers for industrial processes. There are two basic assumptions to simplify the regulator design. One is the Certainty Equivalence Principle where the unknown true parameter values are replaced by the estimated ones; another is the Separation Principle which is the separation of the process parameter estimation and the control signal computation, i.e., the separation of the two loops. Therefore, the detailed structure of a self-tuning control law can be decomposed into an on-line parameter estimation module together with an on-line control law synthesis procedure. A common characteristic of STRs is that they are based on a discrete time input-output model. The discrete polynomial A R M A X (auto regressive and moving average exogenous) model representation for the dynamic plus noise transfer function has been widely used as the main mathematical tool in designing, studying and analyzing STRs. For the on-line mathematical modelling of dynamic systems, there have been various identification methods developed such as Least Squares, Maximum Likelihood, Instrumental Variables, etc [Ljung, 46]. Independent of each method, the exactness and convergence mostly depends on the richness of the exciting signal. The theory of parameter estimation has another impact on the adaptive system. Eariier STR using Recursive Least Square was introduced by [Âstrôm and Wittenmark, 3] who applied a minimum variance strategy to obtain a self-tuning minimum variance regulator. A more general approach was presented by [Qarice and Gawthrop, 12] to include the control of a stable non-minimum phase system by using a generalized output function and employing a cost function which included a penalty on the control signal. Similar self-tuning algorithms were developed for pole-zero placement control schemes [Âstrôm and Wittenmaric, 4]. Confusion may arise with the terms of explicit STR and implicit STR, or equivalently direct STR and indirect STR. In a direct situation, the controller parameters are updated directly, while for a indirect algorithm, controller parameters are derived indirectly via a design procedure. It is however, sometimes possible to re-parameterize the process model such that it is possible to use either a direct or an indirect controller [Âstrôm and Wittenmaric, 5]. The use of a long range predictor in adaptive control known as E H P C (extended horizon predictive control) [Ydstie, 72] was first considered for SISO systems. Developed from the same idea, the G P C (generalized predictive control) strategy was introduced to minimize a finite time horizon quadratic cost function of the function errors and control increments [Qarke, 13]. A s a further development of the SISO stochastic control system design, significant progress has been made to extend the SISO design to the M I M O case [Dugard, 17; Goodwin, 25]. Hao theory has been used to guarantee the asymptotic stability [Polak and Salcudean, 61; Grimble, 27]. The strategy based on Hoc, worst case design increases the robust consideration in adaptive control [Polak and Salcudean, 60]. To overcome the problem that a precise a priori knowledge of the model structure such as system order and time delay is required, different approaches have been launched to find a more robust structure than A R M A representation. Impulse representation of the process [Guimarson and Ljung, 29] and sinusoidal analysis [LaMaire, Valavani, Athans and Stein, 42] have been used to obtain the transfer function estimation . Parameter identification of linear systems via orthogonal functions has been imder study for quite a long time [Zervos, 77]. Among them, an appropriate, simple, convenient orthogonal set is the Laguerre set that has been developed into an adaptive regulator and tested in practice [Zervos and Dumont, 76]. The estimation of model uncertainty is an active area of research in frequency domaiit Many researchers have developed techniques for estimating the probabilistic bound, i.e., the soft bound for model uncertainty [Zhu, 78; Rivera, Pollard, Sterman and Garcia, 62]. Just as robust control, adaptive control is simply another way to diminish the system sensitivity to the modelling uncertainty. Unlike robust control which deals with unstructured imcertainty, adaptive control has focused on the parametric modelling uncertainty or structured uncertainty. Research has shown that regular adaptive controllers can not adequately handle unstructured uncertainty [Rohrs, Valavani, Athans and Stein, 64; Rohrs, 63]. This explains why the behavior of adaptive controllers in the presence of unmodelled dynamics has been a topic of concem in the adaptive control community in recent years [Fu and Dumont, 21; loaimou and Tsakalis, 36; Ortega and Yu, 57]. 1.2 Motivation and Thesis Outline One of the important design criteria in feedback control theory is the method of dealing with system uncertainty. A well known approach for dealing with unmodelled or unstructured uncertainty in control system design is the approach of robust control. In this approach, one aims to design a fixed control law which trades off the performance to achieve robusmess and particularly to preserve stability in the worst case design. To achieve this result, robust control design requires information on the nature of the imcertainty. In recent years, the requirement of uncertainty description becomes another starting point for adaptive control theory. Goodwin [Goodwin and Salgado, 24] states that "There should be a consistency between the kind of information used in robust control design and that used in adaptive control design." Kosut [Kosut, 39] observes that "The model estimator should be designed to be compatible with the same class of uncertainty anticipated by the robust control design rule." However, the parameter estimation procedures in common use in current adaptive controllers yield an estimated nominal model without any measuring of modelling uncertainty. The adaptive control design is then based on the estimated nominal model alone which is inconsistent with the robust control philosophy. Based on the above argument, our ideal approach should be able to combine adaptive control and robust control and to build a bridge between these two design methodologies. The inspiration from this idea becomes the main motivation of this research. In this direction, the estimator is required to provide information on both the estimated nominal function and the associated uncertainty function. The expected performance of the resulting controller should be at least as good as a fixed robust design. There are many ways of specifying the uncertainty about the model or the obtained system function. Obviously the tighter the uncertainty, the more useful the function. A n appealing approach is to describe the system function and the uncertainty in the frequency domain because working directly with the frequency domain data overcomes the need to express the input-output data in terms of a model and could be more natural since many design objectives are most effectively expressed in terms of the frequency response of the system. For many problems, it is in fact much more interesting and revealing to know the quality of the Bode-plot estimate than the quality of some parameters. There is thus motivation to do the estimation in the frequency domain. Another motivation for modelling in the frequency domain is the possibility of estimating the time delay. Most of the commonly used procedures for estimation and adaptive control assume knowledge of the process delay. However, xinknown pure time delays occur frequently in process control. Thus there has been continuous interest in ways of dealing with unknown and variable time delays. This thesis is organized to present a complete adaptive robust control strategy for minimumphase, linear-time invariant systems, and related theoretical supporting material. Section 2.2 will introduce several interesting methods in frequency domain estimation which appear in recent literature. Section 2.3 wiU study two methods in estimating system frequency function and associated imcertainty bounding. Section 2.4 will develop a way of estimating the disturbance spectrum at the expense of collecting more input/output data. Section 2.5 will give a comparison among the methods in this Chapter and a brief conclusion. Two ideas based on frequency domain robust estimation are relegated in Appendix since they still need fiirther study. One is the estimation of pure time delay and another one is the estimation of A R M A model from frequency domain data. Section 3.2 will investigate the stability of adaptive control system. Section 3.3 will derive the main robust stability theorem for the Q A R C (quantitative adaptive robust control) system. Section 3.4 wiU put together the related robust performance criterion for the Q A R C . Section 3.5 will develop and implement the Q A R C for both off-line and on-line design to make use of the information provided by the rolwst estimator. Section 3.6 will discuss non-minimum-phase systems. Section 3.7 wUl give the conclusion for this Chapter. The overall conclusion is given in Chapter 4 and also there are some ideas about the further development of this research. In particular, we will look at how our strategy can benefit fitom parallel computing, an area under extensive study now. 1.3 Contributions of This Thesis Robust considerations and real-time adaptations as parts of a unified controller design strategy are considered to be superior to the traditional robust control design or adaptive control design. To achieve this objective, the robust estimator presented here will provide a goodness measurement—an uncertainty bounding for the current estimate of the system transfer function. Furthermore, the derived estimator can supply the information about the system environmental disturbance and system time delay. These properties contribute a new perspective to the estimation theory and is essential when the estimate is to be used in an adaptive robust control situation. Frequency domain analysis provides the chance to distribute the estimation error by attaching different weight to the different frequency zones. However, a more common practical situation is that the model has a lower complexity than that of the real system. From the point of view of nominal model estimation, the closest match between the estimate and real system can be achieved around the frequency range of interest. It is important to have an accurate estimation in certain frequencies to help determine the trade-off between the performance and robustness in controller design. On the other hand, in the time domain approach, the above objectives are most difficult tasks to accomplish. In the literature, the phrase "adaptive robust control"-ARC and "robust adaptive controF'-RAC have been used more and more in the last decade [Cusumano and Poolla, 14; Ortega and Yu, 58]. The main idea of A R C and R A C is to investigate the robust property of the adaptive control. This is the same as in classical control where phase and gain margin are tested for each lead-lag compensator. To satisfy a given margin of any meaningful robust measure where closed-loop system remains stable, the design procedure has to be a trial and error process. This becomes much more difficult for adaptive control because a trial and error procedure does not fit into on-line controller design. A n original idea to improve this problem is to pass the specific modelling accuracy information into the control design rule. However, so far it has been unsuccessful because the time domain estimator fails to provide modelling uncertainty information. This thesis, by means of a frequency domam approach, will combine the advantages of classical sensitivity, robust synthesis and adaptive control theory. The closed-loop design will place a nontrivial robust property over the adaptive control; precisely speaking, it is an adaptive property that is added to the quantitative robust control. Therefore, the title of this thesis is "quantitative adaptive robust control"-QARC instead of A R C or R A C . In addition to the estimation of frequency function and uncertainties, another part of the research is to adapt an on-line frequency domain quantitative synthesis design method in order to make use of the modelling accuracy information supplied by the robust estimator. If the derived bound is a hard bound on modelling imcertainty, then die derived robust stability theorem and performance criterion provide guarantees for the stability of the Q A R C . The absolute certainty about derived hard bound relies upon absolute certainty in our priors. In order to derive a useful, i.e., small hard bound, we use some prior assumptions such that the unmodelled uncertainty is bounded by a decaying exponential [Middleton and Goodwin, 50; Goodwin and Salgado, 23]. According to Ljung [Ljung, Wahlberg and Hjalmarsson, 46], we can mark the derived hard bound with an imfalsilied label since, with a reasonable number of filters or data sets, the bound holds for the systems under study. The end product of the research is an adaptive robust controller that slowly learns and produces a successively better Q A R C compensator which is able to maintain the best closed-loop performance out of a stable closed-loop system. So far, this research is restricted in dealing with stable linear lumped time-invariant systems. The linearity is required by sinusoidal transfer function analysis where different frequencies pass the system independently of each other. A minimum phase condition is imposed in the design of Q A R C due to the possible pole-zero cancellation. Chapter 2 Frequency Domain Robust Estimator 2.1 Introduction In specifying the area vmder consideration a distinction has been drawn between identification mediods which estimate points on an unparameterized system model and algorithms which estimate parameters of some structured system model. This split between parametric and non-parametric identification allows a clear distinction to be made between methods which incorporate the model validation problem and those that do not. In the case of non-parametric identification one can isolate frequency and impulse response estimation as the basic techniques. Though some structured models like orthogonal filters may be used, our scope is reduced to merely the frequency response estimation. The frequency domain characterization of system dynamics has, like its time domain counterpart, a long and varied history. A s far as control engineering is concerned, however, frequency domain identification gained deep relevance with the development of stability and design methods based upon frequency response measurements. Historically, frequency response estimation began with the technique now known as transfer fiinction analysis. The sinusoidal transfer function analysis is, without reservation, a most robust and practically useful control system identification method [Wellstead, 68]. In part this is due to its intimate relationship with the Nyquist stability criterion as well as the recent development in Hoo theory and robust control synthesis methods. In the field of adaptive control, the criticism of frequency estimation is the lack of the real-time implementing power. As a result, frequency identification has lost popularity in adaptive control because the realtime control synthesis and design tools rely heavily on parametric system models. Our effort in this chapter is devoted to the development of several algorithms which can estimate the systems frequency response in real-time. Later, we will proceed to the estimation of uncertainty bounding function, noise function and dead-time. 2.2 Methods of Frequency Function Estimation 2.2.1 Empirical IVansfer Function Estimation It is difficult to ascertain exactly who first proposed the use of numerical Fourier transforms as a tool in the analysis of random data. However, following Ljung's woric [Ljung, 46], we can obtain an infinite dimensional frequency function describing a given plant through Fourier transforms. In order to be consistent with the literature, we introduce the following notations: • Grie^'^) : true frequency function at frequency w; • GN{e^'^°) '• estimated function for a given input frequency Uo', • GN{e^'^) ' estimated function for multi-frequency input u. where: subscript N means that N pairs of input-output data are used. Consider a linear, stable, time-invariant system with input signal [u(t)} and output {y(t)}, subject to an additive disturbance {v(t)}: oo y{t) = ^9T(kHt-k) k=i + v(t) ; i = l,2,..., (2.1) This system has the transfer function: Grie^'^) = f^gT{k)e-^'^ k=i (2.2) A n immediate and well-known consequence of the linear structure of Eq.(2.2) is that an input sinusoid of a given frequency gives, after a vanishing transient, an output sinusoid of the same frequency, but with an amplitude change and a phase shift that is described to be the transfer function or frequency function: u{t) = Acoswot y{t) = B cos {ojot + V') + transient = IGrie^'^') B = \GT{e''^'')\A (2.3) (2.4) (2.5) (2.6) It may sometimes be difficult to determine ip and B from inspection of {y(t)}, in particular i f the system is subject to disturbances. A commonly applied method is then to correlate the output with cosugt and sinoJot respectively: N N ^c(^) = S y W c o s i y . i t=i ; Is{N) = Y, vit) sin Wot f=i (2.7) Then we take: =arctaii{-/,(iV)//,(iV)} lÔNie''^'') NA/2 (2.8) (2.9) A more detailed account is given in [Ljimg, 46]. In this part, we work with the following Fourier transform: N UN{M = -J=T 1 YNUO,) u{t)e-^'^' (2.10) ^ =-^Y,y(^y~'^' <2.ii) Based on the observations: u{k),yik); fe = l , 2 , . - - , i V and where T, is the time between each observation. B y assuming Wo is one of those rearranging Eq.(2.7), we find that: h{N)-jIsiN) = VNYNiJu;o) with YN{JOJO) defined by Eq.(2.11), then Eq.(2.8, 2.9) can be rearranged as: G Nie"'") = 2VNYN{JO;O)/NA Since input is given by Eq.(2.3), we know as in Eq.(2.10): UNiM) = VNA/2 and slightly Consequently, frequency analysis by the correlation method actually employs the following estimate: GN (e^"° ) = YMUOJ,)IUN{3<^O) In Eq.(2.12), (2.12) is the frequency of the input sinusoid. In a linear system, different frequencies pass the system independently of each others. It is therefore quite natural to extend the frequency analysis estimate Eq.(2.12) to the case of multi-frequency inputs, that is, we introduce the estimate: GN{e^'')=YN{jo:)IUN{j<^) (2.13) with YN{JIJJ) and UN{J<^) defined by Eq.(2.10, 2.11). This estimate results from a direa application of Fourier analysis and will be called the empirical transfer function estimation-ETFE since no assumptions have been imposed odier than linearity of the system [Ljung, 46]. To demonstrate this technique, the following first-order system is used as a plant to produce VO data: G{e^'^) = . " ' ^ ^ ^ eJ'^ - 0.9 (2.14) where e-''^ = ^ is a commonly used notation. In our simulation, the input-output model is: y{t) = 0.9j/(t - 1) -I- 0.1tt(t - 1) (2.15) where [u(t)} is filtered signal with normal distribution N{0,1}. Both [u(t)} and (y(t)} pass through a Harming window before entering the DFT block. The E T F E Eq.(2.13) of the system is shown in Fig 2.1 where the transfer function is represented by gain and phase plots. The problem here is that the variance of the estimate does not decay with N (N=256 in this example) but remains equal to the noise-to-signal ratio at the corresponding frequency. This property makes the empirical estimate useless in most practical cases [Ljung, 46]. It is easy to understand the reason why the variance does not decrease with N . We determine as many independent estimates as we have pairs of data points. The information per estimated quantity does not increase as more data is collected. In other words, we have no feature of data and information compression. This in tum is due to the fact that we have only assumed linearity of the true system. So it leaves two possibilities _41 1 • J ' '—I—'—'—'—' lOO ' ' ' '—-J—*—*—^ 101 102 DFT Figure 2.1: Empirical Transfer Function Estimation to increase the information per estimated quantity. One is to assume that the system's behavior at one frequency is related to that at another, this is the starting point of the proposed work on nominal model fitting in frequency domain. Another is to accumulate the information at each frequency but at different times; this is the use of iterative methods and uncertainty bounding in the estimation to complete the data and information compression. 2.2.2 Transfer Function Estimation via K a l m a n Filtering The Kalman filter is a well-known recursive filter used for optimal state estimation or parameter estimation. We will show the same filtering method can be applied to optimal spectral estimation. The superiority of the Kalman filtering approach in estimating the transfer function is its real time capability and its noise smoothing effect Just like the estimation in ETFE, the computational burden is to perform the harmonic estimation, i.e., to estimate [Û{joJk) and Y{ju>k) ; k = 1,2,...,N] from sampled signal {u(t) and y(t) ; t=l,2,... }. Instead of performing the DFT block data transform, the spectral estimation can be done in a recursive manner. This idea can be traced back to a paper V(t) u(t) G(z) Figure 2.2: Input, Plant, Noise and Output by [Hostetter, 35]. We can show that Kalman filtering is optimal in the sense that noise is rejected. This in tum improves the transfer function estimate. We consider input and output signals as shown in Fig 2.2. Again, G(z) is a linear, stable, and time-invariant system. The input signal is assumed to have at least all harmonic components at Wfc. The method starts with the construction of the state-space equations. Followed by the Fourier expansion of input and output signals: u{t) - + [a{2k)co&ijJkt + Oi{2k -f- l)sinwfci] -|- v„(t) (2.16) k=i KO = ^ + E [/?(2A;)coswfct + (5{2k + 1) smi^kt] + Vy{t) + v{t) Jt=i (2.17) where: Vu{t) and Vy{t) are the truncation error, v{t) is the system disturbance accumulated at output Our problem at first is almost like an ordinary Fourier transform, 0^(0» Pi') i.e. how to estimate o:(-), /?(•). Taking Eq.(2.16) as example, the signal u(t) may equivalently be described as the output of die signal model: x(t + l) = ^x{t) (2.18) u{t) = Exit) + v^{t) where, $ is a ( N - l ) x ( N - l ) matrix and i f is a l x ( N - l ) matrix: ri coswi sinwi -sinwi coswi C0SW2 sinw2 — sin 0^2 C0SW2 cosa;£i._i sma;ii_i -,i,o,i,o,---,i,o and x(t) is a state vector as follows: a(l) a(2)cosu)it + a(3)sinwi< —a(2) sinwif + a(3) coswif a(4) cosa;2t + a(5) sina;2< (2.19) a(7V - 2)COSWN_I< + a{N — l)sina;Ai_if . - a ( i V — 2)sinwN_jt + «(iV - l)cosa;N_it. For y(t), there is a similar description: z ( i + 1) = (2.20) y ( 0 = i f ^(0 + w{t) where: m /3{2) cos wi t + I3{3) sin wi t -/?(2)sina;if + /?(3) coswi^ /?(4)cosw2t + P{5)smu)2t (2.21) z{t) = /?(iV - 2) cosw|:_it + /?(iV - 1) s i n w ^ _ i i .-I3(N - 2)smu>N_it + P{N - l)cosuN_^t, W(t) = Vy{t) + vit) The fact is that U{jijJk) and Y{jiJk) are given by U{3Uk) = a{2k)-ja{2k + l) Y{M) = (3{2k)-Jmk + l) (2.22) A;=l,2,.-.,^-l A n important observation from the above equations is that the Fourier transform is direcdy related to the value of the state vector when t = ~ where p is an integer. Therefore, U{M\t) = X2k{t)-jX2k+lit) YiM\t) Z2kit)-jZ2k+l{i) = (2.23) fc= l , 2 , - - - , f - 1 In applying the Kalman filtering formula to the Fourier analysis problem denoted by Eq.(2.18) and (2.20), it is easier to consider the unknown w{t) as zero mean white measurement noise. Taking Eq.(2.20) as an example, we have the following standard solutions: zit + 1) = ($ - Ktn^)z{t) + Ktvit) Kt = Pt-iH{H'^Pt-iH + R) -1 (2.24) (2.25) (2.26) where R, the covariance matrix of w{t), in this case, is only a scalar. Given initial conditions, the above formulas allow us to implement the optimal state estimator in terms of minimizing the covariance of the state estimation error. A simple form is derived by assuming R=0. Kalman filtering becomes an ordinary state observer problem. Then Kt in Eq.(2.25) is a fixed N - 1 vector observer gain {K}. The state estimate z{t) converges to {z(t)} geometrically according to die eigenvalues of ( $ - KH"^). And fiirther, since is a completely observable pair, these eigenvalues may be arbitrarily placed by choosing {K] properly. A deadbeat observer is obtained by placing all the eigenvalue of ($ - KH'^) at the origin. Then the z{t) converges exactly to z(t) in N - 1 steps and yields exactly 2N+1 D F T components as a simple fimction of the elements of the state vector Eq.(2.24). In the other end of stable range, the characteristic equation of the system matrix # of Eq.(2.20) is N (A - 1) J] ( A | + 2Afc coswfc + 1) = 0 (2.27) k=2 Following the solution to the quadratic equation A? -|- 2\k coswjt + 1 = 0, i.e.. Afc = ^^—2cosu;fc ± y/Acos^Wf. - 4^ = — coswfc ± J sinwfc we can see that all eigenvalues of $ , A;^ = e*-''^*, are equally spaced around the unit circle. If we consider the signal models involving additive measurement noise, the fast convergence rate of the deadbeat estimator implies poor noise smoothing feature. Conversely, slow convergence rate implies good noise rejection but poor tracking of the time-variations of the harmonic components. Eventually, any block-processing spectral estimation method will have similar noise problem to the deadbeat observer because the number of data points is usually quite short when compared to the coherence length of the signal. Consequendy, it is beneficial to consider methods which allow observers to be designed to trade-off these opposing objectives. The Kalman filtering solution is an attempt to compromise between these two extremes of maiginaUy stable Eq.(2.27) and deadbeat observer, and also take into account the minimization of E[z{i) - Since the signal model (2.20) is time invariant, and the covariance of {w(t)} is assumed beforehand, here we only consider utilizing the steady-state time invariant version of the Kalman filter. We introduce two theorems to guarantee the performance of the filter. Theorem 2.1 [Hostetter, 35]: if] is a completely detectable pair, then there exists a unique maximally nonnegative definite solution P to the steady-state Riccatti equation P = $ P-PH{H^PH + R) ^H^P Theorem 2.2 [Bitmead, Tsoi and Paricer, 9]: Consider the same pair (2.28) H] in Theorem 2.1. Let R = r > 0 and let P = el for some scalar e>0. Then defining K = ^PH (H^PH + R)~^ (2.29) the matrix ($ - KH^) has all its eigenvalues strictly inside the unit circle. Now we derive a solution for a state observer by using the concept of Kalman optimal filter. This method leaves us a user-adjustable parameter e. With R>0, it is easy to see: when c - * oo, we have the deadbeat observer; and e 0, we have marginally stable observer. Convergence rate and noise rejection are traded off against each other by choosing different design parameter e. The method finally leaves us a real-time harmonic estimator which contains the result of the Kalman filter formalism: ÛiJCOk\t) = X2k{t) - jhk+l{t) Yijojk\t) = Z2k{t) - jZ2k+i(t) (2.30) fc=l,2,--.,y-l To estimate the transfer function, we simply perform the division as mentioned in E T F E : ^ ( e ^ - | 0 = | ^ (2.31) The index "t" in the left side reflects the iterative or real-time property of the algorithm. To show the performance of the Kalman filter based estimator and the effect of adjusting parameter e, the same 101 100 w=2*PI*x/256 Phase 101 Kalman Filtering Figure 2.3: Frequency Response Estimate Using Dead-beat Observer 101 -T 1 Gain 1—I -1 1 1 1 100 10-1 10-2 -1 I I I I 101 Kalman Filtering Figure 2.4: Frequency Response Estimate Using Kalman Filter 1 1 1 T; Gain DFT With Noise (N/S=0.1) Figure 2.5: Frequency Response Estimate Using DFT with External Disturbance 101 -1 Gain 1—I—I—I—I—I— 100 10-1 J -4 100 101 1 I L. 102 K.F. With Noise (N/S=0.1) Figure 2.6: Frequency Response Estimate Using Kalman Filter with External Disturbance plant Eq.(2.15) from last section is used. Fig 2.3 is the estimate from above procedure with €=1 and r=0.1; in this case the Kalman filter is close to a dead-beat observer. As expected, the estimate is almost the same as the one using D F T which is shovm in Fig 2.1. By reducing e, we can see the noise rejection or smooth feature gained from Kalman filtering. With e=0.001 and r=0.1, the improvement of the estimation is shown in Fig 2.4. Furthermore, we add disturbances to the system Eq.(2.15): y{t) = O.dyit -1) + 0.1«(< - 1) + 0.01e(0 (2.32) where e(t) is a PRBS signal between - 1 and +1. The usual estimate from DFT deteriorates as shown in Fig 2.5. A much better estimate is obtained from the Kalman filter as shown in Fig 2.6. In the graph, the solid lines are the system response; the dashed lines are estimates. A later section of this chapter wiU deal with the problem of division by zero in Eki.(2.31) and show the comparison with other frequency domain methods. 2.3 Estimation of Transfer Function and Uncertainty Bounding Function 2.3.1 A Stochastic Approach i n Describing the Modelling Uncertainty for Continuous-Time Systems There are many problems with the adaptive control algorithms which have been developed to date. In particular, most adaptive control algorithms are not robust to unmodeUed dynamics and non-measurable disturbances, particularly in the absence of a persistentiy exciting input signal. So undermodeling has drawn attention because it is a crucial problem i f we want to make the adaptive control more robust In the following description, the undermodeling uncertainty is strongly related to the input characteristics because the input frequency dependence of the estimated nominal model is not removed. There are many ways to describe the uncertainty. The most convenient way is to describe it in the frequency domain. This allows us to specify the degree of smoothness, or extrapolability of the true system transfer function. In this section, we use the following notation for continuous-time systems: • GrCJw): true system frequency response at frequency w, • GoiJoj): nominal model frequency response at frequency w; • G{joj): estimated model frequency response at frequency w; where: GoU<^) = B{eo,MiA{eo,3u) = [tm,i'm-i,---,to,-a„_i,-a„_2,--->-ao] A method to describe the uncertainty has been recently proposed [Goodwin and Salgado, 23]. It defines the difference between GT{JU)) and Go{j<^)'- AGiJo;) = GriM - (?o(jw) (2.33) as a particular realization of a stochastic process, AG'(jw), in the frequency domain. The class of possible functions, AG(jo;), is described in terms of the second-order statistics o f the associated stochastic process. The mean is taken to be zero and the covariance function to be E{AG(ju;)} = 0 E{AG{M)AG'Uo;2)] = i?K,W2) , u>2), that is: (2.34) (2.35) where the asterisk denotes the complex conjugate. We can associate a time domain process h{t) defined in [0,oo) with AG{ju>) via the Inverse Fourier Transform under the condition: oo oo J j \R{wi,U2)\d^\diJ2<oo (2.36) —oo —oo The process h{i) has a zero-mean and satisfies: E{h{h)h%h)] = r{tut2) (2.37) where r ( t i , t 2 ) and Riuji^u}^) are related by a 2-D Fourier Transform, that is: oo oo '•(<i,^2) = 7 - ^ / / i2(a>a,a;2)e^'^''V'^''^dwidu>2 (27r) J J —oo —oo (2.38) to ensure that h{t) is a real valued process, or equivalently, that AG{ju)) = AG{-ju>), we require that: i?(wi,W2) = i Z ' C - W i , -wa) (2.39) A case of interest is when the frequency domain process is a wide sense stationary process. L i this case, the covariance function R(u;i,oJ2) depends only on the difference in frequencies wi and « 2 , and we define: R(u) = R(ui - W2 ) = R(<^i, W2) (2.40) A particular choice of R{u) which satisfies the requirement in Eq.(2.39) is: The corresponding time domain impulse response is: f(t) = a^e-^' (2.42) This is consistent with the common assumption made in recent literature on adaptive control [Middleton, Goodwin, Mayne and Hill, 50; Goodwin and Salgado, 23] that the urmiodelled dynamics are described by an impulse response which is bounded by a decaying exponential. We express the system in regression form: yit) = <f>{tfGo + r){t) (2.43) where: <f>{tf = [p'^uit), u{t),p^-'y{t), • • •, y{t) p denotes the differential operator and i]{t) denotes the modelling error which is related to the input by the frequency domain transfer function: 7?(0 = A{6o,jw)AGiJu)u{t) (2.44) The non-recursive ordinary least-squares estimate is: T ê = J 4>{t)y{t)dt (2.45) 0 T 0 = 0-00 = P^ J 4>{t)rj(t)dt (2.46) where: -I -1 j <i>{t)<j>^{t)dt i L 0 Our central interest here is in the difference between G(jw) and GT{J<^)- A preliminary expression for the square magnitude of this difference is given in the following result: G-Gi = \AGf + G-G. -2Re (G-Go)AG* (2.47) where: G-Go = dGo \p=3<^ de \tèP 1 j m n m where: À = A{ê,p) Ao = A{9o,p) The proof can be found in [Goodwin and Salgado, 23]. The result given in the above for G-Gi 2 depends upon the particular realization of A G . When applying the formula we remove this dependence by finding the 'average' of the squared error over the class of admissible functions A G : G — GT = R(oj) + E G - G , -2ReE (G-G)AG*' (2.48) The effects of limited model complexity on estimation are quantified based on the assiunption about the covariance function for A G . It is true that the derived difference is not a hard bound. However, it indicates how a posteriori errors are distributed in frequency as a result of the estimation procedure. Along the same direction, some other fimctions for A G have been studied and still result in a soft bound [Middleton and Goodwin, 49]. In the frequency domain, the nominal model Go is difficult to define. So, one way to give the bounding function for A G is to use a stochastic description. If we choose another expression, for example, Laguerre functions for Go, we can get a deterministic description by estimating the bounding for j G j - Go| directiy in the frequency domain. This is done in the next two sections. First, we will try to use another structure to define the nominal model Go and to estimate A G . 2.3.2 Frequency Function and Unstructured Uncertainty Estimation with Laguerre Functions The use of orthononnal functions for time-domain approximations goes back as early as the use of Fourier series. But the use of orthonormal functions for frequency domain approximation has not really been exploited. As we know, die traditional impulse function is an orthonormal function representation in frequency domain. For a causal system, the impulse response may exist over the whole positive time axis. The traditional approach of expanding transfer function in the delay operator to obtain linear in the parameters predictor models leads to approximations of very high order in case of rapid sampling and lai^e dispersion in time constants. Therefore, it is very useful to look at other orthonormal fimctions for obtaining frequency domain approximations. There are several reasons for selecting the Laguerre functions here instead of other orthonormal functions. In particular, they have a convenient network realization, a flexible stmcture, are similar to transient signals and to Fade approximants in identifying time-delay [Zervos, 76]. With Laguerre networics, each filter has longer memory than the delay operator of one sampling step. The expansions then need much lower order to obtain reasonable approximations, and improve the numerical properties of the estimation algorithm. In the time domain, the Laguerre function are described by /i(0 = X ^ ( 7 ^ | ^ r ^ - - ' ^ ' ) p>0, (2.49) i=l,2,---,N where p is the time-scale. These functions form a complete orthonormal set in the time domain [0,oo) and the corresponding Laplace transform for this set is Fi{s) = ^p^tzïÇ. (2.50) i=l,2,---,iV The orthonormal property from the time domain is preserved in the s-domain: N F{s) = Y,riFi{s) t=i (2.51) is the Laplace transform of N /(0 = E ^ ' / ' W (2.52) f(t) is a real and continuous function in L2[0, oo) which is defined as a closed subspace in the square integrable Lebesgue space [Zervos, 76]. We still need a theorem to show in what space the s-domain Laguerre functions exist, or where they form a complete set. First of all, let us restate a theorem from [Francis, 18]. Theorem 2.3 [Francis, 18]: The Fourier transform is a Hilbert space isomorphism from I'2[-oo5 0 o ) onto L 2 . It maps L 2 [ 0 , 0 0 ) onto H 2 and L 2 ( - o o , 0 ) onto H ^ . Where Hilbert space isomorphism means that the mapping is continuous, preserves norms, is injective, and has a continuous inverse. H 2 is the Hilbert space of all functions G(s) which are analytic in Re[s}>0, take values in C " and satisfy the uniform square-integrability condition. H 2 is the orthogonal complement of H 2 . The following relations are useful: L 2 ( - o o , 0 ) = L^[0,oo) (2.53) L2(-oo,oo) = L2(-oo,0)®L2[0,oo) (2.54) L2 = H 2 © H ^ (2.55) By substitute s=jw, this important theorem says in particular that H 2 is just the set of Laplace transforms of signals in L 2 [ 0 , 0 0 ) , i.e., of signals on i > 0 of finite energy. Theorem 2.4 (the existence or completeness theorem): Fk{s) = (frf ) Let the fiinction G{s) € H 2 and ' ^^'^for any p>0, we can express N G{s)^ Urn TrkFkis) N-KX> (2.56) ^—^ Proof: Because the Laguerre functions form a complete orthonormal set in the L 2 [ 0 , 0 0 ) space [Zervos, 76], given a real and continuous function g{t) e L2[0,oo), there exists an integer N (see Ekj.(2.52)) and a real number e > 0, such that: CO / 0 m - mfdt < s (2.57) This completeness is preserved in s-domain by Theorem 2.3 or Parseval theorem [Oppenheim and Schafer, 55], i.e., we have JOO JOO 27r ^ j \G{s) - F{s)Us < e-, -JOO (2.58) s=ji. TV where F(s) = 53 fiFiis) and G(s) is the Laplace transform of g(t). According to Theorem 2.3, t=i G{s) e H 2 , and Fi{s) forms a complete orthonormal set in H 2 space. So the proof is complete. ### The discrete expression is obtained by following bilinear transformation which maps the left half plane onto the unit disc. z{l-a) + (l-a) s = pzil + a)-{l + a) ^ s(l + a) + pjl - a) ^ s{l + a) - p{l - a) (2.59) where \a\ < 1. This is actually the same idea as [Wahlberg, 66], except diat in that paper the transform is completed in two steps. After the transform, we have _ Jk-i (l-Ha)(^-l)/l-a^Y y/H^iz-a) \z-aj (2.60) which is the prototype of our discrete Laguerre base fimction. The bilinear transform preserves die Hoo norm [Wu and Gu, 71] and is a conformai mapping of the whole left plane to unit disc since \a\ < 1. So the orthogonality is preserved. The problem is that the transform does not preserve the H 2 norm. So the actual space, where Eq.(2.60 ) forms a complete orthogonal set, is smaller than H 2 . Here is the explanation: Notating P{z) = WF^oo = N r , i ^ ( 2 ) to correspond Eq.(2.51), so that. / .(l-a) + (l-a)N [Pzil + a)-il + a), = \\Piz)\ suppose: ,/ ^(l-a) + ( l - a ) \ V ^ ( l + « ) - ( ! +a); then FB.2 C H j [Francis, 18], i.e., {^(^) ••p(.)eH J C H2 A n interesting fact is that the H 2 norm is preserved when we drop the zero (z-1) from P(z), i.e.. 11^(^)112 = 2^11^(^)112 where Proof: nGM - E [t^) - B y definition: 00 —00 IT 11^(^)112 = ^ / ei« - 1 de with bilinear transform Eq.(2.59), i.e.. (-1) e^'^(l-a) + ( l - a ) j<*> = P—^rh s—h ^€^«(1 + a) - (1 + a) where <T = p(l — a ) / ( l + a), we have: de = . ^ 7 = jccot 2<Tdu Referring to Eq.(2.59), we have: ei" - 1 de C 2(7cL; I = 12;r y —00 00 27r y VJW —00 00 = 2^11^3)112 - 4a2 a2 + a;2 (731) ### The orthogonality is also preserved for base functions Gi{z). Proof: For any given Gm{z) and Gn{z), die inner product is defined as: {Gm{z),G^iz)) = G,n{e^')Gr.{e-^')de .. 1 ( l + «)^ fil-ae^T-" (l-«^-^T"\,, —It 27rj 2;? 7 (ei^ - a)"* (1 - aei^f -•7! 1 ( l + g f / (1 - g - ^ r - " - ^ -dz 27ri 2p / (^-ar-"+^ Assuming m>n, we have: 2Xm - n)! 2 - 0 ef^'"-" ^ =0 When m=n, we have: C 2p 2p(l - 2—a (1 - az) a2) ### To make die base fimction orthonormal and H 2 norm equal, we choose: " " l + 2p ^ """2 So then, the discrete orthononnal Laguerre base function is derived as: ^/ï^ n-azY-"- {z-a)\z-a ) (2.61) li(t) K z-a u(t) y(t) 1-az z-a l-az z-a 1-az z-a 1-az z-a ^ Figure 2.7: Laguerre Filter Netwwk It is interesting to note that Eq.(2.61) has exactly the same fonn as that in the literature [King and Paraskevopoulos, 38 and 37] where it is derived from the time-domain discrete Laguerretype polynomials [Gottlieb, 26] through DTFT (Discrete Time Fourier Transform). The work here establishes an inner-relationship between continuous and discrete Laguerre functions. Because of that, the physical meaning of the time scale a becomes clearer and, more importantly, we define the space where the discrete Laguerre function belongs. Finally, the discrete version of Theorem 2.4 is derived as follows. To be consistent with literature, we denote K = -y/(l — o^). Corollary 2.5: Let the function G(z) be strictly proper, analytic in \z\ > 1, continuous in \z\ > 1, andGkiz) = ^(^l^y~\ then for any \a\ < 1 N J^TkG.iz) Giz)=Um iV—•oo (2.62) ' The network realization of discrete Laguerre expansion is shown in F i g 2.7. We call Laguerre filter gain. Let the true system to be estimated be vit) = G(z)uit) + v{t) 33 (2.63) It follows from Eq.(2.62) that N oo G{z) = J2rkGkiz)+ '-''Gkiz) k=l k=N+l = Go{z) + G^{z) (2.64) Following the same notation as in the previous section, we call Go(z) the nominal model. Widi the orthononnal property of the Laguerre filters, there should be no projection from G „ ( 2 ) to Go{z), We use the predictor model: N m = J2hGkiz)u{t)+v{t) k=l N = J2 ^Mz) + v{t) k=i (2.65) Using standard parameter identification algorithm which minimizes the output difference in a least squares sense, the estimate fk should converge to r^, i.e., Go{z) = Go{z), when data collection hit) = hinT) ; n > N at least (2.66) are driven by white u(t). The corresponding frequency function estimation is N Goie^'^)=J^fkGkie^'') k=i (2.67) The same example as in previous section Eq.(2.15) is used here as a simulation plant We use a discrete Laguerre network with a=0.7. F i g 2.8 is the simulation result with filter number N=3 and Fig 2.9 is witii N=6. Following the notation of the previous section, the uimiodelled imcertainty is as follows: AG{z) = G{z)-Goiz) = [G{z)-Goiz)]+ [Goiz)-Goiz) = A^G{z) + A.Giz) (2.68) With slight difference from literature, we call A „ ( ? the unstructured uncertainty, and AgG the structured uncertainty. Bodi can be eliminated from the view point of estimation. Structured 101 -I r- -I 1 Gain r—r—I—r- solid line: real system 100 I dashed line: estimation / 10-1 10-2 100 _| 1 1 L. _J 101 Phase -2 - 1 I 1_ 102 -1 1 1 1—I—r- real system 100 102 101 1st order system; N=3 Figure 2.8: Frequency Response Estimate by Laguerre Filters, N=3 101 -1 1 1—r- Gain -| 1—1—T- 100 10-1 10-2 100 101 Phase 100 101 Ind order system; N=6 Figure 2.9: Frequency Response Estimate by Laguerre Filters, N=6 102 uncertainty can be decreased by correcdy identifying the Laguerre gains. Unstructured uncertainty can be reduced by increasing the filter number N and take log a close to the dominating time constant of the system to be approximated. For structured uncertainty, with the proper algorithm and growing data collection with white excitation, fk wiU convei;ge to asymptotically. For unstructured uncertainty, we avoid the discussion of how to choose the ' N ' and 'a' (diis has been discussed thoroughly in the paper [Fu and Dimiont, 20]) and concentrate on the question of how to estimate the unstructured uncertainty in terms of frequency domain description because it is vital information for robust control. In this approach, we follow the direction that the unstructured imcertainty can be described by an impulse response which is bounded by a decaying exponential. In the last section, this assumption is used to bound the covariance function in a stochastic description. Now, we use it in a deterministic description as in recent literature [Middleton and Goodwin, 50; Kreissehneier, 40]. The unstructured uncertainty defined in our problem is: (2.69) For a stable transfer function G{z): N Applying die orthonormality, we have c (2.70) where zi,Z2,- " t^n are poles of G(z), for \a\ < 1, Zk-a 1-azk < 1 v k l < i According to Eq.(2.70), we have: \rk\ < KaP fc-i (2.71) where Ka is a constant and /3 = max 1 - azk By observing the right side of Eq.(2.69), we may derive: |A„G'(^)|< J2 k=N+l z — a\z — a ) K ( \ - a z \ k-\ -a\z-a J k=N+l 1-P K z-a i.e., the unstructured uncertainty can be bounded by the magnitude of a first-order model, which is consistent with our assumption. From Eq.(2.71), another conclusion is l i m jr^l = 0. Therefore it is possible to use a finite-term summation to approximate an infinite-term simi of the upper bound for Eq.(2.69): |A„G(z)| < K lim z-a M+N+l V \rk\ k=N+l (2.72) A practical way to evaluate Eq.(2.72) is to ignore the infinitesimal terms \rk\ "^k > M + N + 1 and take the truncated expression of Eq.(2.72) as our estimation of Eq.(2.72). Replacing by f^, we then have the bounding estimation: M A„G(z) z-a (2.73) k=i Therefore, with more filters added to the network, we can estimate the frequency domain bounding which measures the goodness of the frequency function estimate G{z), If there is no special reminder, for all the simulations, we take M=5, and the results shown arc after 100 samples with white excitation. Fig 2.10 is the estimation of frequency response and bounding for the first order system Eq.(2.15) Gain and Bounding \Q-3 I I 1 1 i I—I—I—I—I lO" 1 I i 1 101 I I j i_ 102 Figure 2.12: Estimation of Frequency Response and Bounding for a 2nd Order System with Laguerre Filters, N=6; M=5 with N=3 and M=5, i.e., three Laguerre filters for the frequency response estimation and five Laguerre filters for the uncertainty boimding estimation. F i g 2.11 is with N=6 and M=5. Even though number ' M ' stays the same, the bounding decreases as N increases. Phase plots are not shown since we concentrate on the boimding on magnitude plots. Referring to F i g 2.8 and 2.9, the bias degree of phase plots is no larger than that of gain plots. Also this time we take a second order system: vit) = l.lyit - 1) - 0.72j/(< - 2) + 0.02«(< - 2) (2.74) as an example. Fig 2.12 is the results with N=6 and M=5; Fig 2.13 is the results with N=10 and M=5. It is clear that the lines for real system, estimate and bounding are hardly separable for N=10. Finally we take Rohrs' model as an example for stiff systems. The continuous expression of the model is as following: ^/ ^ 2 / 229 \ Gain and Bounding 10' F , , , 1 I I I I . 100 I 100 2Q'3 I I 1 I I I I I ' , ' I 1 I 101 1 1 I — I 102 Figure 2.13: Estimation of Frequency Response and Bounding for a 2nd Order System with Laguerre Filters, N=10; M=5 With sampling interval = 0.2 sec and Zero-Order-Hold method, the discrete model derived is as follows: . _i> _ 0.1584(1 + 0.9952z-^) ( l -|- 0M49z-^)z-^ ^ ~ (1 - 0 . 8 1 8 7 z - i ) ( l -1- 0.0917^-1 - 0.0025^-2) 0.1584^-^ + 0.1647^-^ + 0.0071^-^ ~ 1 - 0.9104^-ï + 0.0776Z-2 - 0.002^-3 Here, we take the Laguerre time scale a=0.7 which is the same value for all examples in this sectioa The Laguerre filter number N is 5 and M is 5. F i g 2.14 is the estimation result when running time is at t=4 sec. Fig 2.15 is the estimation result when running time is at t=20 sec; Bounding estimation is improved comparing with Fig 2.14. Throughout the section, the inputs to the system are white noise signals. Different input signals will have different effects on the estimation, e.g., a square wave input will force the estimation to have a better fit i n lower frequencies because of the lower frequency components within the square waves. In a very direct sense, the estimated bounds can be used as a measure of the difference between the estimated frequency function and the true but unknown function. Besides its usefulness for robust Gain and Bounding -I 1 r 1—I—I—I solid line: real system (Rohrs-model) dashed line: estimation 100 dotted lines: botmding 10-1 - 10-2 100 -J —1 i_ I I I L. 102 101 Figure 2.14: Estimation of Frequency Response and Bounding for Rohrs' model at t=4 sec Gain and Bounding 101 -I 1 1—I—I—r- -T 1 1 1 1—I—I—r solid line: real system (Rohrs-model) 10° dashed line: estimation dotted lines: bounding 10-1 - 10-2 100 -1 101 1 1 1 L. 102 Figure 2.15: Estimation of Frequency Response and Bounding for Rohrs' model at t=20 sec Figure 2.16: Estimation of Frequency Response of Rohrs' Model for L=10 (N=10; M=0) control, for the estimator, it can be used to judge i f the estimation is good enough. Although there is no guarantee that the derived bound is hard bound unless the number ' M ' in Eq.(2.73) is large enough, simulation shows that, for all the examples in this section, with only M=5, we always capture the modelling error. If we need stronger confidence, we simply use a larger number ' M ' . For further application to robust control. It is very interesting to look into the interplay between the choice of filter nimiber ' N ' and ' M ' . Suppose is a very good estimation of a stable system G(z), we still take the example of Rohrs' model Eq.(2.75). In Fig 2.16, the Bode plot is shown with L=10 (N=10; M=0) and a=0.7. Comparing the true system with the estimated one, it is reasonable to take G{z) as the real system G(z). Then we set off to estimate the uncertainty bounding by sacrificing some accuracy of system transfer fonction estimation. As before, we use ' N ' to represent the number of Laguerre filters used for system estimation and ' M ' the number used for uncertainty bounding estimatioa But in this case we take the relation: N + M = L; 42 L = IQ 101 ^ -I 1—I—I I I I 11 1 1—I—I I I I t 101 = 1 1 1—1 1 1 1 11 1 1 b 100 10° Ë simulation C 10-1 10-2 100 N=6; _i ( I simulation D 10-1 = \ M=4 I 11111 -I I I 1 1111 102 101 10-2 100 N=8; M = 2 1 1 1 1 1 1111 101 1 t 1 1 1 111 102 Figure 2.17: Magnitude Plot of Frequency Response and Uncertainty Bounding for Rohrs' Model Fig 2.17 is the simulation result with different ' N ' and ' M ' . If we define the set ft as follows: G{z) :|G(.)-G(z)|<lA.G(r)|} It is easy to see from Fig 2.17: fir) C ftc C ÎÎB C For any given robust control law which deals with the uncertainty, the closed-loop system designed for ÇIA will have the best robustness but the worst performance; and vice versa, i.e., one designed for ÇÎD has the best performance but the worst robustness. Therefore, instead of using all ' L ' filters for system estimation, the choice of ' M ' or ' N ' can be used as a tuning parameter between system robustness and system performance. This will supply die system designer with one more important control knob. 2.3.3 Estimation of Frequency Function and Uncertainty Bounding Function via Recursive Digital Fourier TVansform In this section, we attack the problem from a different direction. We stiU confine ourselves to the discrete-time single-input and single-output case. Assume that the true system is subject to Eq.(2.1) same as shown in Fig 2.2, using Eq.(2.10)(2.11), then: YNU<^) = GT (e^") UNU^) + WN{J^) + VN{J'^) (2.76) The term W}v(ja>) is due to the effects of finite-length data. VN{JU;) is defined as: 1 VNUU) ^ = -7= E K n ) e - ^ - " (2.77) We can describe the E T F E as: Since we will not always be working with N-point sequences that begin at n=l, we define the following version of Fourier transform for a sequence of N points ending with time index t: UMiMt) = - ^ J2 t*(n)e-^'^" (2.79) YNiMt) = ^ J2 yin)€-^'^ (2.80) A useful recursive equation for computing f^ArOw|f) and YN{j<^\t) is obtained as: UNiMt) = UNiMt - 1) -H [ « ( 0 - u{t - N)e-^^'^] e-^"^ (2.81) YNiMt) = Y^iMt - 1) + [vit) - vit - N)e-^^'^] e'^"^ (2.82) Eq.(2.81)(2.82) can be derived directly from Eq.(2.79)(2.80). Now we define the N points discrete Fourier transform (DFT) of y(t) and u(t) at the N frequency points. Let: ^k = j<^s ; k = l,2,---,N 44 (2.83) where: = 27r/T, is the sampling frequency. We shall assume unit sampling interval for simplicity. Then: t UN{M\t)^^ 1 YN{ju;k\t) = - ^ •V ^ w(n)£"^ (2.84) ' J2 yi^)^"' (2.85) n=t-N+l -N+1 where: E = (2.87) At frequency point w^, we have: é . ( . - | , ) = «.(e>".) + ^ H . ^ where: As the first step we will develop results that can be used to boimd the effects o f using finite-length data to compute the frequency functioa Again, assimie that the system is subject to Eq.(2.2). Let {y(t)} and {u(t)} be related by the stable system Gr(e^'^). y{t) = GT{e^''Ht) (2.90) The input {u(t)} for i < 0 is unknown, but obeys \u{t)\<Cu Vf (2.91) Theorem 2.6 [Ljung, 45]: Let N UN{JO^) = ^Y.<^>~''^' 1 ^ (2.92) Then: where: with: Proof: We have by definition: ^ k=\ t=l = ^'E9Tik)e-^'-J^u{t)e-^' Now: t=iV-A;+l N < V-'^ t=N-k+l ViV hence: 1 2 < °° —C^Cg N-k So the Eq.(2.95) is true. fi n 11 The remaining problem is that this bounding for is too loose to be truly useful. So at first we assume that we know a coarse bounding function on the magnitude o f the impulse response of the true plant (this assumption can be removed later using IDFT), and make the bounding function a litde tighter. For example, assume: |<7r(0l<5oA' ; 0<X<1, go >0 (2.99) that is aU the poles of grit) are in the open unit disk. Inspect Eq.(2.98) and use our notation with index t, we have: oo = Y,9T(n)E^'[UN{M\t - n) - UNiM\t)] (2.100) n=l Under assumption of Eq.(2.99), we find that give some finite integer M , the magnitude of WN{ji^k\t) is bounded for each Wjt as following: Theorem 2.7: Given the system Eq.(2.94), the bounding can be ejq)ressed as M-l \WN{Mt)\ < E \9T{n)\\UN{3<^k\t - n) - UN{Mt)\ n=l + 7 f ^« E "l^^-WI + -^C„ E N\9T{n)\ n=M * (2.101) n=N for M<N, and M-l \WN{3'^k\t)\ < E \9T{n)\\UN{Mt -n)- UNU'^M n=l oo "^Tm^- E " N\9r{ri)\ (2.102) n=A/ for M>N. Proof: Using the triangle inequality on Ek}.(2.1(X)): M-l \WNiM\t)\ < E \9Tin)\\UN{M\t - n) - UNiM\t)\ n=l oo + E n=M \9T{n)\\UNiM\t - n) - UN{M\t)\ (2.103) with: E \9T{n)\\UN{j'^k\t -n)- UNUiOk\t)\ = n=M N-l E \9Tin)\\UNiM\t -n)- UNiM\t)\ n=M oo + E \9Tin)\\UNiJo^k\t - n) - (2.104) UN{M\t)\ we can rearrange the components as following, \UN{j^k\t - n) - UN(M\t)\ m=t-N-n+l ^.iV m=t-iV+l m=t-JV-n+l ^ • ' ^ m=t-n+l n< N (2.105) and if^Ar(M|i-n)-t^Ar(M|t)| n>N (2.106) So the Theorem is proved by substituting Eq.(2.105) and (2.106) back to (2.104). ### Using assumption Eq.(2.99): E ^ M O i < E*5oA* (2.107) t=M we have, when M<N, 7V-1 n=A/ oo (2.109) n=N when I V ^ , the Eq.(2.108) disappears. Eq.(2.109) becomes: E^SoX- =^ ^ (2.110) n=M So we obtain: Af-l WNiM\t) = Y\gTin)\\UNiJoJk\t-n)-UN{joJk\t)\ n=l 2Cu9o{M\^ - ( M - 1)A^+^ - A^+^) (2.111) ^/iV(l-A)2 for M < N , and M-l WN{3<^k\t) = \9T{n)\\UN{Mt-n)-VN{3'^k\t)\ +VN{1 (2.112) - A) for M>N. From Theorem 2.7, we know \WN{j<^k\t)\<WN{jt^k\t) ; (2.113) A;=l,2,.--,iV The bounding function of Eq.(2.111) and (2.112) can be computed on line by using the current N-point DFT of u(t) along with M - l old N-point D F T of u(t). We note that the second line o f the Eq.(2.111) and (2.112) can be made arbitrarily small by choosing M to be sufficiently large. (In the fiiture, it is possible to use IDFT to get grin) estimate and use it as 5r(w) in Eq.(2.111) and (2.112). ) Continuing the above deduction, now we use a direct search method to update the error botmding function and frequency function estimate by judging i f useful information is present. Here we assume the additional noise variance is known, i.e.: \VN{.j<^k\t)\<<^ yt (2.114) We define the corresponding frequency domain error AGsie^'^''\t) from Eq.(2.88): (2.115) the estimated error: AGN{e='^^\t)\ = (\wN{jo;k\t)\ + here \VNiM\t)\)/\UNiM\t)\ (2.116) is given by Eq.(2.111) and (2.112). So: AGN{e^'^-\t)\ < AGN{e^'"'\t) (2.117) By definition, when |y/v(ja;)i.|t)| - 0, diis error corresponds to the mistructured wicertainty. Otherwise, this error is the uncertainty bound for both unstructured uncertainty and the variance of the noise. In the following we will discuss a technique or algorithm to combine the frequency function estimates and the corresponding uncertainty bounding fimction firom different time intervals, i.e. how to combine all of the past frequency domain information into a recursive way to estimate the frequency fimction and the uncertainty boimding function. The basic idea is that at a given frequency point Wfc we use the value of GN i^^"''\t) which has the smallest corresponding uncertainty bounding function A G i v (c'''^'lO . We define the smallest bounding function at w^: AG(e'''*'*|t) = m i n ( AGjv(e^'^* } (2.118) and the corresponding frequency function estimate at ojk' G(e^-|f) = {GAr(e^-|m) : |^ô(e>"Mm)h|AG(ei-W|} (2.119) We define, for time index t: AG{e^'^''\t) = Gie^'^-lt) - Gr(e^''^*) (2.120) Then, Eq.(2.117) ensures that at time index t: AG{e^'^''\t) \ < AG{e^'^''\t) (2.121) For use in an on-line implementation, the following simple recursive algorithm will be used to compute (7(e-"^*|^) and A(?(e-"^'|f) at a given frequency u>k', if: AÔNie^'^'lt) < AG(e''"*|<-1) (2.122) Figure 2.18: Estimation of Frequency Response and Bounding via RDFT tlien: AGie^'^"\t) = AGN (e^'^'\t) G{e^'"'\t) =GNie''''\i) (2.123) (2.124) else set: AG{e^'^'\t) = A(?(e''^*\t - l) (2.125) (?(e'"'^»\t) = G{e^'''\t - 1) (2.126) Thus, this algorithm can guarantee that the frequency fimction estimate and the corresponding error bounding function are updated only when useful information is learned at a given frequency. The estimation is shown in F i g 2.18. Comparing with Fig 2.1, we can see the improvement for gain and phase plot estimation. But the bounding estimation is very rough because there is still no connection between each frequency. For Fig 2.18, the simulation model is Eq.(2.15). Since we are working with real valued time domain signals, the properties of the DFTs of real valued signals can be used to show that: G(e"^*\t) = G(e''^''-* (2.127) (2.128) for {A; = 1,2, • • •, iV/2}. where denotes complex conjugate and where we have assumed that N is even. This means that the information for the frequency points [k = 1,2,-••,N] is contained in the information for the frequency points {k = 1,2,• ••,N/2}. The remaining problem is that the additive disturbance v(t) is not known in usual cases. Some method must be developed. The difficulty is that there is no efficient way to estimate V(ja;) i n the frequency domain. A practical hypothesis for improving die E T F E with unknovm V(jw) is suggested by Ljung [Ljung, 45], i.e., die disturbance spectixim does not change very much between frequency intervals or in the "width" of the weighting function. Then a weighted average may cancel the effect of V(ja;). So far, it remains die only way to improve the frequency function estimation and is used by all the researchers who deal with frequency domain estimation problems. In the next section, we will develop an algorithm to eventually estimate V(ju;) in the frequency domain. In doing so, we are able to significantly improve the frequency function estimation. 2.4 External Disturbance Spectrum Estimation via Frequency Domain Recursive Least-Squares Algorithm 2.4.1 Frequency Domain Recursive Least-Squares Algorithm The initial problem is how to deal with the potential problem caused by Eq.(2.31) and Eq.(2.89), i.e., ^^(^"''" = TO '^-'^^ if the estimated denominator U{j(jJk\t) happens to be zero at a particular time for a particidar frequency. This leads us to develop a very straightforward frequency domain recursive algorithm which uses the same procedure as its time domain equivalent. First, let us recall the recursive algorithm for real-time parameter estimatioa A general formula is as follows: ê{t) = 9{t - 1) -1- fi{i)i^{t)e{t) (2.130) where 0 is the parameter vector; rj} is the regression vector, e is die prediction error, fi is the adaptation gain. As we know, i f the adaptation rate is much faster tiian the parameter changing rate, the estimation based on Eq.(2.130) is convergent and consistent under the condition that the noise is white [Âstrôm and >\^ttenmark, 5]. Even for the same formula Eq.(2.130), different approach of getting ^(t) has the different conveigent and consistent characteristics. For instance, the adaptation gain ^(t) can be simply chosen as a constant matrix based on a gradient method approach. The choice of the gain is made to balance tracking of time variations enhanced by large /x, versus smoothing of noise, achieved by small fi. In the later case, the initial condition becomes a very crucial factor to the convergence speed. To avoid guessing the initial condition and choosing fi, a more appropriate approach is the R L S (recursive least-squares) algorithm where fi(t) is a gain matrix P(t) which is automatically updated and is closely related to the prediction error covariance matrix and so is the optimal choice in terms of minimum variance and convergence rate. Therefore, we choose to adopt the R L S algorithm in our context and use the version of Eq.(2.130) based on a R L S updating rule later o a Now it is time to come back to our original zero denominator problem about Eq.(2.129) where Y{ju:>k\t) and U(jwk\t)is obtained through either Kalman filtering or RDFT algorithms as derived in previous sections. To prevent Eq.(2.129) from dividing by zero, we update Eq.(2.129) in a recursive manner and rewrite it as YiM\t) = Gie^'^'-\t)UiM\t) (2.131) From now on, the subscript " N " wiU be ignored since it only indicates that our D F T is N-point DFT. A regression model in time domain is yit) = e'^imit) (2.132) Comparing Eq.(2.131) with Eq.(2.132) and referring to the time domain R L S algorithm Eq.(2.130), it is not difficult to derive the following frequency domain version of R L S : G(e^'^»lO = ^ ( e ^ - M ^ - 1) + pUo^kmUoJkmMt) (2.133) where = Y{juk\t) - G'(e^'^*\t - \)U{jUk\t) ti\M\t-l)U\3'^k\t) with forgetting factor A, the adaptation gain updating equation becomes: A + UU(^k\t)fi{M\t - i)UiM\t)j At this point, it becomes clear that we not only solve the problem of zero denominator, but also derive a much better way to estimate frequency function because we arc able to make use of the prediction error e{jwk/t). 2.4.2 External Disturbance Spectrum Estimation It becomes possible that we could actually estimate the disturbance spectrum when the recursive updating equation Eq.(2.133) was derived. After that point, we simply expand the idea a littie further by considering the overall system Fig 2.2 which is described in frequency domain by equation Eq.(2.88), tiiat is: Y{M\t) = GT{e^'''')U{M\t) + W ^ C M I O + V{M\t) (2.134) As defined eariier, the term W{jwk\t) is due to the effect of using finite-length data. A FIR filter of order N can be completely described by giving the value of its frequency function at N points equally spaced around die luiit circle [Bitmead and Anderson, 8]. In other words, we can completely describe a system function by its values at N frequency points i f the system impulse response lengdi is less than N . Furthermore, we have already derived the algoridim to eliminate and bound the term W{jwk\t) which is only related to the input spectrum Eq.(2.100). Therefore, it is reasonable to put the first two terms in Eq.(2.134) togedier and leave only one term V(jwk\t) — the additional noise or external disturbance by definition. As a result, our discussion in this section is based on the following equation: YiMlt) = GNie^'^-MMlt) + V{M\t) (2.135) where So far, in the literature, there still is no effective way to estimate an external noise fimction V(jw) directly from a harmonic analysis method. The most practical way concerning disturbance spectrum estimation is under the condition that the frequency function estimation exists already [Ljung, 45]. From our simulation results which will be shown later, we wiU see that it is impossible to have the correct frequency function estimation without the knowledge of noise function at the same time because V(jw) will not only affect the variance of die frequency function estimation but, more importantly, will cause large bias. In terms of frequency fimction estimation by means of Fourier analysis, a weU-known statement concerning noise is that the variance of the frequency function estimate at a certain frequency is asymptotically given by the noise-to-signal spectnun ratio at that frequency [Ljung, 45]. This statement is based on the assiunption that VQu) is a zero mean variable since v(t) is assumed to have zero mean value. Following the same line in the literature, the weighting fimction, which is often called the frequency window in spectral analysis, is used to trade off the bias and variance of the frequency function estimate in order to obtain a smoothed estimation. Actually, diis windowing technique only redistributes die bias and variance at one frequency to odier frequencies i f the estimation itself is a biased estimatioa It wiU be apparendy disabled i f the noise signal has wide-band frequency components like white noise because the estimation is biased at all the frequencies as shown in Fig 2.22. Clearly, for a frequency domain method, the term V(jw) can not be simply ignored from Eq.(2.135) as it is in time domain methods. Using frequency domain R L S algorithm, it is possible to estimate the noise function together with frequency function. Here we treat the frequency domain noise function V(jw) as a deterministic signal. A n unrealistic assumption here about the noise function V(ju)) is that it does not change as time changes (this condition is man-madely true for the simulations later) or that the adaptation rate is much faster than the power spectrum rate of change. Therefore, Eq.(2.135) can be denoted as: (2.136) where And Eq.(2.133) becomes: èiMt) (2.137) = è{jojk\t - 1 ) + iJiU<^k\t)i^{jojk\t)eU'^k\t) where e ( M | 0 = Y{jwk\t) - ë^{juk\t - i)V'(Mlf) fi{ji^k\t - iMJ^k\t)i^'^{M\t)KM\t 1 + i^^{jo^k\tMM\t - -1) imM\t) with forgetting factor A, the adaptation gain updating equation becomes: KM\t -1) - ix{M\t - iWM\i)i'^iM\tMM\i -1) In the simulation, we take the same model in previous section Eq.(2.15) as a example. In order to demonstrate that the algorithm can be used in a very general framework, we take two extreme cases into consideratioa In the first case, the noise function is a single frequency sinusoidal function. 101 100 w=2*PI*x/256 phase RDFT with external noise G^/S=10) Figure 2.19: Estimation by Eq.(2.133) with Single u Noise 102 -1 r- _j L- ;am -T 1—I—I—n 101 10» 10-1 10-21 100 102 101 w=2*PI*x/256 -T r- -I phase 1—I—r-T—I -2 100 -J I I I I (_ 101 -I I RDFT with external noise G^/S=10) Figure 2.20: Estimation by Eq.(2.137) with Single w Noise I I I I L_ disturbance spectrum estimation 1 1 1 r 1— T i l l •T 1 I I I I I I• real spectnim / t M - f 1 1 1 1 \ i estimation 1 1 1 1 1 1 1 n 1 1 I 100 1 1 1 1 1 1 I1 1—1 1— 11 I1 I1—I— — 101 1 1 1 1 1 1 1 1 1 1—1— w=2*PI*x/256 Figure 2.21: Noise Spectrum Estimation for Single w Noise In the second case, the noise is white which means that noise fimction is a mixture of sinusoidal functions from all the frequencies. We apply the same algorithm to both cases. For the first case, without noise function estimation, thefi^equencyfunction estimation with single fiiequency distmbance is shown in Fig 2.19. It should be mentioned here that the noise to signal ratio is ten which means that the input signal is almost buried by the noise. We can see that the estimation at that frequency has a large bias. After we put the noise function into the prediction function Eq.(2.136), a great improvement is achieved, the frequency function estimation is shown in Fig 2.20 and the disturbance spectrum estimation is shown in Fig 2.21. For the second case, without noise fimaion estimation, the frequency function estimation Fig 2.22 becomes totally erroneous because the bias happens at every frequency point; the whole curve is pushed up by the noise energy. After considering the noise term, the frequency function estimation Fig 2.23 is improved as for die case of a single frequency disturbance. The explanation is quite simple because the algorithm takes into account all N frequency points. The noise spectrum estimation for white noise case is shown in Fig 2.24. -41 100 • ' ' ' — ' ' ' ' ' — — ' ' 101 ' — • ' ' ' ' 102 RDFT with external noise (N/S=10) Figure 2.22: Estimation by Eq.(2.133) with White noise 100 101 RDFT with external noise Qi/S-lO) Figure 2.23: Estimation by Eq.(2.137) with White Noise 102 disturbance sp>ectnim estimation 100 101 102 w=2*PI*x/256 Figure 2.24: Noise Spectrum Estimation for White Noise Qearly, the noise function estimation helps to improve the frequency function estimation, which is similar to the case in time domain that noise model estimated by Generalized R L S algoridun helps to improve the estimation of system model. For the time domain method based on A R M A X model, a whitening filter— noise model has to be specified before hand i f the noise is not white. It is a much more difficult problem for a time domain method, for instance, to specify the stmcture of noise model i f the disturbance energy is only concentrated in a single frequency which is the first case in our simulation mentioned above, i.e., a very narrow band-pass filter has to be constructed within A R M A X model. With our frequency domain method, we can deal with a more general situation, that is, we do not have to know whether the noise energy exists at one frequency or all the frequencies as shown in Fig 2.20 and Fig 2.23. 2.5 Conclusions In this chapter, we have summarized and developed some methods for frequency domain estimation. Besides the frequency fimction estimation and disturbance spectrum estimation, the focus has been on the uncertainty bounding function estimation, which is the key issue in an adaptive robust framework. This chapter started with the E T F E . Then, the DFT block data transform was changed to a "recursive" transform by means of Kalman filtering. The advantage of using Kalman filter is its real time capability and noise smoothing feature. The Kalman filter also provides a way to deal with the finite-length data problem in block D F T since it is able to memorize the past information. Later, as the main results, two methods are developed to estimate uncertainty bounding fimction. The Laguerre orthononnal function approach gives a much tighter uncertainty bounding, e.g.. F i g 2.10, but widi the assumption that the uncertainty is bounded by Eq.(2.72). The RDFT approach gives a sufficient boimding. Fig 2.18, but it is still too loose to have any practical use. In the end of this chapter, we developed a method to estimate disturbance spectrum which is closely related to least squares algorithm. In doing so, the frequency function estimation is significandy improved. It is remarked here that so far, for a stable, linear, lumped, and time-invariant system, we are able to estimate the frequency function, uncertainty bounding function and disturbance spectrum function. Chapter 3 Quantitative Adaptive Robust Control in Frequency Domain 3.1 Introduction In this chapter, we proceed further to adopt a robust design strategy which could make use of the information provided by the frequency domain estimator derived in the previous chapter. To distinguish from other adaptive robust control systems, we call our method Quantitative Adaptive Robust Control (QARC) since we take the estimated uncertainty bounding into consideration through the design method originally from Q F T (quantitative feedback theory) [Horowitz, 32], and also the whole scheme adopts ideas from both adaptive and robust viewpoints. A very close interpretation is to think of S O A C (self-osciUating adaptive control) systems [Astrom and Witteimiark, 5] with the improvement diat multiple frequency points are considered and stability margins are quantified via uncertainty bounding estimation. From an overall perspective point, both currently studied adaptive systems, robust systems and so-called adaptive robust systems have dieir advantages and disadvantages. It is well known that adaptive control can effectively eliminate the stmctured uncertainty caused by parameter variations in the nominal model, but has difficulties handling unstructured uncertainties which are beyond nominal models. Aldiough the objective of adaptive control is to deal with unknown process and environment, a substantial amount of a priori knowledge is needed to guarantee stability. System order, delay time and noise model have to be known, and restrictions such as SPR (stricdy positive real) condition, slow adaptation rate and richness of driving signal have to be imposed. On the other hand, robust control, with tolerance margins covering imstructured uncertainties, can guarantee the stability of the closed loop system. But because it assumes hard bounds on the uncertainty, it ends up with a very conservative design. The direct combination of adaptive control and robust control, so-caUed adaptive robust control or robust adaptive control [Cusumano and Poolla, 14; Kosut, 39], adds robustness to adaptive control. But it does not solve the problem completely since there is stiU no on-line access to unstmctured uncertainties. Therefore, we may still end up with eiUier a very conservative design or an unstable closed-loop system. -{^ parameters Design Regulator y Estimation < U Process Figure 3.25: Block Diagram of a Self-Tuning Regulator At this point, with very little a priori knowledge about the system and the envirorunent, the information provided by our frequency domain estimator makes it possible to implement the kind of adaptive robust controller where the performance-robustness trade-off is made in a quantitative way. Therefore, it is possible to obtain the best closed-loop performance out of a stable closed-loop system. 3,2 Stability of Adaptive Control Most adaptive controllers are model-based control strategies. The discussion below is based on the STR (self-tuning regulator) scheme since M R A C is similar to the direct STR, The block diagram of the whole system is shown in F i g 3.25. In the diagram, the process is described by the A R M A X model: y(t) = A{z-'')y{t) + B{z-'')u{t -d) + e{t) (3.138) and the predictor model for the estimation part is: y{t) = Â{^~')yit) + s{z-'Ht-d) (3.139) The adaptive control system can be thought of as consisting of two loops: an irmer loop, which is an ordinary feedback loop composed of the process and the regulator, and an outer feedback loop H(z-) -1 z e(t+i) t ) - X £\|r(t) Figure 3.26: Error Model Drawn as a Feedback Loop that makes the adaptive system non-linear. The parameters of the regulator are adjusted by the outer loop in such a way that the error between the process output y and the predictor model output y becomes small. So the outer loop is also a regulator loop. The key problem is to determine the adjustment mechanism so that a stable system, which brings the error to zero, is obtained. Stability is a basic requirement in a control system. For adaptive control, analysis is very complicated because of the existence of the outer loop. Much effort has been focused on the issue of outer loop stability, i.e. the parameter updating loop. A n effective way to piu-sue this is to look at the evolution of the adaptive error system [Anderson, Bitmead, Johnson, Kokotovic, Kosut, Mareels, Praly and Rielde, 2] or commoidy the so called error model of the following general form: ait) = -erj;{t)n{s){i;'^{t)9{t)} (3.140) In discrete time, error model is described by an ordinary difference equation: e{t + 1 ) = eit) - ,i>{t)H{z-'){rp'^it)e(t)} (3.141) This equation is depicted as a feedback loop in Fig 3.26. The non-linearity appears as the product tp{t)0(t) in Eq.(3.140) and Eq.(3.141) since, in a closed loop configuration, the values of regression vector ip(t) are affected by the parameter vector 6{t). As we know, adaptive control is based on a simplified system model and noise model. With the existence of the unstructured uncertainties, i.e., model mismatch as usually noted, the behavior of the closed loop system may depend drastically on the command signal and the disturbances. The solution obtained for one command signal may be stable; another command signal may give an unstable solution. A famous example is the Rohrs' model [Rohrs, 63]. In another word, different process models are obtained for different driving signals. Also, although a solution is stable, a perturbed motion may diverge since the parameter estimation is not convergent with unknown disturbance. Therefore, non-linear theories have to be introduced into adaptive control. One of them is averaging theory which can apply directly to the error model Eq.(3.140) and (3.141). Averaging theory permits the analysis of equilibria and local properties around equilibria. The basic idea is that the parameters change much more slowly than the other variables of the system. Under the assumption that the parameters 0{t) are constant and by averaging Eq.(3.140) or (3.141) over a period of time, rp{t) become a time invariant vector ^ . Then Eq.(3.140) or (3.141) as weU as Fig 3.26 becomes a linear time invariant error system e{t +1) = e{t) - etpH{z-'^) {rp'^eit)} (3.142) which can be stabilized by imposing conditions on operator fiinction H{z~^). Fmally, the stability of Eq.(3.142) implies tiie stability of Eq.(3.141) for all 0 < e < Co [Guckenheimer and Holmes, 28]. Qearly, in order to apply averaging theory, it is necessary to use small adaptation gain c. Only by choosing e sufficientiy small, the rate of change of the parameter 9 can be made arbitrarily small. Along the same line, the small gain theorem [Anderson, Bitmead, Johnson, Kokotovic, Kosut, Mareels, Praly and Riedle, 2] gives the key result that the closed-loop system is stable i f E{z~^) is SPR and the unstructured uncertainty is passive. Hopefully, averaging theory may lead to a unification of analysis of adaptive systems. There are unfortunately no good methods to determine analytically how small the gain CQ should be. Another difficulty is that the SPR condition is not directiy interprétable in process model or control model since H[z~^) only arises in the error system Fig 3.26. The formation of H(z~^) depends on the controller design strategies and parameter adaptation algorithms [Ortega, Praly and Landau, 56]. A transfer function H{z~^),m error system Eq.(3.141) based on a standard pole-zero placement adaptive controller is introduced in [Fu, 19]. Furthermore, the lack of a parameter equilibrium will disable the averaging method. With all the troubles associated with parameter model estimation, one advantage of frequency domain estimation is that we do not have a model validation problem. Also it is not necessary to treat the frequency domain estimation loop as a closed-loop system since the estimation loop does not produce an output estimate y. Even with the R L S algorithm in Section 2.4.2, the vectors V'CJW, i) and 6{joj,t) are different from their time domain counterpart — regression and parameter vectors. Moreover, the R L S algorithm is only used as a refining process in the estimation. To simplify the problem, we wiU limit the discussion of Q A R C to the territory of linear theory. Therefore, the robust stability in the next section is considered as the overall stability of the complete Q A R C system. 3,3 Stability of Quantitative Adaptive Robust Control Limited by the scope of this thesis, the discussion of die Q A R C stability is in the domain of linear theory. In the adaptive control, the outer loop stability or the stability of adaptive error system is discussed assuming that the inner loop is stabilized. Q A R C here is merely focused on the stability of the inner loop and takes into account the bounding function of the outer estimation loop. Besides the reasons mentioned in the previous section, it is oiu" intention to avoid the discussion of the non-linear error model in frequency domain since it is still not clear how to derive one and how strong a non-linearity exists in Q A R C . So far, widi the persistent excitation and plant linearity conditions, the frequency domain estimation is proved to be unbiased. Instead of analyzing the "error model", we assume it is good enough to take the "error bounding" into consideration for closed-loop stability. Therefore the stability of die Q A R C becomes the robust stability of the inner loop, that is the regulator feedback loop as in Fig 3.27. Analysis of feedback amplifiers by Nyquist [Nyquist, 53] and Bode [Bode, 10] is one of the cornerstones of feedback theory. Further refinement of the theory leads to design methods that explicitiy take process uncertainty into account. For this purpose, we simply employ the Nyquist stability criterion as the starting point. The Nyquist plot of the discrete system of Fig 3.27 is drawn in Fig 3.28. Stability analysis based on Nyquist plot is a very mature mediod. It also gives a direct —• K(z) u , G(z) Figure 3.27: Inner Loop Block Diagram I m Figure 3.28: Nyquist Plot description of the loop sensitivity. The stability margin a is a distance between the curve crossover point and the -1 point. The stability is defined as: (3.143) where a is used to take care of uncertainty. For SISO systems like that of Fig 3.27, the appropriate notion of smallness for the sensitivity fijnction requires the complex scalar [ l + K{e^'^)G[e^'^)Y^ to have small magnitude, or conversely that 1 + K{i^^)G{e^'^) have large magnitude, for all real finequencies w where the commands, disturbance and plant uncertainty AG{e^'^) are significant. In fact, the sensitivity objectives of SISO feedbaclc systems are commonly stated in the same inequalities of the fonn as E:q.(3.143): a{w)< | l + ii:(e^'^)G(e^'") V W < Wa (3.144) where a(u>) is a positive fimction and tOa specifies the active frequeiKy range. Actually Eq.(3.144) is fundamental for the optimal sensitivity theory such as QFT and I M C (Internal Model Control) [Morari, 51]. The basic idea can be readily extended to M I M O problems through tiie use of matrix norms. Selecting the spectral norms as a measure of matrix size 11 + K (e-''^) G (e''^) j , for example: = v^A„.n{[l + K G r [ l + KG]} (3.145) )u(w) is called minimum singular values The corresponding feedback requirements become: (3.146) Going fiirther along this line wiU lead to the synthesis method [Doyle and Stein, 16]. Since our approach is limited to SISO systems, we take Eq.(3.144) as om" initial objective. From our estimation of the unstructured uncertainty, open loop transfer functions will appear as a Nyquist band of Fig 3.29 instead of a Nyquist curve. In terms of frequency domain estimation, we actually have the firequency function estimate at each frequency: (3.147) and the unstructured uncertainty bounding estimate at each frequency: AG{w) = | A G ( e ' ' " ) Gu{e^'')-GLie^'^) which is usually called additive uncertainty description since it appears as G{e''') e ( ? ( e ^ ' ^ ) ± A G ( e ' ' " ) Figure 3.29: Nyquist Band If we define: G{^'^) -v\AG{€P'^)^LG{eP'^) (3.149) where Gu is the tightest hard upper bound (discussion of hard and soft bounds can be found in [Ljung, Wahlberg and Hjalmarsson, 45; Wahlberg and Ljung, 67]), we have the following theorem: Theorem 3.1 (robust stability dieorem for Q A R Q : Assuming that all plants G{e^'^) G H2 in the family, set or band: J7 = |G(eJ-)-G(eJ")|<AG( (3.150) Then the closed-loop system is robustly stable with the controller K (e^'^ ) G H2 ifand only y the inverse sensitivity function for the plant Gu {'^'^) defined by Eq.(3.149) satisfies the following conditions: K{e^'^)Gu{e^'^) not encircling the (-1,0) point and l + K{e^'^)Gu{^'^)\>Q (3.151) Proof: The idea is very straightforward. For the plant family defined by Eq.(3.150), we will have the biggest contour K(e^'^)Gu (c-''^) (see Fig 3.29) where controller K{e^'^) design is assumed available. Then implied by Nyquist Criterion, the stability of K [e^'^)Gu{e^'') will guarantee the stability of the whole family |('(e-'''')||G(ej«)_G(ei")|<AG(u<)}* Therefore, we have proven the theorem. U It II Let us clarify what it means that Theorem 3.1 is not only sufficient for robust stability but also necessary. If condition Eq.(3.151) is violated, then for the set ft defined by Eq.(3.150), there must exist a plant G(e'"^) G ft (flunk of Gt/(e^'^)) for which the closed loop system with the controller K[e^'^) is imstable. However, i f the set ft is a hard bound obtained by approximating the original uncertainty regions, and therefore contains plants not present in the original imcertainty set, then condition Eq.(3.151) is generally only sufficient for the original uncertainty set. Nevertheless, the necessity implies that at least for the set ft defined by E;q.(3.150), with an equal sign i n Eq.(3.151), this is the tightest robust stability condition that can be derived. Theorem 3.1 can be stated in a more general way by defining G{ei'^) - G{e^'^) < A^G(u;) G(eJ'^) (3.152) Where AmG{w) is referred to as a multiplicative uncertainty description and is related to additive uncertainty description as follows: (3.153) G(ei'^) Then the radius, that is the half of width of the band at each frequency in F i g 3.29 can be described as K(e-''*')(?(&''*') AmG{ij). The inequality of Eq.(3.144) in this case becomes: K{e^'^)G{e^'^) A „ G ( w ) < 1 + K{e^'^)G{e^'') or K{eJ'^)G{ei'^) A^G(u;)<l 1 + K{ei'^)G{e^'^) Vw (3.154) With the notation of the complementary sensitivity function for the estimated frequency function G{e^'^) as: (3.155) 1 + K{e^^)G{e:i'^) we have a general version of Theorem 3.1: Corollary 3.2: Assuming that all plants G(e^'^) in the family, set or band fir, fim = { G{e^-) : < AmGiu;) G{e}^) (3.156) have the same number of clockwise encirclements of the point (-1,0) by their Nyquist contour (think of the same number of RHP poles in a model), and that a particular controller K{^^) stabilizes the estimated plant G[e^'^). Then the closed-loop system is robust stable with the controller .K'(e-''^) // and only if the complementary sensitivity function 'r(w) for the estimated plant G{e^'^) satisfies the following bound: T:(u)AmG{u) < 1 (3.157) where: t(u;)A^G(w) ^ sup { t ( w ) A ^ G ( a ; ) } According to Corollary 3.2, robust stability imposes a bound on the oo-norm of the complementary sensitivity fimction T(w) weighted by the midtiplicative uncertainty boimding AmG{u)). 3.4 Robust Performance of Quantitative Adaptive Robust Control 3.4.1 Quantitative Feedback Theory A n important reason for using feedback in the design of control systems is to reduce the sensitivity of the system response to the uncertainties about the system. The fundamental principles of such a design have been studied in several papers [Horowitz, 31 and 32]. A frequency domain method, the so-called QFT, has been developed for designing a linear, time invariant, SISO feedback system, in the simplified configuration of Fig 3.30, to achieve a prescribed degree of insensitivity of output to Figure 3,30: Block Diagram of 2-Degree-of Freedom Control 3 4 5 Time (sees) Figure 3.31: Bounds on Response in Time Domain the uncertainties AG(e-''^) about the process. With a pre-filter F{e^'^) in the configuration. Fig 3,30 is sometimes called Two-Degree-of Freedom Control since r and v have different effects on e: e = •.r — (3,158) Here / i (e-''^) can be selected for stability and good disturbance rejection, and then the pre-filter i^(e"^) can be chosen independently for good set-point tracking. The most direct form of the specifications on the output y(t) would be in the time domain, such as bounds on the step response shown in Fig 3.31. For the QFT design, however, it is assumed that such time domain sensitivity requirements have been somehow, by trial and error, transferred into 'corresponding' frequency domain specifications [Horowitz, 31] in the form of the bounds on T = KG 1 + KG (3.159) This translation from time domain to frequency domain requirements is one of the weaknesses of QFT. These frequency domain bounds are then used as the starting point of their design. Q F T is a controller design method that meets the closed-loop tolerance: (3.160) with the process constraints: e (3.161) fia where 0.^^ is represented as a rectangular shaped template in Nichols Chart Notice that, Eq.(3.161) gives both gain and phase uncertainty, since QFT originally is used to deal with the parameter 6 variations in the process G{€^'^,6), i.e. structured uncertainty by definition. Such numerical specifications of uncertainty templates Eq.(3.161) and the closed-loop tolerance Eq.(3.160) are called "quantitative synthesis". The actual design is performed by means of a Nichols chart. It has been shown [Horowitz and Sidi, 33] diat, with any Go(e''^) G fi^ chosen as the nominal plant, the above quantitative synthesis leads to bounds b(a;) on a nominal loop transfer function i(e-"^) = K{e^'^)Go{€^'^) for each w e (0,oo). For example, fi^^,, the entire uncertainty range of Goie^"^'), is shown as template A B C D in F i g 3.32; and the closed-loop tolerance specification at a; = Wo is derived as Aln|r(e^'^«)| = A l n 1 + i(eJ'^°) < 5dB Since InL = lnK + In G, the pattern outiined by A B C D may be translated (but not rotated) on the Nichols Chart, the amount of translation being given by the value of /ir(e"^°). For example, a trial design of i(e-''^°) corresponding to the template of Goie^^") is at A"B"C"D", where the extreme 40 0 0;25 30 o;5 • .. 20 •ri 10 ?3 6 0 12 -10 D -20 20 -30 -40 -1 -350° -300° -250° -200° -150° L. -100° -50^ ! : 40 0" degrees Figure 3.32: Bounds on Loop Transfer Functions values of liir(eJ'^») are C " ( - l d B ) and A"(-6dB). If for condition A " is chosen to be (-3dB, arg - 7 0 ° ) , then it is guaranteed that A In T{e^'^') \ < hdB over the entire uncertainty range of the plant Also, it happens that, i f argX^(e-''^°) = —70°, then -3dB is the smallest magnitude of LA{e^'^') which satisfies the 5dB specification for A In r(e''^<')|. Any larger magnitude is satisfactory, but represents, of course, overdesign at that frequency. The manipulation of the w = template may be repeated along a new vertical line, and a corresponding new minimum of LA [e^^") found. Sufficient points are obtained in this maimer to permit drawing a continuous curve of the lower bound on LA {e^'^") as b(u>o) shown in Fig 3.32. The entire process may be repeated at other frequencies. In the high frequency range. LA (e^'^) \w>h has to be outside b(w>h). In fact, b(c<;>h) guarantees a certain minimum amount of damping in the response to a disturbance v(t) since the usual constraint on the damping factor can therefore be transferred into a constraint on the peak value of In |r(e-''^) I [Horowitz and Sidi, 33]. Suppose this happens to be 6dB in the present example as shown in Fig 3.32. Therefore from QFT design, we actually derive an infinite number of loop transfer functions which meet both specifications of Eq.(3.161) and Eq.(3.160). A unique Xx(e''^) is obtained by imposing extra conditions. For instance, let i^(e-"^) be assigned a fixed excess number "n" of poles over zeros, so that l i m LA (e''^) = kol{e^'^ ~ ! ) " • t^^" ^ realistic optimization criterion is the minimization of kg. Such an optimum LA (e-''^) is unique and the rather lengthy design procedure is presented in [Horowitz, 31; Gera and Horowitz, 22]. So far the QFT design problem can be only solved graphically; no analytical method is available; and the final residt is in niunerical form which is adequate for research purpose, but is a shortcoming for practical design, in as much as rational function approximation remains to be done. 3.4.2 Robust Performance Criteria for Quantitative Adaptive Robust Control In terms of loop transfer function shaping, Q A R C follows die same procedure as QFT. The improvement here is that we try to simplify the procedure and to find the bounds bO*) analytically instead of moving templates around graphically in the Nichols Chart. This allows the quantitative feedback synthesis to be realized in real-time. First, we begin with using a simple and rational method of converting the time domain sensitivity specifications of F i g 3.31 into frequency domain criteria. These general criteria together with disturbance rejection consideration are fiirdier reduced to simple conditions on the sensitivity function and complementary sensitivity function. Let us suppose that the response y(t) to an input r(t) is required to remain within the bounds: a{t) < y{t) < Kt) i > 0 (3.162) as shown in Fig 3.31, for all processes G[e^'^) in a given closed and bounded set ft as Eq.(3.150). Defining: and The œndition Eq.(3.162) then can be stated as [y{t) - m{t)f < Am2(f) V € ft (3.165) When one considers the fact that the value of a function x(t) at a single instant t is determined, in principle, by the values of X (e-''^) for aU w, the difficulty of finding an exact equivalent to Eq.(3.165) in die frequency domain can be understood. As an alternative to Eq.(3.165), consider the slighfly weaker condition: p p - f [y{t) - m{t)fdt < - f Am^{t)dt PJ PJ 0 (3.166) 0 This condition can also be described in discrete time: I J2 P t=i bit) - mit)f < ^ ^ ' " ' ( O P t=i (3-167) It is to be noted that above conditions apply for every p > 0, that is, Eq.(3.166) or (3.167) ensures that Eq.(3.165) is satisfied on the average over every finite mterval {0,p). In view of this property, it is argued diat Eq.(3.166) or (3.167) is a reasonable altemative to Eq.(3.165) as a robust performance design specification. Referring to Theorem 2.3 or Parseval's theorem, we can derive the following theorem for continuous-time systems. Theorem 3.3: Let <B(s) = £ { / ? ( t ) } , A{s) = £ { ^ ( t ) } and H{s) = ^{s)/A{s) (3.168) H{s) represent a causal system, continuous in s. A necessary and sufficient condition for p p J 0^{t)dt < j 6^{t)dt 0 V 0 < i) < 00 0 is that for all w where the equality sign may hold for some but not for all w. (3.169) Proof: Consider first the case p = oo. Applying Parseval's theorem to Eq.(3.169) and using Eq.(3.168), gives oo y A ( j a ; ) ( l - ^(ju;)/r*(ia;))A*(iw)(L; > 0 (3.171) —oo Clearly, Eq.(3.170) is sufficient for Eq.(3.171). It is also necessary because i f inequality Eq.(3.170) is reversed at some u?o. then by continuity it is reversed in some neighborhood of and a ^(<) exists which will reverse the inequality of Eq.(3.171) and Eq.(3.169). To proof die case p < oo, let ^p(f) denote the function i^{t) = 6(t) yt<p ëj,{t) = 0 Vi > Let Pp{t) be related to Sp{t) by ^p{s) = His)Apis) (3.172) where ^p{s) = &{Pp{t)} and Ap{s) = £>{Sp{t)}. By die causality of H{s), Pp{t) = m yt<p (3.173) Now, oo p Js^t)dt = jêl(t)dt 0 0 Using Eq.(3.172) and we already proved for p = oo, we have oo oo J6l{t)dt> 0 j0'p{t)dt 0 and by Eq.(3.173), oo p oo j Pl{t)dt = If3'it)dt+ j (3l{t)dt 0 0 p p (3.174) which is Eq.(3.169). Again, Eq.(3.170) is not only sufficient for Eq.(3.169) but is also necessary. If the inequality Eq.(3.170) is reversed, a 6p{t) exists that will reverse Eq.(3.174) sufficiently to override E:q.(3.175) and dius reverse Eq.(3.169). ### Since ^(t) is arbitrary, it can be identified with A m ( i ) . According to Theorem 3.3, an equivalent condition to Eq.(3.166) in frequency domain is as follows: Vw \Y{jtj) - M{ju;)\ < \AM{ju)\ where Y{s) = S,{y{t)}, M{s) = Z{m{t)} and AM(s) = 2{Am(t)}. Same procedure applies to Eq.(3.167) and gives Yie^"^) - M(e''^)| < |AM(e''^)| Vw (3.176) It is not difficult to show that a sufficient condition in the frequency domain corresponding to Eq.(3.176) and Eq.(3.167) is given by: < max AM{eJ'^) (3.177) Substituting + K{ei'^)G{ei'^) (3.178) and expressing the magnitude in d B , we obtain: max G(e>")€n Af(eJ'^) 1+ K{&'^)G{ô'^) -1 AM(e?'^) M(ei'^) (3.179) Which is the general robust performance design criterion in the frequency domain for the time domain specification Eq.(3.166). The simplification is seen by comparing with Eq.(3.160). If die pre-filter F(e-"^) is chosen to be: ^ then Eq.(3.179) simplifies to: max G(eJ-)en > ii(ei-) (3.180) where '^^^'^^ " 1 + ir(ei-)G(ei<^) is the sensitivity function. In a well-posed problem, one can expect F(e^'^), as defined in Eq.(3.180), to be a realizable transfer function; its response to the input r(t) is m(t), which is defined in Eq.(3.163), as the mean of admissible responses and is itself an admissible response. Therefore, i f the specifications Eq.(3.165) can be met by a realizable system, then Eq.(3.180) defines a realizable transfer fimction and can be designed beforehand. In addition to the specification Eq.(3.165) on the response to the reference input r(x), it is usual to require a certain minimum amount of damping in the response to the disturbance v(t), govemed by the sensitivity function S(e^'^). The specification can be expressed by a condition of the form: max max 5(e''^) < B (3.182) or more specifically max Sie^'^) <B(u) (3.183) Combining Eq.(3.181) and Eq.(3.183), we obtain the overall robust performance design criterion: max S(e^'^) < min{A(w),B(w)} Vw (3.184) Following Theorem 3.1, we have the following criterion: Criterion 3.4 (robust performance criterion for Q A R C ) : Assuming that all plants G(e''^) 6 H2 in the family, set or band: and that a particular controller A'(e"^) e H2 stabilizes the estimated plant G(e^'^) defined by Eq.(3.147). With prefilter chosen as Eq.(3.180), the closed-loop system will meet the performance specification Eq.(3.167) if the sensitivity function satisfies the following condition: In order to understand the relation between the robust stability and robust performance criteria, referring to geometric arguments in Fig 3.29, and the frequency function estimate Eq.(3.147), we find that: (3.187) Using multiplicative imcertainty description of Eq.(3.153), Eq.(3.187) then becomes: |1 + KG\ > 1 + KG KG A^G ^Genr, (3.188) With (3.189) 1 + K{ei'^)G{ei'^) and referring to Eq.(3.155), Eq.(3.188) can be rewritten as: < 1- ^mGiw) VG(e^'^) en. (3.190) or Vw max 1- AmGiu) (3.191) Using Ekj.(3.186) and defining: AB(a;) = — min{A(u;),B(u;)} Vw (3.192) then, substitute Eq.(3.191). we have: T(e^'^) AmG{c>j) + 5(e''") ABiu) < 1 Vw (3.193) Now we have derived the inherent relationship between robust performance and robust stability as well as nominal performance, which is stated in the following corollary. Corollary 3.5: Assuming that all plants G(e-"^) in the family, set or band ft„.- have the same number of clockwise encirclements of the point (-1,0) by their Nyquist contour (think of the same number of RHP poles in a model), and that a particular controller /ir(e-''^) stabilizes the estimated plant G{e^'^). Then the closed-loop system will meet the performance criterion max 5(e^'") I < min {A(w), B(u;)} Vw (3.195) if and only if the sensitivity function S^e^^) and the complementary sensitivity function T(e-''*') for the estimated plant G(e-"^) satisfy T(e''^) AmG{io)+ 5(e^'^) A5(u;) < 1 Vw (3.196) In Eq.(3.196), AmG(w) is related to uncertainty bounding estimation and AB(u>) is related to closed-loop tolerance and disturbance rejection. From Eq.(3.196), we can see that robust performance criterion implies robust stability and nominal performance. The interdependence of 5(e-'") and S{e^^) +i{e^'^) = 1 (3.197) makes it a challenge to meet Eq.(3.196). For instance, with a given robust performance specification, improving the performance by increasing open loop gain K{e^'^)G{e^'^) , which decreases S{e^") A 5 ( w ) , will worsen the robusmess, i.e., increases ï ( e ^ ' ^ ) A^G(c<;), and pushes the system closer to the point of instability for some G(e^^) £ fi^. Therefore, the trade-off presented by Eq.(3.196) is readily appreciated. Robust performance is achieved simply by satisfying botii nominal performance 5(e''^) AB{u}) < 1 and robust stability T(e-''^) AmG{ijj) < 1 with some margin; diat is, i f 5(e>'^) A5(a;) < /?(w) and T(eJ'^) A „ G ( w ) < 1 - y9(w) widi /?(w) < 1, dien robust performance is automatically guaranteed. 3.5 Implementation of Quantitative Adaptive Robust Control 3.5.1 Off-Line Design: Open Loop Shaping In this section, we still use X{e^^) to notate the discretized system X{e^^'') for simplicity although there is no information available beyond ùj>ir (T, = 1 is assumed in this thesis). The objective is to shape the loop transfer fiinction i(e-''^) to meet the robust performance Criterion 3.4 which implies stability from Corollary 3.5. The design procedure is best presented by the following example. We begin with the specifications for Eq.(3.163) and (3.164): m{t) = 1 - e-"'^^' (3.198) ' Am{t) = 0.75e -0.75« (3.199) so ^-0.75 AM(e-) = -0.75 1 - e -0.75 (3.200) (3.201) gju> _ e-0.75 and e-''^ = z as usual. With the pre-filter chosen as: -0.75 R{Z) ^ - ^-0-75 (3.202) the A(a;) in Criterion 3.4 (section 3.4.2) becomes: A(a,) = AM{z) M{z) (l_e-o.r5)(,_i) + 1 - 2e-o-75 (3.203) which is shown in Fig 3.33 in dB. Then the B(u>) in Eq.(3.183), the constraint on the response to the disturbance, is simply assigned as a constant: (3.204) On a second-order model, this would correspond to a damping factor ^=0.3 in the response to a disturbance at the output. The performance criterion Eq.(3.186) in dB becomes: 0.5276 (eJ'^ - l ) 0.4724eJ'^ -|- 0.0553 dB 20*logA(w) 10 Figure 3.33: Bound A(w) versus w In order to use die Nichols Chart to deal with system uncertainty, we define: r(^) = 1 L{z) (3.206) where L{z) = K(z)G{z), the loop transfer fimction; then, Eq.(3.205) becomes: max G(ei")efi„ 1 + r(e-''^) < min dB 0.5276 (e''^ - l ) 0.4724eJ'^ + 0.0553 ,6dB} Vw (3.207) dB Therefore we can foUow the design procedure in Q F T to design F . Fig 3.34 shows the Nichols Chart loci b(') derived for T in Eq.(3.207). The difference from die regular Q F T design is diat we obtain an upper boimd b(«) instead of a lower bound in F i g 3.32. Note the upper boimd is derived as soon as Eq.(3.207) is available, which is a much easier procedure than that in QFT design. For example, die locus b(0.15) at w=0.15 is obtained as follows: from Eq.(3.207), for u;=0.15, r(eJ0.i5) max 1 + r(eio-i5) < - 1 6 . 6 dB (3.208) dB therefore, r(e-'°-^^) must lie on or below locus b(0.15) which is exactly the contour of constant l n | r / ( l + r ) | = -16.6 on die Nichols Chart. Thus, i f aigT(e^°-^^) can be no larger than -17 dB, i.e. point " O " on Fig 3.34. = - 9 4 ° . dien \T{e^'^-^^)\^g Since \L{e^'^)= -|r(e"^) it is not difficult to see that bounds b(w) for i(e''*') can be derived by turning the whole picture 180° around 0 dB axis, which is shown in Fig 3.35. Examining the contours of constant magnitude of the closed-loop fiinction on the Nichols Chart, we note that i f argX(e'°-i5) = - 9 4 " , tiien \L{e^^-^^) can be no smaller tiian 17 dB, i.e. point " O " on F i g 3.35. So locus b(0.15) is obtained as the lower boundary of the admissible region for L (e-'°-^^). Comparing to the traditional QFT method of using templates, we get a lower boundary with much less effort. The disturbance condition Eq.(3.204) gives rise to the closed contour marked locus b(>1.5) on Fig 3.35. For w>1.5, the loci b(a;) are contained outside contour b(>l,5). Hence, in view of Eq.(3.205), the design of Ujw) is based on the loci b(w) for w<1.5 and on the closed contour b(>1.5) for a;>1.5. -350" -300" -250" -200" -150" -100" -50" 0" Degrees Figure 335: Design Loci b(') and L(e^'^) on the Nichols Chart As we mentioned at Section 3.4.1, we are still facing the problem of how to determine a imique loop transfer function X(e^'^). To deal with diis problem, first, we simply revise some optimum properties of i(e-"^) finom [Horowitz and Sidi, 33; Horowitz, 32]: Boundaries on X(e-''^) have die property diat d\L{€^'^)\/doj < 0; Boundaries on L[e^'^) have portions of both positive slope and of negative slope. A n optimum £(e-''^) first crosses boundaries where the slopes are positive, followed by the crossing of boundary where the slopes are negative. The above properties, plus practice and intuition, lead to the L(jw) in Fig 3.35. For instance, it should be recognized that the phase at zero frequency is fixed at -90° to meet the zero steady-state error requirement 3.5.2 On-Line Design: Controller Design In this section, we will develop an on-line method to use the plant frequency function estimation (?(e-'") and its uncertainty bounding estimation AG{oj) or Am(7(w) to carry out the controller design, i.e., K(€^'^). As a bridge to the goal, we first seek the way to derive L[e^'^) which is the open-loop transfer function corresponding to the estimate G{e^'^): Studying R g 3.35 carefiilly, we can see that the bounds b(w) for a;<1.13 are lower bounds for the uncertainty range at each frequency (picture the template " B " outside the locus b(0.76) in F i g 3.35), i.e.: X(e^'*') = Lcie^"") V w < 1.13 (3.209) where the subscript ' L ' means that i£,(e-''^) should have the smallest magnitude for the entire uncertainty range at frequency w. Furthermore, we can see that the bounds b(a») foru;>1.13 become upper bounds (picture the template " A " outside the locus b(1.13) in Fig 3.35) since the entire imcertainty range at w should be outside the contour b(w) for c<;>1.13, i.e.: V w > 1.13 (3.210) where subscript ' U ' means that Lu (e-''^) should have the largest magnitude for the entire uncertainty range at frequency w. To be sure tiie bounds E:q.(3.209) and (3.210) hold for all tiie fimctions in il^n defined by Eq.(3.152), the following relation has to be true: > L{e^'^) - AL{e^'^) V w < 1.13 < L{e^'^) + ALie^'^) V w > 1.13 where Ai(e-"^) = A'(e-"^)AG(a;). Using multiplicative uncertainty description E:q.(3.153), the above equation becomes: (1 - A,„G'(w))|2(e^'^)| > \L{e^'^)\ (1-f-A„G'(w)) <|X(e^'^)| > 1.13 Va;<1.13 Therefore, an optimal loop transfer function is derived for 2(e-''^): tit) (3.212) where L{e^'^) is die designed loop transfer function from off-line design. ATOG(W) is from robust estimation. So eventually, Lopt(e^'^) can be manipulated on-lme. There will be a straight line in Nichols Chart F i g 3.35 for Xop/(e-"^) around The length of this line depends on the size of the uncertainty bounding = 1.13. Since this discontinuity is physically unrealizable, the straight line will be smoothed by a rational function fitting approximation: N. with wi and u>2 chosen to cover die discontinuous point, we change Ng to match the size of AmG{cjs). The relation is as follows: Ns = Integer / l o g [ ( l - 2 A G ^ K ) ) | Z + , (ei-^ ) / X + , (e^-^ ) ^ log [|^(ei-OU.=i] J (3.213) which can be also calculated on-line. Finally, the real L T F (loop transfer function) is derived as (here we stiU use X(e-''^') to denote it): X[e^'^" ) = L(e^'^" )/(I - AmGiu^k)) X(e^'^») = X + , ( e ^ - O i . ( ^ " * ) L{e^'^^) = X(e^'^')/(1 + AmGM) V a;^ < c^i V w i < wfc < (3.214) V u ; , > u;^ where X+,(e^'^^)i,(e^-») « X ( c ^ ' ^ i / ( 1 + A^G'(a;2)) Fig 3.36 and Fig 3.37 are only for demonstration, where we assume that the dotted line is the designed L T F i(e-''^) and the solid line is the corresponding optimal L T F Lopt{e^'^)- We can see that the real L T F î ( e J ' ^ ) is foUowing Lopt{e^'^) closely. The frequency function of controller K{e^'^'') at each frequency a>fc is calculated from: (3.215) G(eJ'^*) By observing Eq.(3.215), minimum phase assumption has to be imposed on the process. So far, we come to the same end result as QFT, that is, a controller with a numerical form which is difficult to use in a practical adaptive context. To deal with this problem, we introduce an auxiliary hypothesis about the real controller, that is to assign the controller an order "n". The problem is then to fit a parametric frequency response model to the numerical frequency response points. The price for this rational function approximation is then the risk of introducing a bias. A least-squares type fitting can be used by assembling a set of " N " numerical frequency response points and solving for the unknown coefficients of a "n"th order transfer function model. Clearly, to expect a solution, " N " has to be greater than "n". Also, " n " is the number of pairs of lead-lag compensators, so keeping a minimum "n" has a very practical sense. The rational function approximation procedure is illustrated as follows: considering our controller model as ^' ' 1 - aiz-^ - a2Z-^ o„_i^-("-i) - a „ z - " (3.216) where: 0 = [ai,a2, - • • •,<l'n-l,0,n,bQ,bx,b2, - • • ,bn-i,bnf^ Using this structure, we can write: K{e,z-') = rP''iz-')0 where: i,{z-') = [z-'K{0,z-'),---,z-^K{0,z--'),l,z-\---,z-^f since K(0,z~^) and z~^ are complex numbers and 0 are real numbers, we have: Im{K{0,€-^'^')) = / m ( V ' ^ ( e - ^ ' ^ ' ) ) ^ (3.217) Thus, if we know the value of /v e"-''^* ) for some known <jjk, we can find two sets of linear equations in the parameters 0. However, in practice, we only have numerical frequency function K{e^'^'') at frequency points ojk for k = l,2,3,---,N (3.219) and N > n-\- 1 So we define the 2Nx(2n+l) matrix widi £'(e-''^*) as the variables: • $(w) = Re{e-^'^'K{e^'^^)} • • • Re{e-J'"^^K{e^'^')} 1 Re{e-^'^^} • • • Re{e-^"^'} i?e{e--"^'^7i:(eJ'^'')} ••• Re{e-^'"^''K(e^'^'')} 1 Re{e-^'^''} ••• Re{e-^'^''} /m{e-J'^» K (eJ'^' )} • • • / m { e - - ' " " ' K{e^'^' )} 0 Jm{e--"^'} • • • Im{e-^"^' } and the 2 N x l matrix: i2e{/i'(eJ'^')} Re{K{e^''^)} (3.221) Im{K {e^'^^)} Then, we have a well known result: 0 = r$^$i (3.222) Therefore, this rational controller transfer function approximates the controller fix)m the frequency domain curve fitting in a minimum variance sense. We have presented a complete design routine for Q A R C and the final resvdt is a controller model. The proposed frequency domain robust estimator provides the information about system uncertainties that Q A R C uses in its controller design in a real time manner. The price we pay for this novel capability is the extensive frequency domain calculation. 3.5.3 Q A R C Design Procedures In short, the whole design procedure can be descrit)ed by the following steps. A. Off-line design: a. Tune domain specification for closed-loop: Specify step response m(t) for unit step reference r(t); • b. Specify allowed tolerance bound Am(t); Associated frequency domain specification: (here r(t)-»i2(e-''^), m(t)-^M (e-''^) and Am(t)^AM(eJ'^)) • Choose pre-filter as F{z) = Specify the Hoo bound for sensitivity function: max Sie^'^) < B V w; Derive the up-bound iH2) of sensitivity function: max Siei^) < min | ^ ^ ^ C i , 51 V w; Derive the bounding L T F L{é''^) from Nichols Chart. (To guarantee the existence of Z(e''^), it is better to make it a realizable rational function at this stage); Fmd the firequency point where L{€^'^) changes from lower boimd to upper bound, (for the example, this frequency is 1.13); • Try out frequency points wi and u>2 with ioi < u)2, construct: and make sure the fast decline rate at w,; B. On-line design (simplified version*): a. Design of smooth function L,(e"^'') and real L T F 2(é"^*): Calculate A „ G ( w J ; *By simplified version, we mean we will use AmG{u>3) instead of AmG{uk) in Eq.(3.212) and Eq.(3.214), which makes easy to realize. Take the optimal L T F as: Calculate N.: Ng = Integer A^g - 2AG.(a;.))|X^,,(M)/X.^,.(iu'2) log • b. [\L,{M)\N.=I] Calculate the real L T F : Design of controller K(e^'^): • Calculate the numerical form (N points) K [e^'^'- ) : Finally, derive the «th-order rational fimcdon form K (0, e~^'^'' ) through L S curve fitting. 3.5.4 Simulations and Discussions For the examples, we take the algoridun in Section 2.3.2 as the estimation part. The closed-loop system configuration is the same as shown in Fig 3.30. Example 3.1 The input-ouçut model is same as Eq.(2.15) except for g is variable: y(t) = 0.9y(^ - 1) -H su{t - 1) (3.223) As first example, we only show the quantitative robust property, i.e., we fix g = 0.1; we then use the estimated frequency function and uncertainty bounding function as the process constraints; after that, we show diat the Q A R C closed-loop performance is widiin the specifications, Ekj.(3.167). The process constraints is shown in Fig 3.38 which is the same picture as Fig 2.10 where we use three Laguerre filters (N=3) for frequency function estimation and five filters (M=5) for the uncertainty bounding estimation. The assumed controller model is: Chapter 3: Quantitative Adaptive Robust Control in Frequency Domain Gain and Bounding solid line: real system (1st order) ^TT^it'tirTiïtTr^^ dashed line: estimation 10-1 dotted line: boimding J I—L Figure 3.38: Estimation of Frequency Response and Bounding with Laguerre Filter N=3; M=5 Here, we take n=l (controller order) and Nu=26 (frequency points, 'N^' is used here for number of frequency points we take, and N is the number of Laguerre filters for frequency function estimation). Going through the whole design procedure from Section 3.5.1 to Section 3.5.2. First, we impose the closed-loop specification as Eq.(3.198) and (3.199). Then, the pre-filter is obtained from Eq.(3.202). The boundary L{e^'^) is derived in Fig 3.35 and A'(e"^'=) is calculated from Eq.(3.215). Finally, we estimate controller parameter from Eq.(3.222): = 1.0, 6o = 8.583 6i = -8.215 With pre-filter F{z) and controUer K{e,z-'^) as: F{z) = 0.4724Z + 0.0552 (z - 0.4724) 8.583 - 8.215z-^ 1-2-1 (3.225) 1 y(t) 0.8 - / 0.6 - 0.4 • 0.2 - 0 - -0.2 ./a(t) - -0.4 I 0 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 2 4 6 8 10 12 14 16 1 1 t 18 1 20 Figure 3.39: QuantitaUve Robust Property of QARC the process output is shown in F i g 3.39. h i Fig 3.39, referring Eq.(3.198) and Eq.(3.199), we notice that the closed-loop response happens to lie between the two lines (not necessarily): b{t) = m(t) + Am(t) a{t) = m{t) - Am{t) Obviously, the closed-loop specification is fulfilled in this case. Example 3.2 We take the same design procedure as Example 3.1, but use it in an on-line situation where we change the plant parameter g in Eq.(3.223) from 0.1 to 0.2 at t=60sec. We use eight Laguerre filters (L=8), so that there are eight Laguerre filter gains to be estimated. The first three (N=3) are used for frequency function estimation and the other five (M=5) for uncertainty bounding estimation. In die Laguerre gain R L S estimation algorithm, the diagonal element in matrix P(0) is 10^ and the forgetting factor A is 0.95 as the initial setting. During the first ten steps, the loop is open; the process input is white excitation. The Q A R C is turned on at t=10sec and the loop is closed. The process output is shown in F i g 3.40. In order to study the performance of the robust estimator process gain doubled r(t) / y(t) 05 QARC turns on -05 20 40 60 80 100 120 140 Figure 3.40: Adaptive Property of QARC in the closed-loop situation, the estimated frequency function is shown in Fig 3.41 and die estimated boimding is shown in Fig 3.42. It is understandable diat the bounding estimation needs some time to adjust itself after the sudden process dynamical change just same as the case in frequency function estimation. This effect is also seen at the output during 60sec to 120sec. Another observation is diat the set-point change helps the estimatioa Example 3.3 version with In this example, we apply the Q A R C strategy to Rohrs* model Eq.(2.75). The discrete = 0.2 sec of the model is as foUows: , _ 0.1584(1 -t- 0.9952-?-^) ( l + 0.0449^"^)^-^ ^ ~ (1 (1 - 0.81872-i)(l 0.81872-1 )(1 -h 0.0917^-1 - 0.0025Z-2) 0.0025Z-2' (3.227) The robust estimator is implemented by Laguerre filter method. In this case, we use ten Laguerre filters (L=10), among them, five are used for frequency fimction estimation (N=5) and five used for uncertainty bounding estimation (M=5). The initial setting is as for the previous example and the closed-loop specification is the same as in Example 3.1. In this example, we try to show the effect of the controller order on the closed-loop performance. With a second-order controller, the process Figure 3.41: Frequency Function Estimation in a 3-D representation Figure 3.42: Uncertainty Bounding Estimation in a 3-D representation Figure 3.43: Performance with a second order controller 1.5 1; ,'! 1 1 1 1 1 ; • r(t) ; \ M ? 1 t ;'\ i i ! ; 1; t 0.5 1 ' ' 1 ! \' 'I i -0.5 1— 20 ; y(t) i \ 1 1 ( I I 1 ' 1 1/ N=5; M=5 I 40 1 60 1 80 1 100 1 120 140 160 Figure 3.44: Performance with a fifth order controller 180 200 output is shown in Fig 3.43. With a fifth-order controller, the process output is shown in Fig 3.44. We can see that the closed-loop response is closer to the specification with a higher order controller. The explanation is that the frequency domain curve fitting Eq.(3.222) has a smaller variance when the controller order gets higher. So far, we only deal with stable, linear and minimum phase system. However, the Q A R C has been shown to be a feasible design strategy. It is important to have a reasonable close loop specification in the first place, i.e., for a higher-order process, the specification Eq.(3.198) and (3.199) have to be a litde more complicated in order to have an acceptable control action. Accordingly, die assumed controller order 'n* also must increases. 100 101 102 100 101 102 Figure 3.45: Magnitude Plot of Frequency Response and Uncertainty Bounding for Rohrs' Model Example 3.4 In this example, we use the same plant as the last example, but widi different choice of N and M . The estimation part for each case is shown in Fig 3.45 (same as Fig 2.17). The closed-loop response for each case is shown in F i g 3.46. For a measure of both performance and stability, the 0.5 - 03 simulation A simulation B N=2; M=8 N=4; M=6 0 10 0.5 15 20 10 0.5 simulation C 20 simulation D N=8; M=2 N=6; M=4 10 15 15 20 10 15 20 Figure 3.46: Closed-Loop Step Response from QARC 0.5 0.5 - simulation A N=2; M=8 10 -I 0.5 15 simulation B N=4; M=6 JJ 20 10 15 —t r- 05 - simulation C 20 1- simulation D N=6; M=4 N=8; M=2 0 -I 10 15 20 _l 1_ 10 Figure 3.47: Closed-Loop Step Response from Regular ARC 15 20 Rohrs-model 1st order model >l< 0.5 N=8; M=2 -0.5 20 40 60 80 100 120 140 160 180 200 Figure 3.48: Switch From Rohrs' model to 1st order system difference is not very significant because of the overall robust stability and performance consideration of the Q A R C design. In another words, the Q A R C design not only cuts the L T F gain for the high frequencies to guarantee the stability, but also increases the L T F gain in low frequencies according to the size of imcertainty. This improves the performance of the system, which is superior to regular robust control methods. Also, it is very easy to reduce Q A R C to regular robust control. We simply take the last equation in Eq.(3.212) (3.214) and use it for the whole design. This feature can be clearly shown, as in the Nyquist plots Fig 3.36, by pushing the jump point a», to the very beginning. The associated closed-loop performance is shown in Fig 3.47. Comparing with Q A R C in H g 3.46, the performance becomes worse as uncertainty gets bigger. As seen from both Fig 3.46 and 3.47, the system performance gets better as N increases and M decreases. This eventually matches the discussion about the interrelationship between N and M in Section 2.3.2. Example 3 J In this example, we will simulate the situation when the controlled plant is switched fiiom one model to another. The purpose is to test the adaptive robust property of the Q A R C design. Fig 3.48 and Fig 3.49 contain the results when the plant is switched finom Rohrs' model to a 1st ; '. Rohrs-model > |< 1st order model !! 11' '<(• 05 N=2; M=8 -05 0 20 40 60 80 100 120 140 160 180 200 Figure 3.49: Switcli From Rohrs' model to 1st order system widi N=2 and M=8 15 1 I. Isl order model > |< 1 1 1 140 160 180 Rohrs-model :a: !'•' 1 1 i! i / t1 \ • ' iiii ;^ 05 1 1 1 1 \ 1 1 ' 1 1 1 1 ! 1 1 : \ 11 \i N=8; M=2 20 1 1 40 60 1 80 L 100 120 Figure 3.50: Switch From 1st order system to Rohrs' model 200 1.5 order model Rohrs-model >l< N=2; M=8 -0.5 0 20 40 60 80 100 120 140 160 180 200 Figure 3.51: Switch From 1st order system to Rohrs' model with N=2 and M=8 order model. Fig 3.50 and Fig 3.51 contain the results when the plant is switched from a Ist-order model to Rohrs' model. For both examples, for fixed L=M+N, the simulation with N=2 and M=8 is more tolerant than that with N=8 and M=2, i.e., the system robustness gets better as M mcreases and N decreases. This matches the discussion about Ûie interrelationship between N and M in Section 2.3.2., i.e., the choice of ' M * or ' N ' can be used as a tuning parameter between system robustness and system performance. 3.6 Discussion about Non-Minimum-Phase Systems Since Q A R C follows almost the same design procedure based on QFT, instability arises when the controlled system is non-minimum-phase (refer to Eq.(3.215)). A method was proposed in [Horowitz and Sidi, 34] to deal with the problem in QFT, i.e. to separate the N M P (non-minimum-phase) part, the imstable zeros from the system. In order to keep the remaining M P (minimum-phase) part with the same magnitude characteristics as the original N M P system, an aU-pass filter is used. Then the controUer design is merely based on the M P part The controller itself is stable and closed-loop Stability is guaranteed as far as the phase margin is larger than the phase-lag represented by the all-pass filter. The problem is that the bandwidth of the closed-loop system is very hmited. In the context of Q A R C , an additional problem is that we are not supposed to know the locations of system poles and zeroes. Inspired by the above idea, we now discuss a way to deal with N M P problems in Q A R C . The idea is straightforward i f we combine the aU-pass filter in [Horowitz and Sidi, 34] and Smith's Principle [Marshall, 48]. If the system G(z) has one or more unstable zeros at a i , 02 ,• • •, a , , let them be shown explicitly by writing G{z) = Gx{z){ai -f- z){a2 + z)---{a, + z) g q = Gm{z)A{z) (3.228) where: Gm{z) = Gi{z) Xl {I + Qiz) IS a M P system; .=1 ^i^) = ft T^fS «•=1 is an aU-pass filter. The phase of A(z) decrease monotonicaUy as frequency a;, increases and the magnitude \A{z)\ = 1 for all Wfc. In other words, the only difference between the original N M P plant G{z) and the M P Gm{z) is the phase-lag produced by the all-pass filter A{z). This leads to the same situation as timedelay systems where we can use Smith's Principle [Marshall, 48]. More specifically, we denote: A{z) = e-J'^('^*) (3.229) where: </>(wfc) is the phase-lag of A(z) at w , . Based on Smith's method, the controller is derived for G(z): ^ ' l + Km{ei'-'')Gm{e^->'){l-e-i'^M) ^ where: 7i„(e-''^'') is the controller designed for Gm{z) Example 3.6 We still use die Rohrs' model Eq.(2.75) as an example. This time, die sample interval Ta = 0.1 sec and the discrete model is: 0.037(1+ 1 . 8 2 1 9 z - i ) ( l + 0 . 1 1 6 4 . - 1 ) . - ! (1 - 0 . 9 0 4 8 . - i ) ( l - 0.4374.-1 + 0.0498.-2) ^^ ^ which appears to tiave an unstable zero at a=1.8219. The design procedure is as follows. From the estimation, we have full knowledge about the transfer function G(e^'^'') at u>k. The corresponding M P plant, which has the same magnitude as G(z), is: _ ""^^ 0.037(z-i + 1.8219) (1 + 0 . 1 1 6 4 2 - i ) z - i (1 - 0.9048^-1 )(1 - 0.4374Z-1 + 0.0498z-2) If we know the phase difference <f>{u>k) between N M P ^^'^^^^ and M P Gm (e-"^*), we can calculate the M P transfer fimction G„»(e"^*) from G(é"^*) by Eq.(3.229): Gm (e^"" ) = G{e^'^- ) e''^^'^*) (3.233) In our frequency domain estimation, some approximation methods have been developed to estimate the phase shift between a N M P plant and corresponding M P plant based on Bode relations or Bode integrals (see Appendix A ) . In this example, we assimie that we know: ^•^(^0 = i ± I : Ë 2 1 9 ^ eJO.K., +1.8219 (3 234) ( ^ and use it to calculate G^Ce""*). Then, E:q.(3.214) and Eq.(3.215) is used to calculate Kmie^'^"). After that, the controller function for the original N M P plant is derived by substituting Gm{e^'^'') and Km (e-''^' ) to Eq.(3.230). Finally, the rational function approximation is completed by following the same procedure shown in Section 3.5.2. Fig 3.52 is the closed-loop response with the controller order set to be n=2. Fig 3.53 is the response with the order n=5. The advantage of the method is that we can treat the unstable zeros and time-delay in a same way. That is because, in the frequency domain, an unstable zero modelled by an aU-pass filter has the same effect on the system characteristics as time-delay; both can be represented by <p{ju>k)- So eventually, we do not have to deal with unstable zeros and time-delay separately. When applying the Smith's Principle to time-delay systems, the stability of N M P system is equal to that of the corresponding M P system if we can accurately model the system G(e^'^''), phase-lag 4>{joJk) and controller K[B,e~^^''). In the presence of system mismatch, the problem 'reduces' to the nature of the characteristic equation: 200 Figure 3.52: Performance with a second order controller -1 T" -1 1 r- 1.5 i<t) 0.5 y(t) -05 - -ii 0 . 20 . 40 -J 60 80 100 120 I 140 J 160 Figure 3.53: Performance widi a fifdi order controller L_ 180 200 In this thesis, we still restrict Q A R C within the limit of M P systems since first we have to improve the estimation accuracy of the phase shift caused by N M P dynamics. Also in a nonparametric context, we cannot derive a rational characteristic equation like Eq.(3.235). Moreover sensitivity to mismatch is a very difficult subject i f we consider all the mismatches in a N M P closed-loop system. 3.7 Conclusions In this chapter, we have derived the robust stability theorem and robust performance criterion for Q A R C . To apply the theorem and criterion, a method is developed from QFT. In order to implement the Q A R C in real time, the time domain sensitivity specification is changed to a frequency domain sensitivity specification, i.e., from Eq.(3.167) to Eq.(3.176). B y means of the Nichols Chart, we are able to realize a complete adaptive robust control strategy which can make use of the imcertainty bounding information provided by the robust estimator in the previous chapter. The quantitative robust property and adaptive property have been demonstrated by five examples. We should mention here that the quantitative robust property is met only when the estimator converges. This is because, in the learning period after the sudden dynamical change, the estimator needs certain time to follow the change as shown in Fig 3.41 and F i g 3.42. Compared to regular robust control, one advantage of the Q A R C design is that it does not sacrifice too much performance for a large stability margin. The trade-off can be clearly seen from Nyquist plot Fig 3.36, i.e., we increase the lower frequency magnitude to compensate the higher frequency magnitude lose. So, with large uncertainty, the system response is not necessarily very slow as shown in Fig 3.46 comparing with Fig 3.47. Also the robustness is preserved as shown in Fig 3.48, 3.49, 3.50 and Fig 3.51. For a complicated system, the problem is the fitting error between a series of frequency points and a structured frequency function. One way to deal with it is to increase the controller order as shown in Fig 3.43 and Fig 3.44. However, this error is not included in the imcertainty bounding function. A drawback of Q A R C is the lengthy procedure of the off-line design on the Nichols Chart and the on-line block L S fitting of rational controller fimction. On one side, the higher-order controller means better approximation of frequency fimction; on the other side, it means we have to have more complicated specifications and more time in real-time calculation. Chapter 4 Conclusions 4.1 Summary of Thesis In this thesis, we have developed a frequency domain adaptive robust controller which we call quantitative adaptive robust controller. The results of the thesis can be summarized into two parts: the frequency domain robust estimation in chapter two and Q A R C design and related theoretical supporting material in chapter three. In frequency domain robust estimation, we have described and demonstrated several frequency domain methods which provide not only the system frequency function estimation but an imstructured uncertainty bounding function and external disturbance fimction estimation. Different orthonormal functions, such as sinusoidal function and Laguerre functions, are used to obtain die results. Compared with the time domain methods based on a parametrized model, we need less a priori knowledge about the process. As a main result, the estimated uncertainty bounding can capture the modelling error. Besides its usefiilness for robust control, it can also be used to measure the goodness of the estimation. In the Q A R C design, we have developed a real-time frequency domain synthesis design method which makes use of the modelling accuracy information supplied by the frequency domain robust estimation. The final product of this research is an adaptive robust controller which satisfies the prespecified closed-loop performance out of a stable closed-loop system. One advantage of the Q A R C design is that it tries to sacrifice less performance to meet the stability condition. At this point, our research has been conducted only on SISO systems and the process is assumed to be stable, linear, minimum phase and time-invariant system. Looking into recent developments in adaptive robust control, e.g., the 1991 European Control Conference, the focus is still on the estimation of uncertainty bounds and the implementation of adaptive robust control. However, most of the time, the knowledge about the shape of the true frequency response is still required to estimate transfer fiinction error bounds [Wahlberg and Ljung, 67]; and the full use of the bound is not in the frequency domain implementation of adaptive robust control [Wittenmark and Kâllén, 70]. Adaptive robust control schemes making use of the system uncertainty bounding with no or litfle knowledge about the system are not very common yet. Hopefully, the Q A R C developed in this thesis could enrich the dialog between the previously separate fields of identification and robust control. 4.2 Suggestions for Further Research In frequency domain robust estimation, as discussed in appendix, one direction for future research is the possibility to estimate the time delay, A R M A model and associated structured uncertainty bounding function. The reason we leave it in the appendix is that it is still a very intuitive approach and much theoretical analysis and improvement are still needed. Ordy then, wiU we be able to have a complete frequency domain estimation package. Our research in the Q A R C design is still at an early stage, so that we only choose a prefixed noise attenuation factor for all the frequencies. As one direction of future research in this part, we will explore the possibility to take into consideration extemal disturbance estimation to add frequencydependent weighting to the sensitivity function. The difficulty we can foresee is how to implement it in real-time. It is very important to look ahead at the possibility of extending the method to M I M O situation. In the estimation part, we already have multidimensional DFT and Kalman filtering, M I M O Laguerre filters and LS algoridun. What is needed is a M I M O uncertainty bounding estimation. In the Q A R C design, for the robust stability theorem and robust performance criterion in chapter three, we need to change from transfer function, imcertainty fimction, sensitivity function and complementary fimction to transfer matrix, uncertainty matrix, sensitivity matrix and complementary matrix correspondingly. Theoretically speaking, all the theorems and criteria still hold. The H2 framework changes to Hoa frameworic. The loop transfer fimction shaping changes to the shaping of the singular values. Then, the problem remains to be that a real-time M I M O synthesis method is not available. There is lot of woric to be done in this direction. Another problem standing in the way of extension to M I M O systems is the complexity of the frequency calculation. It is even a potential problem in the SISO case i f the system has fast changing dynamics and environment. So anotlier direction for fiiture woric is to increase the speed of the algorithm. This research should thus benefit from the newly developed parallel process techniques. The highly parallel structure becomes evident from the beginning of chapter two. The DFT Eq.(2.10), Kalman filter Eq.(2.24) and R D F T Eq.(2.81) can be described as single input and ' N ' output matrices, i.e., they represent banks of *N' single input-single output filters which run in parallel. This parallelism is further exploited in imcertainty bounding estimation Eq.(2.74) and Eq.(2.101), and in frequency function and external disturbance spectrum estimation Eq.(2.137) because the function is estimated independendy at each frequency . Thus the frequency domain estimation is actually decomposed into ' N ' decoupled estimators operating in separate frequency bands. The parallel structure is also depicted in the Q A R C design in Eq.(3.214). The only non-parallel operation is the final L S rational function approximation which can be also implemented by a L S hardware processor. Therefore in order to deal with the speed problem in the future implementation, it is worth looking into the use of parallel processor. References [I] IEEE Special Issue, IEEE Trans, on AC, volume AC-16, D e c , 1971. [2] B . D . O. Anderson, R. R. Bitmead, C. R. Johnson, P. V . Kokotovic Jr., R. L . Kosut, I. M . Y . Mareels, L . Praly, and B . D . Riedle. Stability of Adaptive Systems: Passivity and Averaging Analysis. MIT Press, 1986. [3] K . J. Astrom and B . Wittermiark. On self tuning regulators. Automatica, 9, 1973. [4] K . J. Astrom and B . Wittenmaric. Self-tuning controllers based on pole-zero placement. lEE Proc, 127, no. 3, 1980. [5] K . J. Astrom and B . Wittermiaiic. Adaptive Control. Addison-Wesley Publishing Company, 1989. [6] M . Athans. A tutorial on the Iqg/ltr method. In American Control Conference, Seatde, Washington, 1986. [7] J. D. BirdweU. Evolution of a design mediodology for Iqg/ltr. IEEE Control System Magazine, April, 1989. [8] R. R. Bitmead and B . D. O. Anderson. Adaptive frequency sampling filters. IEEE Trans, on Circuits and Systems, 28, no. 6, 1981. [9] R. R. Bitmead, A . C. Tsoi, and P. J. Paricer. A kalman filtering approach to short-time fourier analysis. IEEE Trans, on ASSP, 34, no. 6, 1986. [10] H . W. Bode. Relations between attenuation and phase in feedback amplifier design. Bell Syst. Tech. J., 19:421^54, 1940. [II] H . W. Bode. Network Analysis and Feedback Amplifier Design. Princeton, NJ: Van Nostrand, 1945. [12] D. W. Qarke and P. J. G a w ^ o p . A self-mning controller. Proc. Inst. Elec. Eng., 122:929934, 1975. [13] D . W . Qarke, C. Mohtadi, and P. S. Tuffs. Generalized predictive control, part i : The basic algorithm. Automatica, 23, no. 2, 1987. [14] S. J. Cusumano and K . PooUa. Adaptive robust control: A new approach. In American Control Conf., Atianta, Georgia, June, 1988. [15] J. Doyle. Analysis of feedback system with structured uncertainties. lEEProc, 29, Pt. D , no. 6. 1982. [16] J. C . Doyle and 1981 G . Stein. Multivariable feedback design: Concepts for a classical/modem syntiiesis. IEEE Trans, on AC, AC-26(2), 1980. [17] L . Dugard, G . C. Goodwin, and X . Xianya. The role of the interactor matrix in multivariable stochastic systems. Automatica, 20:701-709, 1984. [18] B . A . Francis. A Course in Hoo Control Theory. Springer-Veriag, 1987. [19] Y . Fu. Robust adaptive control. PhD tiiesis. Dept. of E E , U B C , 1990. [20] Y . Fu and G . Dumont. A n optimum choice of time scale for discrete time laguerre network, to appear, IEEE Trans, on A C . [21] Y . Fu and G . Dumont. Robust discrete str design by using intermittent adaption. In IFAC Workshop on Robust Adaptive Control, Newcastie, August, 1988. [22] A . Gera and I. Horowitz. Optimization of tiie loop transfer function. INT. J. CONTROL, 31(2):389-398, 1980. [23] G . C . Goodwin and M . E . Salgado. A new paradigm for estimating restricted complexity models of dynamic systems. Technical report. The Univ. of Newcastie, 1988. [24] G . C. Goodwin and M . E . Salgado. Information and intelligent control. In International Workshop on Adaptive control Strategies for Industrial Use, Alberta, Canada, June, 1988. [25] G . C. Goodwin and K . S. Sin. Adaptive Filtering Prediction and Control. Prentice-Hall Inc., N.J., 1984. [26] M . J. Gottiieb. Polynomials orthogonal on a finite or enumerable set of points. Am. J. Math., 60:453-458, 1938. [27] M . J. Grimble. Eco robust controller for self-tuning control applications, part i i : Self-tuners and stability. INT. J. CONTROL, 46, no. 5, 1987. [28] J. Guckenheimer and P. Holmes. Nonlinear Oscillations, Dynamical Systems, and Bifurcations of Vector Fields. Springer-Veriag, 1983. [29] S. Gunnarsson and L . Ljung. Frequency domain tracking characteristics of adaptive algorithms. IEEE Trans, on AC, 37, no. 7. 1989. [30] J. W. Helton and A . Sideris. Frequency response algoridmi for Hoo optimization with time domain constraints. IEEE Trans, on AC, 34, no. 4, 1989. [31] I. Horowitz. Optimum loop transfer function in single-loop minimum feedback system. INT. J. CONTROL, 18(1):97-113, 1973. [32] I. Horowitz. Quantitative feedback dieory. lEE Proc, 129, Pt D , no. 6, 1982. [33] I. Horowitz and M . Sidi. Synthesis of feedback systems with large plant ignorance for prescribed time-domain tolerances. INT. J. CONTROL, 16:287-309, 1972. [34] I. Horowitz and M . Sidi. Optimum synthesis of non-minimum phase feedback systems with plant uncertainty. INT. J. CONTROL, 27:361-386, 1978. [35] G . H . Hostetter. Recursive discrete fourier transformation. IEEE Trans, on ASSP, 28, no. 2, 1980. [36] P. A . loannou and K . S. Tsakalis. A robust direct adaptive control. IEEE Trans, on AC, AC-31, 1986. [37] R. E. King and P. N . Paraskevopoulos. Digital laguerre filters. Circuit Theory and Applications, 5:81-92, 1977. [38] R. E . King and P. N . Paraskevopoulos. Parametric identification of discrete-time siso systems. INT. J. CONTROL, 30:1023-1029, 1979. [39] R. L . KosuL Adaptive robust control via transfer function imcertainty estimation. In American Control Conf, Atianta, Georgia, June, 1988. [40] G . Kreisselmeier. A robust indirect adaptive control approach. INT. J. CONTROL, 43:161175, 1988. [41] K . R. Krishnan and A . Cmickshanks. Frequency-domain design of feedback systems for specified insensitivity of time-domain. INT. J. Control, 25(4), 1977. [42] R. O. LaMaire, L . Valavani, M . Athans, and G . Stein. A fiiequency-domain estimator for use in adaptive control system. In American Control Conference, Minneapolis, M N , 1987. [43] I. D . Landau. A survey of model reference adaptive techniques. Automatica, 10, 1974. [44] I. D . Landau and B . Courtiol. Design of multivariable adaptive model following control systems. Automatica, Sept., 1974. [45] L Ljung. System Identification - Theory for the User. Prentice-Hall, Ehglewood Qiffs, New Jersey, 1987. [46] L . Ljung, B . Wahlberg, and H . Hjalmarssoa Model quality: the roles of prior knowledge and data information. In Proceedings of the 30th CDC, Brighton, England, 1991. [47] G . J. Macfarlane. Frequency-Response Methods in Control Systems. IEEE Press, 1979. [48] J. E . Marshall. Control of Time-Delay Systems. Peter Peregrinus Ltd, 1979. [49] R. H . Middleton and G . C. Goodwin. Digital Control and Estimation —A Unified Approach. Prentice-Hall, Inc. 1990. [50] R. H . Middleton, G . C. Goodwin, D. Q. Mayne, and D . J. H i l l . Design issue in adaptive control. IEEE Tram, on AC, 33, no. 1, 1988. [51] M . Morari and E . Zafiriou. Robust Process Control. Prentice-Hall, 1989. [52] K . S. Narendra, Y . H . Ym, and L . S. Valavani. Stable adaptive controller design. IEEE Trans, on AC, AC-25, 1980. [53] H . Nyquist. Regeneration theory. Bell Syst. Tech. J., 11:126-147, 1932. [54] K . Ogata. Modern Control Engineering. Prentice-Hall, Inc., Englewood cliffs, N J , 1970. [55] A . V . Oppenheim and R. W. Schafer. Digital Signal Processing. Prentice-Hall, Inc. 1975. [56] R. Ortega, L . Praly, and I. D. Landau. Robustness of discrete-time direct adaptive controllers. IEEE Trans, on AC, 30, no. 12, 1985. [57] R. Ortega and T. Yu. Theoretical results on robustness of direct adaptive controllers: A survey. In 10th IFAC World Congress, Munich. Germany, 1987. [58] R. Ortega and T. Yu. Robustoess of adaptive controllers - a survey. Automatica, 25, no. 5. 1989. [59] P. C. Parks. Lyapunov redesign of model reference adaptive control systems. IEEE Trans, on AC, AC-11. 1966. [60] E . Polak and S. E. Salcudean. Adaptive control of arma plants using worst-case design by semi-mfinite optimization. IEEE Trans on AC, 32. no. 5. 1987. [61] E . Polak and S. E. Salcudean. On the design of linear multivariable feedback system via constrained nondifferentiable optimization in Hoo spaces. IEEE Trans on AC, 2A, no. 3. 1989. [62] D. E. Rivera. J. F. Pollard. L . E . Sterman. and C. E . Garcia. A n industrial perspective on control-relevant identification. In Proceedings ofACC, San Diego, California. 1990. [63] C. E . Rohrs. Adaptive control in the presence of unmodelled dynamics. PhD thesis, DepL of E E , MIT, 1982. [64] C. E . Rohrs, L . Valavani, M . Athans, and G . Stein. Robustness of adaptive control algorithms in the presence of unmodelled dynamics. In Proc. 21st IEEE Corf, on CDC, Orlando, F L . , 1982. [65] M . G . Safonov and M . Adians. Gain and phase margin for multiloop Iqg regulators. IEEE Trans, on AC, AC-22(2), 1977. [66] G . Stein and M . Athans. The Iqg/ltr procedure for multivariable feedback control design. Technical Report, UDS-P-I384, May, 1984. [67] B . Wahlberg. System identification using laguerre models. IEEE Trans, on AC, 36, no. 5,1991. [68] B . Wahlbeiig and L . Ljung. On estimation of transfer function error bounds. In European Control Conference, Grenoble, France, 1991. [69] P. E . WeUstead. Non-parametric methods of system identification. Automatica, 17, no. 1, 1981. [70] P. E . Whitaker, J. Yamrom, and A . Kezer. Design of model-reference adaptive control systems for aircraft. Technical report. Report R-164, Instrumentation Lab. M I T , Cambridge, Mass., 1958. [71] B . Wittenmaik and P. KâUén. Identification and design for robust adaptive control. In European Control Conference, Grenoble, France, 1991. [72] N . E . Wu and G . Gu. Discrete fourier transform and Hoo approximation. IEEE Trans, on AC, 35, no. 9, 1990. [73] B . E . Ydstie. Extended horizon adaptive control. In IFAC 9th Triennial World Congress, Budepest, Hungory, 1984. [74] D . C. Youla, J. J. Bongiomo, and H . A . Jabr. Modem wiener-hopf design of optimal controllers; part i : the single-input-output case. IEEE Trans, on AC, AC-21(1), 1976. [75] G . Zames. Feedback and optimal sensitivity: model reference transformations, multiplicative seminorms, and approximate inverses. IEEE Trans, on AC, AC-26(2), 1981. [76] C . Zervos and G . Dumont. Deterministic adaptive control based on laguerre series representation. INT. J. CONTROL, 48, no. 6, 1988. [77] C. C. ZJSTVOS. Adaptive control based on orthonormal series representation. PhD thesis, DepL of E E , U B C , 1988. [78] Y . C. Zhu. Estimation of transfer functions: asymptotic theory and a boimd of model uncertainty. INT. J. CONTROL, 49:2241-2258, 1989. Appendix A Estimation of Pure Time Delay There are mediods in signal processing to estimate time delay by correlation techniques in the frequency domain. But they are not tailored for use in control systems because of additional plant dynamics. For our purpose, we have to devise a method to estimate the time delay. From classical control theory. Starting from continuous time domain, we know that for a minimum phase system the magnitude and phase-angle characteristics are direcfly related. This means that i f the magnitude curve of a system is specified over the entire frequency range from zero to infinity, then the phaseangle curve is uniquely determined and vice versa [Marshall, 48]. The so-called Bode Relationship [Bode, 11] are oo dlog\GiJo:)\ log du U?-|-Wi du (A.236) — CO where u = log / In the digital signal processing area, the similar relation for minimum phase system is called Hilbert Transform [Oppenheim and Schafer, 55] that has the form ^G{e^^^) = log|G'(e^-)|cot (j^)d^ (A.237) — IT where cot /a»-u;i\_ sln(a>-a>i) 2 J ~ 1 - cos (w - W i ) So, i f we know the amplitude plot in a Bode plot, the phase plot can be detennined exclusively. This observation gives us a good suggestion that the time delay can be distinguished from the other part of the plant by assumption that the other part of the plant has the minimum phase feature since the pure time delay only affects the phase value of the system: = 1 Iç-j'^r, ^ (A.238) (^239) where r j is the delay time. Here we assume the Bode relationship is used. The difficulty is how to describe the relation between amplitude and phase values in the frequency domain with a simpler relationship. Because our objective is to estimate the time delay rather than the exact function between amplitude and phase, some approximation methods for a particular model is shown in this sectioa Eq.(A.236) relates phase shift at frequency wi to the derivative of the gain characteristic on a logarithmic frequency scale. For example in high frequency, it is possible to derive an approximate relation for calculating phase characteristics: iGiJu) = wdlog\G{jw)\ 2 dlogu (A.240) From intuition, it is true in high frequency range for minimum phase system. Unfortunately, it gives a lot of bias in other part of the spectrum. Also with the limit of Nyquist frequency in our robust estimator, we can not get full information for the higher frequency range. Therefore, we have to look for the approximate relations in lower frequency. We are going to try two method. Since neither has a very sound theoretical background, we did not include this part in the thesis, instead here in this Appendix. Method One is a direct interpretation deducted from Eq.(A.236). Using series expansion, we can write: (A.241) Substituting Eq.(A.241) back to Eq.(A.236), we can still say that the phase characteristic is an odd function of the frequency because the integration is along u) rather than On the other hand, many process can be simply represented by simple model, low order with real poles and dead time. In industry process, most plants can be particularly described by a first order model plus dead time. The first order model without time delay can be described by (A.242) where T is the time constant The phase characteristic is LG{jij) = - arctanruj = - (A.243) The result is consistent with the interpretation from Bode relationships that the phase characteristic is an odd function of the frequency. Considering the fact that the process may have the integrator and dead lime, we can write a phase characteristic as a function of the frequency IG{3<^) = Co + C i w + Czu? + Csa;^ + • • • Referring to Eq.(A.243), i f (A.244) < 1, it is reasonable to use first three terms as an acceptable approximation in the lower frequency range Z G O w ) a Co + C i w + Caw^ (A.245) Where Co is the phase shift caused by the integrator, it is zero i f the process is first order plus dead time. C i is the part caused by dead time and process dynamics. C3 is only affected by process dynamics. Actually, referring to Eq.(A.239) and (A.243), C i can be written as C i = r + rd (A.246) From the phase characteristics of the frequency response, C o , C i and C3 can be estimated immediately from Eq.(A.245). Then, r can be carried out from C3 and, finally, r j is obtained from Eq.(A.246). Method Two is used to estimate phase characteristics from gain characteristics for the case that the delay-free part is minimum phase. Then the phase shift caused by time delay is the difference between measured phase shift and estimated phase shift. Still from delay-free system Eq.(A.242), we have: dlog\G(joj)\ du dlGUu) doj O;T^ 1 -I- w2r2 ^ T 1 -I- w2r2 (A.248) If w r < 1, we obtain the relation: dlGiJu) (A.247) 1 dlog\G{joj)\ v'-21og|C(jw)| du a;<0.5 plant I plant I I plant m Method One Td = 1.0496 Td = 1.8006 Td = 6.4894 Method Two Td = 0.90738 Td = 1.40739 Td = 7.07542 Figure A.54: Time Delay Estimation Therefore, refer to Eq.(A.239): V-2log|G(ju>) To eliminate the effect of noise, we average (A.250) in the lower frequency range: M dZG(M) - , \ dlog G C M ) (A.251) ^^ = - M 2 . A:=l where dlGiM) = IGiJojk+i) - iGiJUk) dlog Gijuk) =log (7(^+1 ) -log G{ju;k) For demonstration, we take the following two plants as examples, and we calculate the frequency response direcdy from them without going through the frequency domain estimation algoridun. plant I = plant II = s+ 1 1 (- + 1)^ plant III = (^ + 1)' The estimation is shown in Fig A.54. Now, the dead time TJ, of the process is derived by manipulating the phase characteristics. The accuracy of dead time estimation depends on die complexity of the process. However, it will be good enough for die control purpose because we usually need the complement (trade off) between system dynamics and dead time, e.g. for the process 1 (1 + sf it is reasonable to assume part of process phase shift caused by dead time. Also the estimated T and Td can be very useful for the controllers based on a first order plus dead time frame. For the use in the discrete time, a phase shift e~^^^' has to be added. We have already got the delay free controller parameters in Chapter 3, it is not difficult to adopt it to the delayed case by using the Smith Predictor where Kd represents the controller for the time-delay system and d = ^. model the plant G from our frequency domain estimation. The problem left is to Appendix B Estimation of A R M A Model from Frequency Domain Estimation In order to compare witii time domain estimation metliod or to get furtlier from our frequency domain robust estimator, we may fit a parametric frequency response model to the much large frequency domain data. By doing this we can reduce the variance of frequency domain estimation since the variance does not decay with increasing frequency point N as in E T F E case for example. The price for the variance reduction is then the risk of introducing the bias, that is the unstructured uncertainty by our definition. The method is similar to the one in Sec.3.5.3; a least square fitting can be used by assembling a set of N measured frequency response points and solving for the unknown coefficients of a transfer function model. The procedure can be best illustrated by following example. We consider our nominal model: hz-^ 1 - az-i (B.254) where: Bo = [a b]^', the parameter vector. d can be estimated from Appendix A ; z = eJ'". This first order plus dead time model is commonly used in the industry process. Using this structure, we can write: (B.255) whert: xp = [z-'^Go z-(<^+^y Since Go{Oo,z~^) is a complex number, we have: Re{Go{eo,z-''))=Re{xp'^) (B.256) lTn{Go{eo,z-'')) (B.257) 122 =Im{rp'^)eo Thus, if we know the complex value of Go{Oo,c~^'^)forsome known w, we can find two linear equations in the parameters Oo- This is why a first order nominal model can be matched perfecdy by a single sinusoid input i f the tme system is a first order system too. Usually the true system has a more complicated structure than assumed model, we can only model it in a certain frequency range with small variance. The variance in other frequencies will be relatively large. A formal description has been given in section 2.3.1. Now from the information supplied by recursive frequency function estimation of section 2.3, we can actually get the recursive nominal model estimate with the variance in a L S sense. Assume that we know the frequency response of plant at frequencies o^t for k = 1,2,3,...,N/2, we define the N x 2 matrix: Re{€-^^^ Go (00, e-^'^' )} Re{e-i'^^/' Go (Oo, e"-''^'^/^ )} ^{Go{eo,e-^'^^)} = igeje'-'^'^+i)'^'} Re{e-i(<'+i} (B.258) and the N x l matrix: Re {Go (00, e-^'-^)} Re{Go{0o,e-^'^-^^)} ${Go(^o,e--"^^)} = Im{Go{0o,e-^'-^)} (B.259) llm{Goi0o,e-'^-^^)} j However, in practice we ordy have our frequency function estimate, e.g., Eq.(2.67) with which to estimate die parameters. If we use Go{€^'^) instead of Go(^o,e~-"^) in Eq.(B.258) and (B.261), we have: # = *|G'o(e^'^*)} (B.260) and: (B.262) where: 6 = à b Then a well known result is: 0 = -1 (B.263) Thus, we can compute the value of the estimated nominal model Go ((9,e-^'^)fork=1.2,3,...,N/2. If we define our unstructured uncertainty as: A G , = \GO(è,e-^*^') - G\ (B.264) the bounding for structured and unstructured uncertainty can be written as: A G = AG„ + A G , (B.265) So, using the triangle inequality, we can show for frequency Wjt. and at time index n: Go(^,e-J'^") -Gr(e''"^) < \Go(ê,e-^'^-^ - G| + |G - GrCe''^*) < AGu + A G , = AG (B.266) Here, we simply apply Eq.(B.263) to frequency estimate of Section 2.2.2, we have a w 0.90197 and b « 0.09748. The information provided in this section, such as nominal model, and pure time delay, will be usefiil in the time-delay system control by our quantitative adaptive robust controller. A s a remarie for this section, we note that from the extension of our frequency domain robust estimator, our method can provide full information that no other methodology can. However, the price is the extensive frequency domain calculations.
- Library Home /
- Search Collections /
- Open Collections /
- Browse Collections /
- UBC Theses and Dissertations /
- Quantitative adaptive robust control based on a robust...
Open Collections
UBC Theses and Dissertations
Featured Collection
UBC Theses and Dissertations
Quantitative adaptive robust control based on a robust frequency domain estimator Wang, Hongtao 1992-12-31
pdf
Page Metadata
Item Metadata
Title | Quantitative adaptive robust control based on a robust frequency domain estimator |
Creator |
Wang, Hongtao |
Date | 1992 |
Date Issued | 2008-12-18T21:17:11Z |
Description | It is widely appreciated that robust control and adaptive control focus on dealing with system uncertainty but emphasize different aspects. Adaptive control is efficient in eliminating structured uncertainty but has difficulty handling unstructured uncertainty. On the other hand, robust control guarantees the stability for the worst case but ends up with a very conservative design. This thesis will be mainly concerned with the development of a frequency domain quantitative adaptive robust control strategy, which can be viewed as a good combination of adaptive control and robust control to obtain the best closed-loop performance of a stable closed-loop system. The research will consist of two parts. Part one is to develop a new estimation method called frequency domain robust estimator which can improve the frequency function estimate and quantify system uncertainties in the frequency domain Part Two is to implement the quantitative adaptive robust controller based on the information obtained from the robust estimator. |
Extent | 5008389 bytes |
Genre |
Thesis/Dissertation |
Type |
Text |
File Format | application/pdf |
Language | eng |
Collection |
Retrospective Theses and Dissertations, 1919-2007 |
Series | UBC Retrospective Theses Digitization Project |
Date Available | 2008-12-18 |
Provider | Vancouver : University of British Columbia Library |
Rights | For non-commercial purposes only, such as research, private study and education. Additional conditions apply, see Terms of Use https://open.library.ubc.ca/terms_of_use. |
DOI | 10.14288/1.0064784 |
URI | http://hdl.handle.net/2429/3136 |
Degree |
Doctor of Philosophy - PhD |
Program |
Electrical and Computer Engineering |
Affiliation |
Applied Science, Faculty of Electrical and Computer Engineering, Department of |
Degree Grantor | University of British Columbia |
Graduation Date | 1992-11 |
Campus |
UBCV |
Scholarly Level | Graduate |
Aggregated Source Repository | DSpace |
Download
- Media
- ubc_1992_fall_wang_hongtao.pdf [ 4.78MB ]
- Metadata
- JSON: 1.0064784.json
- JSON-LD: 1.0064784+ld.json
- RDF/XML (Pretty): 1.0064784.xml
- RDF/JSON: 1.0064784+rdf.json
- Turtle: 1.0064784+rdf-turtle.txt
- N-Triples: 1.0064784+rdf-ntriples.txt
- Original Record: 1.0064784 +original-record.json
- Full Text
- 1.0064784.txt
- Citation
- 1.0064784.ris
Full Text
Cite
Citation Scheme:
Usage Statistics
Country | Views | Downloads |
---|---|---|
China | 16 | 29 |
Russia | 7 | 1 |
United States | 6 | 1 |
Luxembourg | 3 | 0 |
Poland | 3 | 0 |
Canada | 3 | 0 |
Japan | 2 | 0 |
India | 1 | 0 |
France | 1 | 0 |
City | Views | Downloads |
---|---|---|
Beijing | 13 | 18 |
Unknown | 7 | 18 |
Saint Petersburg | 5 | 0 |
Vancouver | 3 | 0 |
Strassen | 3 | 0 |
Shenzhen | 3 | 9 |
Ashburn | 3 | 0 |
Tokyo | 2 | 0 |
Washington | 2 | 0 |
University Park | 1 | 0 |
{[{ mDataHeader[type] }]} | {[{ month[type] }]} | {[{ tData[type] }]} |
Share
Embed
Customize your widget with the following options, then copy and paste the code below into the HTML
of your page to embed this item in your website.
<div id="ubcOpenCollectionsWidgetDisplay">
<script id="ubcOpenCollectionsWidget"
src="{[{embed.src}]}"
data-item="{[{embed.item}]}"
data-collection="{[{embed.collection}]}"
data-metadata="{[{embed.showMetadata}]}"
data-width="{[{embed.width}]}"
async >
</script>
</div>
Our image viewer uses the IIIF 2.0 standard.
To load this item in other compatible viewers, use this url:
http://iiif.library.ubc.ca/presentation/dsp.831.1-0064784/manifest