LINEAR SYSTEMS IDENTIFICATION AND OPTIMIZATION WITH APPLICATION TO ADAPTIVE CONTROL ADEL ABDEL-RAOUF HANAFY B.Sc, U n i v e r s i t y o f C a i r o , C a i r o , Egypt, 1965 A THESIS SUBMITTED IN PARTIAL FULFILMENT OF THE REQUIREMENTS FOR THE DEGREE OF DOCTOR OF PHILOSOPHY i n the Department of Electrical Engineering We a c c e p t t h i s t h e s i s as conforming t o the required THE standards UNIVERSITY OF BRITISH COLUMBIA J u l y , 1972 In p r e s e n t i n g t h i s t h e s i s in p a r t i a l f u l f i l m e n t o f the r e q u i r e m e n t s an advanced degree at the U n i v e r s i t y of B r i t i s h C o l u m b i a , the L i b r a r y shal1 make i t f r e e l y available for I f u r t h e r agree t h a t p e r m i s s i o n f o r e x t e n s i v e I agree r e f e r e n c e and copying of t h i s representatives. It is understood that copying or thesis t h e s i s f o r f i n a n c i a l g a i n s h a l l not be a l l o w e d w i t h o u t my written permission. The U n i v e r s i t y o f B r i t i s h Columbia Vancouver 8, Canada Date 26-1-7 Z or publication of t h i s Depa rtment that study. f o r s c h o l a r l y purposes may be g r a n t e d by the Head of my Department by h i s for ABSTRACT T h i s t h e s i s i s concerned w i t h the problem of i d e n t i f y i n g c o n t r o l l i n g l i n e a r continuous for systems. Algorithms r e a l - t i m e d i g i t a l computations are developed the c o n t r o l l e r and The the which a r e feasible f o r the d e s i g n o f Due g e n e r a l i z e d e q u a t i o n e r r o r i s shown to be a p p l i c a b l e t o a to the imposed s t r u c t u r e of the e s t i m a t o r the of h i g h o r d e r m a t r i c e s i s avoided. Examples i l l u s t r a t e ness of the e s t i m a t o r f o r a v a r i e t y of cases n o i s e as w e l l as measurements n o i s e . the effective- dealing with quantization In some c a s e s , this d i g i t a l identi- matrix. a l g o r i t h m f o r computing the g e n e r a l i z e d i n v e r s e of a m a t r i x i s developed. T h i s a l g o r i t h m e l i m i n a t e s the need f o r Gram-Schmidt g o n a l i z a t i o n and i t s a s s o c i a t e d ( i n the i n t e r e s t of a c c u r a c y ) g o n a l i z a t i o n as new v e c t o r s are time-varying ortho- reortho- introduced. A two-stage e s t i m a t o r i s developed and para- manipulation f i e r r e q u i r e s the computation of the g e n e r a l i z e d i n v e r s e of a A simple both identifier. mean-square method of r a p i d d i g i t a l e s t i m a t i o n of l i n e a r system meters. and for estimating time-invariant parameters i n l i n e a r systems. During the second stage, the t i m e - v a r y i n g parameters are c o n s i d e r e d as unknown c o n t r o l i n p u t s to a l i n e a r subsystem of known dynamics. identifier F o r the f i r s t s t a g e , the digital i s shown t o be e f f e c t i v e i n i d e n t i f y i n g the t i m e - i n v a r i a n t parameters. Numerous examples i l l u s t r a t e the e f f e c t i v e n e s s o f the method. To d e s i g n a feedback c o n t r o l l e r a method of s u c c e s s i v e a p p r o x i mations f o r s o l v i n g the two constant gain matrices p o i n t boundary v a l u e problem f o r optimum i s developed. The method i s shown to be t i o n a l l y e q u i v a l e n t to a d e f l e c t e d g r a d i e n t method. ii Convergence comptacan always be a c h i e v e d by c h o i c e o f a s c a l a r s t e p - s i z e parameter. l i n e approximate method i s developed whose parameters must be i d e n t i f i e d . An on- which appears s u i t a b l e f o r systems The two p o i n t boundary v a l u e problem i s r e p l a c e d by an a l g e b r a i c problem, the s o l u t i o n o f which g i v e s a optimal constant gains. The sub- problem o f t r a j e c t o r y s e n s i t i v i t y r e d u c t i o n by augmented s t a t e feedback i s shown to be w e l l - p o s e d s t r u c t u r e of the g a i n m a t r i c e s i s taken. The i f constrained simple s t r u c t u r e o f constant gain matrices i s considered. Based on the developed i d e n t i f i e r and c o n t r o l l e r , one p o s s i b l e s t r a t e g y f o r o p t i m a l a d a p t i v e c o n t r o l which i s p a r t i c u l a r l y from an e n g i n e e r i n g p o i n t of view i s s t u d i e d . The attractive s t r a t e g y i s to i d e n t i f y system parameters at the end of an o b s e r v a t i o n i n t e r v a l then to use the parameters to d e r i v e an " o p t i m a l " c o n t r o l f o r a sub- sequent c o n t r o l i n t e r v a l . illustrate and Two examples are c o n s i d e r e d i n o r d e r t o the e f f e c t i v e n e s s of t h i s s t r a t e g y . iii TABLE OF CONTENTS Page ABSTRACT i i TABLE OF CONTENTS iv LIST OF TABLES vi LIST OF FIGURES v i i ACKNOWLEDGEMENT ix 1. INTRODUCTION 1 2. RAPID DIGITAL IDENTIFICATION OF LINEAR SYSTEMS 2.1 2.2 2.3 2.4 2.5 2.6 2.7 2.8 3. 8 9 14 17 18 27 30 2.7.1 2.7.2 2.7.3 30 36 37 Correlation o f Noise Samples The Error Covariance Matrix Estimation of the Confidence Interval Length Discussion 42 RECURSIVE COMPUTATION OF THE GENERALIZED INVERSE OF A MATRIX. 3.1 3.2 3.3 3.4 3.5 4. Introduction Problem Formulation Recursive and Block Data Processing Measurement Requirements Applications Comparison with a Limited Memory F i l t e r Some S t a t i s t i c a l Properties of the D i g i t a l I d e n t i f i e r . . . Introduction Properties of the Pseudo-Inverse Recursive Evaluation of R^ Recursive Evaluation of A Examples 44 44 46 49 50 1 TWO-STAGE ESTIMATION OF TIME-INVARIANT AND TIME-VARYING PARAMETERS IN LINEAR SYSTEMS 4.1 4.2 4.3 4.4 4.5 4.6 4.7 Introduction.. Problem Formulation Estimation of A . Algorithm 1 Estimation of a P e r i o d i c \ Algorithm I I Estimation of a P e r i o d i c A . Algorithm I I I Two-Stage Estimation of Af and A Examples 56 57 59 62 65 67 69 4.7.1 4.7.2 70 71 v v v Example 1 Example 2 Page 4.7.3 4.7.4 4.7.5 4.7.6 5. 71 72 73 73 THE Introduction Problem F o r m u l a t i o n U n i f i e d Treatment o f D i f f e r e n t F o r m u l a t i o n s G r a d i e n t Descent and S u c c e s s i v e Approximations Optimum Constant G a i n : T r a j e c t o r y S e n s i t i v i t y Reduction. On-Line Computation o f Gains Examples 76 78 82 84 88 90 93 OPTIMAL ADAPTIVE CONTROL 6.1 6.2 6.3 Introduction Optimum C o n t r o l D u r i n g the F i r s t Optimization Intervals Examples 6.3.2 6.3.3 6.3.4 7. 3 4 5 6 OPTIMUM CONSTANT FEEDBACK GAINS FOR LINEAR SYSTEMS AND PROBLEM OF TRAJECTORY SENSITIVITY REDUCTION 5.1 5.2 5.3 5.4 5.5 5.6 5.7 6. Example Example Example Example 102 and Subsequent 102 104 Example 1 Example 2 Discussion . CONCLUSIONS 105 107 110 116 APPENDIX I 119 APPENDIX I I .' . 122 APPENDIX III.- 125 APPENDIX IV 127 APPENDIX V 130 REFERENCES • 132 v LIST OF TABLES Pagi T a b l e (4.1) 75 T a b l e (4.2) 75 T a b l e (4.3) 75 T a b l e (4.4) 75 LIST OF FIGURES Figure Page 1.1 The Proposed General Optimal Adaptive System 2.1 Example 1. Estimation of a^ with Quantization Noise 21 2.2 Example 1. Estimation of 21 2.3 Example 1. Estimation of a^ with Quantization and Unbiased Measurement Noise 22 Example 1. Estimation of a^ with Quantization and Unbiased Measurement Noise 22 Example 1. Estimation of a-^ with Quantization and Biased Measurement Noise 23 Example 1. Estimation of a^ with Quantization and Biased Measurement Noise 23 2.7 Example 2. Estimation of a^ with Quantization noise 24 2.8 Example 2. Estimation of a^ with Quantization Noise 24 2.9 Example 2... Estimation of. b^ with. Quantization. Noise.. 25. 2.10 Example 2. Estimation of a^ with Quantization and Unbiased Measurement Noise 25 Example 2. Estimation of a^ with Quantization and Unbiased Measurement Noise 26 Example 2. Estimation of b^ with Quantization and Unbiased Measurement Noise 26 Example 2. Estimation of a., with Quantization and Unbiased Measurement Noise Using a Standard Limited Memory F i l t e r . . . 28 Example 2. Estimation of a~ with Quantization and'Unbiased' Measurement Noise Using a Standard Limited Memory F i l t e r . . . 28 Example 2. Estimation of b with Quantization and Unbiased Measurement Noise Using a Standard Limited Memory F i l t e r . . . 29 4.1 The Recursive Observation-Estimation Intervals 64 5.1 Example 1. Control as a Function of Time For: (a) Method. (b) Constant Gain and Augmented States. (c) Constant Gain and Non-Augmented States 2.4 2.5 2.6 2.11 2.12 2.13 2.14 2.15 5.2 2 with Quantization Noise 1 Lee's 98 Example 1. Trajectory Senstivity Function z^ For: (a) Lee's Method. (b) Constant Gain and Augmented States, (c) Constant Gain and Non-Augmented States 99 Figure Page Example 1. Trajectory S e n s i t i v i t y Function z.^ f o r : (a) Lee's Method. (b) Constant Gain and Augmented States. (c) Constant Gain and Non-Augmented States 100 Example 2. Block Diagram of a Sychronous Generator Voltage Control System 101 6.1-a The I n i t i a l i z a t i o n Interval and the Operation Period 103 6.1-b The Successive Observation Intervals 103 6.1-c The Successive Optimization Intervals (X, = r) 103 6.1-d The Successive Optimization Intervals (I = 2r) 103 6.2 Example 1. Noise I d e n t i f i c a t i o n of Parameter a^ with Quantization., 112 Example 1. Noise I d e n t i f i c a t i o n of Parameter a^ with Quantization 6.4 Example 1. The Optimum Piecewise Constant Feedback Gain K^... 113 6.5 Example 1. The Optimum Piecewise Constant Feedback Gain K^... 113 6.6 Example 2. Noise I d e n t i f i c a t i o n of Parameter a Example 2. Noise I d e n t i f i c a t i o n of Parameter a^ with Quantization Example 2. Noise I d e n t i f i c a t i o n of Parameter a^ with Quantization 6.9 Example 2. The Optimum Piecewise Constant Feedback Gain K^... 6.10 Example 2. The Optimum Piecewise Constant Feedback Gain 6.11 Example 2. The- Optimum Piecewise Constant Feedback Gain K,.-..- 115 5.3 5.4 6.3 6.7 6.8 112 with Quantization 114 114 114 viii 115 l^... ACKNOWLEDGEMENT I am indebted to my supervisor, Dr. E. V. Bohn, for h i s continuous guidance and assistance throughout the period of t h i s research. Thanks are due to Dr. V. G. Haussmann f o r h i s valuable comments on the d r a f t . i s duely Reading the draft by Dr. M. S. Davies and Dr. G. F. Schrack appreciated. 1 also wish to thank Miss Linda Morris for her patience i n typing the manuscript. The f i n a n c i a l support by the National Research Council of Canada, grant 67-3134 (1968-1971), the Scholarship by the B. C. Hydro Company (1970), and the U. B. C. Graduate Scholarship (1971-1972) are appreciated. ix 1 1. INTRODUCTION An optimal adaptive control i s a control that w i l l cause a plant to operate i n an "optimum" manner without a P r i o r i knowledge about the constraints on e i t h e r the magnitude or the rate of v a r i a t i o n of the plant uncertainties. These uncertainties may be i n the plant parameters, i n i t i a l conditions, and/or noisy measurements. One possible mathematical design approach for an adaptive control i s to use s t a t i s t i c a l methods. Plant uncertainties are then described by p r o b a b i l i t y density Such s t a t i s t i c a l methods are generally not computationally functions. f e a s i b l e even when a priori knowledge of the s t a t i s t i c s of the plant uncertainties e x i s t . In most p r a c t i c a l s i t u a t i o n s , s i m p l i f y i n g assumptions must be and these detract considerably from the usefulness introduced of a mathematically precise s t a t i s t i c a l approach. An a l t e r n a t i v e approach i s the optimally sensitive control [1] where no a p r i o r i knowledge of s t a t i s t i c a l parameters i s required. However, the range of plant parameter and/or i n i t i a l condition variations i s restricted. This thesis i s concerned with the development of an optimally s e n s i t i v e feedback control which does not require detailed s t a t i s t i c a l data for i t s r e a l i z a t i o n . However, i n order to increase the range of optimal control action, parameter estimation i s used. Algorithms which are f e a s i b l e for real-time d i g i t a l computations are developed for the design of both the c o n t r o l l e r and the Figure (1.1) i l l u s t r a t e s identifier. the general structure of the optimal adaptive control system which has been developed i n t h i s t h e s i s . The -n LEAST Sq. I.C ESTIMATOR A- NOISE 11 PLANT OUTPUT DIGITAL IDENTIFIER FEEDBACK Y- - 1 SENSITIVITY L_ M?Z?£Z. ^ GAINS ^ ADJUSTABLE C = 3 MODEL ~—7> /7/.7F/P 11 r 11 'L -SCOMPOSITE T ^JDENTIFI^ _p FIG. (I.I) THE PROPOSED GENERAL OPTIMAL ADAPTIVE SYSTEM, 3 blocks drawn i n dashed l i n e s represent elements which may not be required i n many p r a c t i c a l a p p l i c a t i o n s . generally be divided into This optimal adaptive control system can two main parts; the estimator and the c o n t r o l l e r . A d e t a i l e d discussion of the components which comprise these two parts together with the thesis layout i s considered In the following two sections. 1.1 The Estimator The estimator i s used to obtain a "best" estimate for the unknown plant parameters and i n i t i a l conditions. I t i s also used to obtain a "best" estimate f o r the noisy and/or the unmeasured plant state v a r i a b l e s . These estimates are needed i n order to construct a reasonably accurate deterministic l i n e a r model of the plant. A deterministic l i n e a r model i s necessary i f use i s to be made of the powerful tools that have been developed from optimal control theory. this thesis. Linear continuous models are considered i n The choice of l i n e a r models i s made because of t h e i r r e l a t i v e s i m p l i c i t y from the instrumentation and computation point of view. This s i m p l i c i t y i s largely due to the high degree of development of the theory of l i n e a r systems. Even when the plant i s nonlinear, however, the use of q u a s i l i n e a r i z a t i o n techniques can result i n a l i n e a r i z e d model of the plant. Over a fixed period of time-, such, a linear... model is. of ten. adequate. to. describe the dynamic behaviour of a nonlinear plant. The estimator i s composed of four major elements (see Figure (1.1)). These four elements are: a- I n i t i a l Conditions Estimator This estimator i s a standard l i n e a r least-squares estimator which minimizes 4 a quadratic error function based on a system of l i n e a r algebraic equations. It i s necessary to obtain a good estimate of the plant i n i t i a l conditions since the values of the c o n t r o l l e r constant feedback gains depend on i n i t i a l conditions. This estimator i s used whenever the plant states are conta- minated by noise and/or some of these states are unmeasurable. A linear least-squares estimator has been chosen due to i t s s i m p l i c i t y . Also, no a p r i o r i knowledge of noise s t a t i s t i c s i s needed. b - Filter The f i l t e r gives a continuous estimation of the plant states. I f some of the states are unmeasurable then a Luenberger observer i s a good choice to reconstruct the unavailable states. On the other hand a Kalman f i l t e r can be used for best estimation of the noisy states. c - Digital Identifier A continuous-time parameter i d e n t i f i c a t i o n method that has been extensively investigated i s the parameter-tracking technique on a steepest descent parameter adjustment. [2,3,4], These This i s based parameter-tracking i d e n t i f i e r s have been largely developed for use with analog computers. P r a c t i c a l on-line i d e n t i f i c a t i o n , however, requires the use of d i g i t a l computers and the d i s c r e t e version of parameter-tracking several d i f f i c u l t i e s . introduces The optimum gain or optimum step-size depends on the input signals and i s d i f f i c u l t to determine on-line [5]. In Chapter 2 a d i g i t a l i d e n t i f i e r i s developed which gives rapid and accurate estimation of the unknown, time-invariant parameters of the plant. The i d e n t i f i e r uses data obtained from augmented input controls and output states which r e s u l t s i n p a r a l l e l data processing. 5 A d i s c r e t i z e d version of the continuous augmented model i s used. simulation r e s u l t s have shown that the i d e n t i f i e r developed has very good rates of convergence. Digital i n this thesis This d i g i t a l i d e n t i f i e r requires the computation of a matrix inverse which may not exist i n some cases. In such cases a generalized-inverse of a matrix i s required. A number 7 of methods have been proposed f o r finding the generalized inverse of a matrix [6, 7, 8], These methods generally consist of some combination of a Gram-Schmidt orthogonalization process, the expansion of products of matrices, and the solution of sets of l i n e a r algebraic equations. In Chapter 3 an algorithm i s developed which replaces these operations and computes- the generalized inverse from simple recursive formulas. Inherent i n the method i s a simple check on the l i n e a r inde- pendency of the column vectors of the matrix. d - Composite I d e n t i f i e r A two-stage estimator i s developed i n Chapter 4. In the f i r s t stage the d i g i t a l i d e n t i f i e r i s used to obtain an estimate f o r the time-invariant parameters. In the second stage the unknown time-varying parameters are considered as c o n t r o l inputs and a l i n e a r state regulator quadratic cost function dynamic optimization problem i s formulated. associated two-point The solution of the boundary-value problem f o r the optimum control" r e s u l t s i n an estimate for the time-varying parameters. This two-stage i d e n t i f i c a t i o n algorithm i s s u i t a b l e for the cases where only a few of the unknown parameters of the plant are time-varying. complicated techniques developed The use of the for the i d e n t i f i c a t i o n of time-varying parameters i s not j u s t i f i e d f o r such cases. The usual approach taken to i d e n t i f y time-varying parameters i s to consider them as unknown augmented states described by d i f f e r e n t i a l equations of known form and then use some form of state estimation [9]. transforms have been used [10] Alternatively, integral to extend the c a p a b i l i t y of parameter- tracking system to track time-varying parameters. The implementation d i f f i c u l t i e s for such techniques make t h e i r use questionable i f not completely u n p r a c t i c a l . 1.2 The Controller In Chapter 5 the design of a c o n t r o l l e r with optimum constant feedback gains i s discussed. A method of successive approximations i s developed and i t i s shown to be equivalent to a deflected gradient descent approach with a variable weighting matrix. This property i s u t i l i z e d to develop an algorithm f o r o f f - l i n e design studies which r e s u l t s i n rapid convergence. An approximate method i s developed to replace the nonlinear two-point boundary-value problem by an algebraic optimization problem. This method appears to be suitable for on-line computation of suboptimum constant feedback gains (see e i n Figure (1.1)). Trajectory s e n s i t i v i t y function- minimization i s an important requirement... The i n troduction of the s e n s i t i v i t y functions i n the feedback structure of the adaptive c o n t r o l l e r makes i t possible to up-date the constant feedback gains less frequently and makes the plant c o n t r o l l e r less s e n s i t i v e to modelling e r r o r s . The s e n s i t i v i t y functions are generated by the sensi- t i v i t y model shown by f i n Figure (1.1). The parameters of the s e n s i t i v i t y model are up-dated each time a new estimate of the unknown plant parameters i s made. In Chapter 6 a recursive observation-optimization strategy i s described. During a f i n i t e observation I n t e r v a l , measurements of the state are taken and used to determine estimates of the plant parameters. The l i n e a r model parameters are adjusted according to these new estimates. The l i n e a r model i s used to compute a control strategy which i s applied i n a subsequent control i n t e r v a l . This scheme seems p a r t i c u l a r l y a t t r a c t i v e for adaptive systems. In summary, the objective of this thesis i s to investigate p r a c t i c a l methods f o r designing an optimal adaptive control system based on a reasonably large domain of plant uncertainties. Emphasis i s placed on techniques which can be implemented on real-time d i g i t a l process control computers. To reduce computational requirements, the control strategy should not require frequent up-dating. This objective has been accomplished by developing a rapid d i g i t a l i d e n t i f i e r f o r continuoustime systems and by developing simple algorithms for the computation optimally s e n s i t i v e feedback gains. of 8 2. RAPID DIGITAL IDENTIFICATION OF LINEAR SYSTEMS 2.1 Introduction D i g i t a l control of adaptive systems i s f e a s i b l e only if parameter i d e n t i f i c a t i o n and evaluation of control laws can be p e r formed i n r e a l time. To reduce the complexity of the computations, s i m p l i f y i n g assumptions must be made. One f e a s i b l e method i s to identify the system as an approximate l i n e a r system within a given time i n t e r v a l and then use a l i n e a r optimal control [5,11], For such a control policy to be e f f e c t i v e i n an adaptive system the i d e n t i f i c a t i o n should be rapid. A continuous-time parameter i d e n t i f i c a t i o n method that has been extensively investigated i s the parameter tracking technique [2,3,4]. This i s based on a steepest descent parameter adjustment of a quadratic error function or a Liapunov function. The error can be taken to be a response or a generalized equation e r r o r . The generalized error approach i s attractive since i t results i n an exact continuous-time steepest descent adjustment. Furthermore, by use of state variable f i l t e r s , i t allows considerable f l e x i b i l i t y i n the generation of additional input-output additional- equation e r r o r s . signals and These signals play an important:- role- in the estimator that w i l l be discussed i n this Chapter. They are considered to be the time response of augmented input controls and output states. These states can be chosen so that the quadratic error function has a unique minimum i n parameter space which then allows very rapid i d e n t i fication. P r a c t i c a l o n - l i n e i d e n t i f i c a t i o n , however, requires the use of d i g i t a l computers and the discrete version of parameter tracking introduces several d i f f i c u l t i e s . 9 The optimum gain or step-size depends on the input signals and i s d i f f i c u l t to determine on-line [5]. Furthermore, a steepest descent search f o r a minimum does not generally r e s u l t i n rapid i d e n t i fication. A popular method f o r on-line i d e n t i f i c a t i o n by d i g i t a l computers i s the mean-square technique. In the l i n e a r system case a l l state variables must be measured and the estimated parameters are biased. Furthermore, i f rapid i d e n t i f i c a t i o n i s attempted, a r i s e with i l l - c o n d i t i o n e d matrices. d i f f i c u l t i e s can A d i g i t a l estimation method i s developed here which can overcome these d i f f i c u l t i e s by use of augmented input controls and output states. Due to the imposed structure of the estimator the manipulation of high-order matrices, which i s generally the case i n the mean-square approach, i s avoided. 2.2 Problem Formulation Consider a completely observable and controllable continuous- time l i n e a r system x = Ax + Bu + w (2.1) where x i s an n-state vector, u i s an m-control vector and w i s additive plant noise. Let T be the sampling i n t e r v a l and l e t t , t ^ , t^, a set of sampling instants. Q be For the moment, to ease the n o t a t i o n a l burden, l e t z^(k) = z^(t^) be the sampled measurement at t ^ where z ^ t ) = x(t) + v(t) and where v ( t ) i s additive measurement noise. (2.2) (The estimator to be developed does not require that a l l states be measured as implied by (2.2). The general measurement requirements are discussed i n Section (2.4) Let Z i+l - ( t ) f o\^ '> d U £ l ( t ) i \ ™ * > = + W = l \ ^ fc) d ^ <'> 2 3 The vectors z^, u^, (£ = 2), are defined as augmented states and augmented controls, respectively. I t i s seen from (2.1), (2.2) and (2.3) that z„ = Az„ + Bu„ + w. % I I I A A • where u^ = u, w^ = t = x = t k k + 1 (2.4) v - Av + w. Integrating (2.4) over the i n t e r v a l and l e t t i n g Z (k) = z ( t > , U (k) - u ( t ) , y i e l d s Z (k+1) = AZ^(k) + BU (k) + N (k) £ £ k £ £ £ £ k (2.5) £ where A = I + TA, B = TB and V k ) " [k t ( < ) k+1 5 z T " £< z £ t ) ) + B ( u k Jl ( T ) " £ u ( t k } ) + V T ) ] d T (2 ' 6) The discrete set of equations forms the basis f o r a l l further discussion. To ease the notational burden a l l n(n+m) elements of A and B are considered unknown. The problem i s to estimate A and B from a set of measured Z^(j)» U-j^(j) and computed (from (2.3)) augmented states and controls. Let z (k) 4 A V - - (2...Z) c = k) b' n where prime denotes transpose and where a^ and b^ are the j - t h n and m row vectors of A and B respectively. These rows are ordered into an 2 (n + nm) vector c as shown. Equation (2.5) can then be written as z (k+i) = e (k)c + N (k) £ £ 4 (2.8) where Yj(k) 0 0 ( Y l ( k ) ) (2.9) Bn •Y k) i( (Due to the frequent occurrence of matrices of block-diagonal form i t i s convenient to introduce a s p e c i a l notation. The f i r s t subscript B i s used to indicate that a block diagonal matrix i s to be formed from the vector or matrix within the bracket. The second subscript n indicates 2 the number of blocks. Note that 9 (k) i s an nx(n + nm) matrix). JO Consider the sampling instants k-N+1, k, so that N blocks of data are available and l e t N(k-l) Z(k) Z £ ^N T (2.10) c + N(k-N) Z (k-N+1) J where z (k) e(k-i) x K N , 6(k) = • mm " (2.11) , Z(k) e(k-N) n+m n+m< > Z k The matrix 9(k) i s n(n+m) x n(n+m) and the vector Z(k) i s n(n+m) x 1. The matrix 6^ i s of dimension nN(n+m) x n(n+m) as can be seen from (2.11). It i s seen from (2.7), (2.9) and (2.11) that the measured and computed augmented states and controls are used to form the matrix 8(k) and the vector Z(k) at the k-th sampling instant. These quantitites are evaluated for N sampling instants and are used to form the matrix 9^ and the vector Z^. The problem i s to determine an estimate of the parameter vector c given 9^ and Z^, where these quantities are related by (2.10). The simplest approach i s to define the optimum estimate to be the estimate which minimizes the quadratic form (Z^ - 8 ^ 0 ) ' ( Z ^ - Q^c). The solution i s the minimum mean-square (also c a l l e d least-squares) estimate [6] c = P 6*Z. N N N M (2.12) T where h 1 = <w (2.13) Since the basic structure of a l l l i n e a r mean-square estimators i s the same i t i s of i n t e r e s t to compare the estimator (2.12) with previously proposed estimators [5, 12-16]. approach. Mayne has used the Kalman f i l t e r This requires that (2.10) be interpreted as a measurement equation with 6^ as a deterministic measurement matrix. This i s the case only i f there i s no measurement noise and i f a l l states are measured. A further and very s i g n i f i c a n t difference i s i n the d e f i n i t i o n of 0^ and N(k). These i n v e s t i g a t o r s take 9 (k-l) N (k-1) x 1 , 0 00 N (k-n-m) 1 (References case. (see [2, 13-15], furthermore, (2.14) N(k) 0 (k-n-m) treat only the scalar measurement Their r e s u l t s do not d i r e c t l y apply to the above case). Since (2.5)) E [ Z ( j + l ) N ^ ( j ) ] = E [ N ( j ) N j ( j ) ] •* 0 1 1 i t follows that the estimate i s biased. viously [13]. (2.15) L This has been pointed out pre- In order to eliminate the bias, s u f f i c i e n t time must elapse between samples. This can make rapid i d e n t i f i c a t i o n d i f f i c u l t . Even when the bias i s small enough to give acceptable i d e n t i f i c a t i o n , further d i f f i c u l t i e s can a r i s e with data i n the form [2.14]. I f a small sampling i n t e r v a l i s chosen to achieve rapid i d e n t i f i c a t i o n there may not be a s i g n i f i c a n t difference between adjacent rows of 6(k) r e s u l t i n g i n an i l l - c o n d i t i o n e d matrix. These d i f f i c u l t i e s tend to be avoided i f data i n the form (2.10) i s used. Some insight into the estimation bias associated with (2.12) can be obtained by considering the noise samples N ^ ( j ) , N^(k), (p = 1, ..., n + m; j ^ k) , as independent. I t follows from this assumption and (2.5) that E[Z (k)N'(k)] = 0, (Z, 0 Consequently 6.. at sampling _L (2.16) time t, i s independent of N (k) and the estimate ic p (2.12) i s unbiased when N = 1. noise i s an important p = 1, .... n + m) From (2.6) i t i s seen that quantization component of N^(k) making an analytic i n v e s t i g a t i o n into the dependency of the noise samples i n t r a c t a b l e . observations can, however, be made. Some general Successive integrations performed on w (T) not only reduce the e f f e c t i v e amplitude of the noise but also reduces the c o r r e l a t i o n between noise samples N^(k), (£ ~ 1 , ..., n + m). The truth of this statement for a certain class of noise samples i s shown i n Section (2.7). Consequently, for a given k, the quantization noise terms of (2.6) for d i f f e r e n t I tend to be uncorrelated. Con- sequently, l i t t l e i f any d i f f i c u l t y i s to be expected with estimation bias. A l l examples t r i e d v e r i f i e d this expectation. The e s s e n t i a l difference between the two estimators can be seen from (2.15) and (2.16). If £ = 1, sequential data processing must be used. I f T i s chosen large to overcome i l l - c o n d i t i o n i n g then quantization noise can become excessive. The combination of heavy noise and estimation bias (2.15) can r e s u l t i n f a i l u r e of the conventional estimator to give acceptable estimation. These d i f f i c u l t i e s tend to be avoided i f augmented states and augmented controls are chosen (£ = n + m) so that p a r a l l e l data processing can be used. H e u r i s t i c engineering j u s t i f i c a t i o n s can also be given for the use of augmented states and controls i n l i n e a r system parameter i d e n t i f i cation. I f combined parameter estimation and optimal control i s attempted signal-to-noise problems can a r i s e [12], Optimal control attempts to keep state deviations from a zero reference state small which decreases the signal-to-noise r a t i o . signal-to-noise r a t i o . Integrators are low-pass f i l t e r s which improve The data used i s then less noisy. Such a s i t u a - t i o n has been reported by Bakke and he suggests using a hybrid i d e n t i f i cation method [16]. B a s i c a l l y h i s hybrid method i s to use data i n the form (2.14) and then replace Y^, 6^ by X^ y 82* Two d i f f e r e n t are taken and a decision i s then, made on the best estimate. estimates However, a f u l l set of augmented states i s not considered, nor i s the data combined i n a mean-square sense. 2.3 Recursive and Block Data Processing The estimator (2.12) operates i n a block data processing mode. A new estimate i s made a f t e r every N a d d i t i o n a l data blocks. Other processing modes can be used such as the recursive method suggested by K i s h i where old data blocks are deleted and new data blocks introduced. - Let, N-R 9 - .> ' N 9 N-R N-R ZZ z- X T= N R (2.17) zR It follows from (2.13) that P N' =N-R R R P + 9 9 (2.18) Suppose that a l l quantities i n (2.12) are known and i t i s desired to determine ^ from these quantities (This amount to determining the best estimate when R data blocks are deleted). Making use of ( 2 . 1 8 ) i t follows from ( 2 . 1 2 ) that N-R N-R - ^ 9 " W Z -N-R " ( 2 1 9 ) By p a r t i t i o n i n g ( 2 . 1 2 ) i n the manner indicated by ( 2 . 1 7 ) , i t i s seen that °N = W R N N-R N-R + P 6 Z ( 2 ' 2 0 ) substituting ( 2 . 1 9 ) into ( 2 . 2 0 ) y i e l d s C VR" ( I " V R V ^ N V R V " (2 ' > 21 With the aid of the i d e n t i t i e s < "N R R " N R =N R X P < 6 6 ) ~ VRV" X lp 1 e = P 1 + 6 ( " I ~N R R ( I P e 0 W i ^ ) _ l p N R R e e ( 2 ' 2 2 ) ( 2 . 2 1 ) can be reduced to the standard from for recursive estimation C VR = *N N R + P 9 ( I " R N R e P 9 ) _ 1 ( V " R> Z N ( 2 ' 2 3 ) applying the matrix inversion lemma to ( 2 . 1 8 ) y i e l d s N-R = N P P + Vi ( I " W R ^ V N ( 2 ' 2 4 ) Equation ( 2 . 2 3 ) gives the best estimate i f the l a s t R data blocks are deleted and ( 2 . 2 4 ) gives the associated P-matrix. If R new data blocks are sampled, a new Z and new values of c^ and P^ can be computed. and 6 D are available The equations giving these values are obtainable from ( 2 . 2 3 ) and ( 2 . 2 4 ) by interchanging N with N - R, replacing 6_ by K -8_ and leaving 0' unchanged. K is. With these changes ( 2 . 2 3 ) and ( 2 . 2 4 ) gives the Kalman f i l t e r algorithm (since the measurement matrix i s a random variable the Kalman f i l t e r i s not optimal). proposed i n References The estimators [ 5 , 1 2 - 1 5 ] are based on a scalar measurement. The bracketed quantity i n (2.24) i s then a scalar and no matrix inversion i s required. However, many scalar samples are required to estimate the parameters within a given confidence interval. This can r e s u l t i n slow identification. The matrix to be inverted i n (2.24) i s of dimension Rn(n + x Rn(n + m). m) This i s a consequence of the use of a vector measurement of augmented states. In e f f e c t , rapid i d e n t i f i c a t i o n i s attempted by p a r a l l e l processing of augmented state samples. the The p r i c e paid for p a r a l l e l processing i s i n the high order matrix to be inverted i n the Kalman f i l t e r algorithm. An a l t e r n a t i v e recursive approach i s therefore desirable. It follows from (2.9), and n+m (2.14) that k-1 SL=1 J=k-N where n+m k-1 Q = £ £ N V j ) - Y i ( j ) ( 2 ' 2 6 ) 1=1 j=k-N consequently P N = %\n ' (2 27) Equation (2.27) gives the i n t e r e s t i n g r e s u l t that only an (n + m) x (n + m) matrix need be inverted to find V^. n 2 + nm parameters, an (n 2 + nm) x (n 2 + nm) (In general, to determine matrix has to be inverted). The structure of (2.27) allows a recursive algorithm to be developed which gives iV I n 11 terms of Q"''. Consider (2.26) and N = Y„(k I - N), S, = Q_ - V.V' ' 1 N 11 T S n+m = Q N-1 1 define S. = S. - V. V'., j J-1 J j' (2.28) applying the matrix inversion lemma y i e l d s s- = 1 (Ql) - v^p- - Q J + Q" v c i + v ^ r ^ o " 1 1 1 1 1 1 (2.29) Note that the quantities i n brackets on the right hand side of (2.29) are scalars. for finding The set of equations (2.29) define a recursive algorithm given Repeating the algorithm R times allows to be expressed i n terms of (see (2.27)). P^^ This procedure replaces (2.24) i n the data deletion step. Substituting (2.18) into the f i r s t i d e n t i t y i n (2.22) y i e l d s P W N-R R = 6 " WP"' 1 <' > 2 30 which allows (2.23) to be expressed i n the form ^N-R " ^N + P N-R R VN 6 ( "V The data deletion i s performed using (2.29) and (2.31). version i s required. . ' (2 31) No matrix i n - Insertion of new data i s performed using (2.29) i n reverse and by changing (2.31) to the appropriate form. 2.4 Measurement Requirements 1 Due to the use of state variable f i l t e r s not a l l states need be measured as implied by (2.2). Most systems of the form (2.1) can be decomposed into subsystems which have a phase variable structure. Suppose, for example, that a subsystem i s represented by x x 2 x = a x 3 x 3 = X 3 + a x 2 2 + a l X ] L + ... 3 = x 2 (2.32) To i d e n t i f y (2.32) i t i s s u f f i c i e n t t o generate the complete s t a t e s of one the augmented models. The f i r s t augmented model t h a t can be g e n e r a t e d c o m p l e t e l y i s the t h i r d model (Jl = 3) . C Z ^ o 3 = ^ o f t o It i x d i T ^ 2 0 x„ d r , d x 0 1 3 Z t 1 1 d r dx„ 1 1 2 dx, l (2.33) 2 i s seen from (2.33) t h a t o n l y x^ need be measured. augmented model (l > 3) can be g e n e r a t e d u s i n g 2.5 x o •^ x o i 2 1 f o x„ d r , d x 2 o d x T h i s i s seen by t a k i n g Any a d d i t i o n a l (2.3) and (2.33). Applications Example 1. C o n s i d e r the l i n e a r system x.^ = X 2 a - | X + a x 2 2 that only simplification. X x X 4 x 6 _8 x = = 0, t ^ = 16., and take I = 4 (see x l " 3 • 5 x Q 5» 8 x = (2.34) i n some augmented system i s o f the form X 3 where x^ = x^, g 0.0 2 "2 = X , (°) two elements need be i d e n t i f i e d r e s u l t s l x 1.0 x The f u l l x + u, x^(o) l X where a^ = -2.0, a^ = -3.0, t The f a c t 1 x x 5 X 7 7' " " a l " l " U U + 2 (2.35) U 3 _4_ U Equation (2.35) i s d i s c r e t i z e d and expressed appropriate d e f i n i t i o n of Z^(k), ^ ( k ) & an c ' i n the form of (2.8) by = ( p a » The estimation of c i s performed by a block processing mode with R = 1. ( i f no old data i s deleted, the estimation algorithm i s e s s e n t i a l l y the Kalman f i l t e r algorithm). The following sample i n t e r v a l s and c o n t r o l are used: (1) T = 0.1, (2) T = 0.2, (3) T = 0.4, (4) T = 1.0 with the corresponding values of N: (1) N = 21, (2) N = 11, (3) N = 6, t (4) N = 3. The control u(t) = 0.1 + s i n 3t + s i n lOt i s taken. Figures (2.1) and (2.2) give the r e s u l t s with zero additive measurement noise (n^ = o = n ) . 2 However, d i s c r e t i z a t i o n of (2.35) introduces quantization errors which can be considered as additive correlated quantization noise. It i s seen from the figures that i d e n t i f i c a t i o n i s rapid even for large quantization noise. Figure (2.3) and (2.4) i l l u s t r a t e the r e s u l t s f o r state dependent noise of the form x,(k) n.(k) = a. • ( - i ) ng., ( i = 1, 2) m (2.36) where ng^ i s a pseudo-random computer generated noise sequence which has a maximum amplitude made. of n . The choice a. = a„ = 0.1, E(ng.) =0. i s m 1 2 l Figures (2.5) and (2.6) i l l u s t r a t e the results for the choice = 0.1, a 2 = 0.25, E (ng^) = 0.1, E (ng ) =0.14. 2 I t i s seen from these r e s u l t s that the estimator i s i n s e n s i t i v e to correlated quantization noise, state dependent measurement noise and to bias i n the measure ment noise. Example 2. Consider the l i n e a r system (See Reference [3]) X X i 2 = --2.0, where = x 2 = a l X l a 2 = , x3^(0) = 1.0 + a x + bjU , X ( 0 ) = 0.0 2 2 -1.0, b 2 1.0, and take Z = 4. 1 system i s of the form • X 2 X 4 X 6 x l x 2 x 3 x 4 " l" X 5 x 6 _ 2_ X 7 X 8 = X u 8, U • 3 4 l U 2 U 3 + b a (2.37) l ^ 4 _ U • • where x_ = x,, x a u 5 r = x,, x_ = x . 6' 7 8 n The following sampling i n t e r v a l s are used: (1) T = 0.1, (2)- T =• 1.0-. (1) N = 21, (2) N = 3. The-cor-responding- values of N are: The value of R i s taken as R = 1 . (1) The i n i t i a l parameter estimate i s c = 0 and u(t) = 0.1 + s i n t + s i n 3 t . Figures (2.7), (2.8) and (2.9) give the results f o r quantization noise only. Figures (2.10), (2.11) and (2.12) give the results f o r quantization noise and additive state dependent noise where « 1 = a 2 = 0.1, E ^ g ^ = 0 , i = 1, 2. Comparing these time responses with those obtained i n Reference [3] f o r the noise-free case, i t i s evident that i d e n t i f i c a t i o n i s more rapid. The speed of response of parameter tracking systems seems to be adversly effected by the i n t e r a c t i o n between parameter generating subsystems which i s an inherent feature of the steepest descent search procedure used. Such subsystems are not required f o r the method developed here. Consequently, i n the noise free case, one-step i d e n t i f i c a t i o n i s achieved. FIG. (2.2) EXAMPLE 1. ESTIMATION OF WITH QUANTIZATION NOISE 0.4 FIG. (2.4) EXAMPLE 1. ESTIMATION OF a WITH QUANTIZATION AND 2 UNBIASED MEASUREMENT NOISE. 0.0 2 4 t 6 8 10 12 U FIG. (2.61 EXAMPLE!. ESTIMATION OF a BIASED . MEASUREMENT NOISE^ 2 16 IS 20 WITH QUANTIZATION AND 0.4 r 1 1 1 1 1 1 2 3 4 5 ®-\ N \ V 5 0 7 9 70 ' FIG.(2.7) EXAMPLE 2. ESTIMATION OF a, WITH QUANTIZATION NOISE. -i 1 1 2 1 3 1 4 1 5 1 6 i 1 7 FIG. (2.8) EXAMPLE 2. ESTIMATION OF a 8 2 1 9 1— 10 WITH QUANTIZATION NOISE. 1.6 1.2 or* 0.4 0.0 I 2 3 10 -0.4L FIG. (2.9) EXAMPLE 2. ESTIMATION OF b. WITH QUANTIZATION NOISE. FIG. (2.10) EXAMPLE 2. ESTIMATION OF a, WITH QUANTIZATION AND UNBIASED MEASUREMENT NOISE. 04. 2.0 . -2.4 . FIG. (2.11) EXAMPLE 2. ESTIMATION OF UNBIASED MEASUREMENT NOISE. WITH QUANTIZATION AND 2.4 2.0 1-6 1.2 0.8 Or** 0.4 0.0 w FIG. (2.12) EXAMPLE 2. ESTIMATION OF b, WITH QUANTIZATION AND UNBIASED MEASUREMENT NOISE. 27 2.6 Comparison with a Limited Memory F i l t e r In this section a comparison i s made between the proposed d i g i t a l i d e n t i f i e r using augmented states (I limited memory i d e n t i f i e r ( I = 1). fier (or f i l t e r ) > 1) and the standard The standard limited memory i d e n t i - i s of the type discussed i n Section 2.3. The f i l t e r operates i n a block data processing mode where R old data blocks are deleted and R new data blocks introduced r e s u l t i n g i n a limited memory f i l t e r . (finite) Example 2 i s used to i l l u s t r a t e the superiority of the augmented state i d e n t i f i e r f o r various values of T and N (R = 1). Figures(2.13), (2.14) and (2.15) show the results for the i d e n t i f i c a t i o n of the parameters a^, a^ and b^, respectively. standard limited memory f i l t e r f a i l s I t i s evident that the to track the parameters f o r values T = 1., N = 3, R = l and exhibit poor tracking accuracy for values T = 0.1, N = 21, R = 1. The f a i l u r e of the standard limited memory f i l t e r to give acceptable estimation of the parameters i s due to the fact that i t i s not an optimal f i l t e r f o r the proposed i d e n t i f i c a t i o n problem. Consequently, excessive quantization noise and c o r r e l a t i o n between samples, as i s the case here, can have a serious e f f e c t on i t s performance. In conclusion i t can be stated that the use of augmented states and controls r e s u l t s i n a more rapid and accurate i d e n t i f i c a t i o n of continuous time l i n e a r system parameters than can be achieved by either parameter tracking systems or the standard limited memory f i l t e r . FIG. (2.13) EXAMPLE 2. ESTIMATION OF a, WITH QUANTIZATION AND UNBIASED MEASUREMENT NOISE USING A STANDARD LIMITED MEMORY FILTER. FIG (2.14) EXAMPLE 2. ESTIMATION OF o WITH QUANTIZATION AND UNBIASED MEASUREMENT NOISE USING A STANDARD LIMITED MEMORY FILTER. 2 29 28. 2.4 . 2.0. 16. FIG. (2.15) EXAMPLE 2. ESTIMATION OF b, WITH QUANTIZATION AND UNBIASED MEASUREMENT NOISE USING A STANDARD LIMITED MEMORY FILTER. 30 2.7 Some S t a t i s t i c a l Properties of the D i g i t a l I d e n t i f i e r In the previous section, the performance of the d i g i t a l i d e n t i f i e r was evaluated by comparing estimated values with known values. a p r a c t i c a l s i t u a t i o n , however, the parameters are not known and means must be used to establish confidence i n the estimates. In statistical In reference [5] K i s h i made use of extensive studies made by Linnik [17] to determine confidence i n t e r v a l s f o r parameter estimates i n the case of a simple single-input single-output system. The relations given by Linnik w i l l be used i n this section to obtain approximate confidence l i m i t s for the parameter estimates of the d i g i t a l i d e n t i f i e r . However, Linnik's results can be applied only i f the d i f f e r e n t components of the noise vector N (see (2.10)) can be represented by a s t a t i s t i c a l noise model. task i s not an easy one. Such a However, i t w i l l be shown that under certain simplifying assumptions i t can be argued that the d i f f e r e n t noise components w i l l be uncorrelated. 2.7.1 C o r r e l a t i o n of Noise Samples In Section (2.2), i t has been assumed that the noise samples of the d i f f e r e n t models are uncorrelated. to be true. Hence, equation (2.16) i s assumed In the following, i t w i l l be shown that i f the noise samples are uncorrelated then the samples of the integrated noise are also uncorrelated. For s i m p l i c i t y , we w i l l consider a single noise component in order to prove the v a l i d i t y of the previous statement. Consider the noise component n(t) which can be represented by the following s t a t i s t i c a l model E[n(t)] = 0 (2.38) E[n(t)n(x)] = a e - 0 1 !^' If a i s of i n f i n i t e value then n(t) w i l l be a type of white-noise component. (In the white-noise case the amplitude a w i l l be equal to i n - finity) . For a p h y s i c a l r e a l i z a t i o n of the noise component n(t) i t w i l l be assumed that a i s a very large but f i n i t e quantity and i s a very large but f i n i t e quantity and that a i s of f i n i t e that a' value. Let N (t) x I t follows from (2.38) and = n(T )dx 1 (2.39) 1 (2.39) that E[N (t)n( )] = E[n(x)n(T )]dx 1 o 1 1 1 T 1 = a / o To evaluate (i) e- l - l' a T dx, 1 T the above i n t e g r a l we w i l l consider two (2.40) cases: Case 1, x ^ t E[N (t)n(x)] 1 = a £ = (ii) t 1 Case 2, a a [e e~ a ( T -a(x-t) _ T l _ dx^ ) e -ax o sc x «_ t / \ i r E[N,(t)n(x)J = a J 1 o T-TUT T e -a(x-x.) , 1 dr, 1 , ,t -a(x..-x) , + a f e 1 dx, x 1 a -ax ax , a ax -ax -at, = — e [e — 1J H e [e -e J a a r a = — [ a r o 2 r -ax - e -a(t-x), - e ] , (2.42) 0 32 For fixed time t and a very large value of a, equations (2.41) and (2.42) can be approximated by E [ N f=0 , x > t I (t)n(x)] J ' " ** U ^ , 0 < x < l V W n (2.43) J As seen from (2.43), the noise components n(t) and N^(t) correlated become un- as a -*- . 00 Using equation (2.41) or (2.42) for the s p e c i a l case where x = t yields E[N (t)n(t)] 1 = f [1 = — The following e" ] at , for fixed t and a ->• °° non-dimensional r a t i o can be used as a measure of the degree of c o r r e l a t i o n of N^(t) E. (t)q(.t)] IMr and n(t): r = - r , for fixed t, a -> « (2.44) at Equation (2.44) shows that f o r fixed t and a large value of a the between N^(t) correlation and n(t) samples i s much less than the c o r r e l a t i o n between the samples of n ( t ) . This shows c l e a r l y the advantage of using the grators i n order to generate the augmented states and controls. inte- It i s also i n t e r e s t i n g to study the degree of c o r r e l a t i o n between the samples of N (t). x Equation (2.39) y i e l d s E[N (t)N (x)] = 1 1 E[N (t)n(x )]dx 1 1 1 (2.45) As before, we have two cases: (i) Case 1, x > t ElN.WN.d)] = 1 1 f]2 - e a o a l T - e -ax, . « [ 2 T I l 1 a ~ a ( " T l T ) dx, ] 1 -a(x-x,) . - £ a a o a -ax , -at -a(x-t) = — j [2at + e +e - e a r o = (ii) 2at a -1] , f o r fixed t and a -*• (2.46) 00 Case 2, o < x < t E[N l (t)N (x)] = / i o E[N (x)n(x )] dx, + / l l l x T E[N (x)n(x.)] dx, l 1 1 1 x a -ax -a(x-x ) , , = J — L2 - e 1 - e 1 ] dx, + o a 1 f r o 1 1 t a -a(x,-x) -ax, -, , / — Le 1 - e 1] dx, x a 1 r r n ~ = —[2ax + e a a a x , ~ +e a t - e -a(t-x) -1] 2ax = , For the case x = t equation f o r fixed t and a •> 00 (2.46) or (2.47) y i e l d s E[N- (t)] = ~ a [at + e 2 = — a Equations (2.47) t a t - 1] (2.48) (2.48) and (2.38) y i e l d E[N (t)] —4 a l tE[n^(t)] 2 a (2.49) Equations (2.46), (2.47), (2.48) and (2.49) show that the c o r r e l a t i o n between the samples of the noise component N^(t) vanishes as a -»• . 00 Equation (2.49) shows that f o r fixed t and large value of a the correl a t i o n between the samples of N^(t) i s much less than that between the samples of n ( t ) . The same analysis can be d i r e c t l y extended to ^ ( t ) where N (t) 2 = N ( x ) dx 1 1 (2.50) 1 Since the same approach i s to be used i t i s s u f f i c i e n t to summarize the corresponding r e s u l t s E[N (t)N (t)] = ^ t 2 E [ N (2.51) 2 1 2 2 at ( t ) ] 3 ~ fo~ ( 2 Efr-Cty^ct)], = t E[n (t)] 5 2 ) (2.53) 9 3 " 2 E[N (t)] 2 ~ t E[n (t)] 4 2 Equations (2.51), (2.52), (2.53) and same conclusions (2.54) 1 3 a t (2.54) show that the that have been reached for N^(t) hold for N ( t ) . 2 Suc- cessive integrations of s t a t i s t i c a l l y independent noise samples r e s u l t i n s t a t i s t i c a l l y independent noise components. Also, successive i n t e - grations of correlated noise samples r e s u l t i n less correlated noise components. The above s i m p l i f i e d mathematical treatment shows that the independence property assumed by equation (2.16) i s j u s t i f i a b l e for many p r a c t i c a l situations where n(t) i s usually taken as white noise. Equation (2.16) can now be used to prove the unbiasness property of the digital identifier. Equation (2.10) has the form Z ( k + 1) = 6 ( k ) c + N(k) N (2.55) N Where the argument value denotes the sampling instant at which the various vectors and matrices are evaluated. Combining (2.55) and (2.12): E[c] = E [ P ( k ) 9 ^ ( k ) Z ( k + 1)] N N = E [ P ( k ) e ^ ( k ) ( G ( k ) c + N(k))] N N = c + E[P (k)e^(k)N(k)] (2.56) N The elements of P ( k ) and 8 (k) depend on Z (k) and U (k), A = 1, M n+m. Equation (2,55) shows that Z^(k) i s only related to N(k - 1). Since N(k) i s s t a t i s t i c a l l y independent above analysis, then ^ ( k ) and 6jj(k) a r of N(k - 1) as shown by the e not s t a t i s t i c a l l y dependent on N(k) and the following relations hold E[P (k)8^(k)N(k)] = 0 (2.57) N E[£] = c (2.58) Equation (2.58) shows that the estimate obtained by using the d i g i t a l i d e n t i f i e r i s unbiased. The results given i n figures (2.1)-(2.4) and (2.7)-(2.12) would not be obtainable i f (2.16) d i d not hold approximately. Equation (2.58) holds when 7^ = (9^9^) ^. to by K i s h i as the observable case. This case has been referred The case where the inverse (9^9^) ^ does not exist i s referred to as the non-observable case. In the non- observable case an estimate i s obtained by use of the equation c - 9* Z where (2.59) N i s the generalized inverse of 6 ^ . That i s , e^e x = X (2.60) N It follows from ( 2 . 1 6 ) 2.7.2 and ( 2 . 6 0 ) that the estimator ( 2 . 5 9 ) i s unbiased. The Error Covariance Matrix Let V = 6a N - ^ (2.61) be the error vector representing the difference between a predicted measurement vector and the measurement vector. For the observable case the above equation can be written i n the form (see ( 2 . 1 2 ) ) = V Using equation ( 2 . 5 5 ) - I ) Z N ( 2 ' 6 2 ) ( 2 ' 6 3 ) yields ( E For the non-observable <Wi N N N " P 6 I ) S case the expression f o r V i s v = (e e^ N (2.64) - I)N where use i s made of the property 6*B = N N N N fl (2.65) The error convariance i s defined by Q w = E[W»] (2.66) Making use of equation (2.63) y i e l d s Q W = E [ ( W N " I ) M ' ( N N N 9 P " E I ) ] ( 2 * 6 7 ) For the non-observable case using (2.64) y i e l d s Q 2.7.3 W = E [ ( 6 " ^ ' ( V N N N 6 " I ) ] ( ' 2 6 8 ) Estimation of the Confidence-Interval Length Only the observable case w i l l be studied i n t h i s subsection since f o r the non-observable case the length of the confidence i n t e r v a l w i l l depend on the degree of redundancy i n the measurements which, would change from one system to the other. = V N N N " ( 6 P 6 I ) Equation (2.63) y i e l d s S ( 2 > 6 3 ) The fol*lowing relations' hold f o r the observable case:" R(8 ) = R(6^) = R(6^6 ) = R(P ) = n(n + m) N N (2.69) N where R(*) stands f o r the rank of the matrix written between the brackets. Let the matrix Y„ be defined as N Y N - N N P 6 ( 2 - 7 0 ) Using the known lemma R(AB) 4 min(R(A), R(B)) (2.71) where A and B are any a r b i t r a r y matrices, i t follows from (2.70) that, R <V = n ( n + m N N ) Y N ) (2.72) Equation (2.70) y i e l d s 3 N = ( 9 9 ( 2 ' 7 3 ) Hence R(8jJ) = n(n + m) < min(R(6^9 ) , R(Y^)) N < min(n(n.+m), R(Y )) (2.74) Hence R ( Y ) > n(n + m) (2.75) N It follows from (2.72) and (2.75) that R(Y ) N Let the matrix N = R(P <3N) = < N N + M n > < ' > 2 76 be defined as "HHW - (2 77) Then W9 M M N N It follows from (2.78) that = 9.. (2.78) N i s a projection matrix. Using (2.78) and (2.71) y i e l d n(n + m) = R(8 ) N < min(R(W ) , R ( e ^ ) N < min(R(W ), n(n + m)) = N (2.79) Hence. R(W ) N ^ n(n + m) (2.80) Applying (2.71) to (2.77) i t i s seen that R(W ) N 4 n(n + m) (2.81) It follows from (2.80) and (2.81) that R(W ) N = n(n + m) (2.82) Considering now the computation of the number of degrees of freedom of the random vectors c and V, that i s (as defined by Linnik [17]) the number of s t a t i s t i c a l l y independent components of c and V. Starting by computing the degrees of freedom of the random vector c we w i l l make use of Theorem 2.3.1 given by Linnik [17] which states that i f the vector U i s expressed as U = AX (2.83) and the matrix A i s of dimension r * £ with r < Z and R(A) = r then the random vector U i s of r degrees of freedom. From equation (2.12) £ = P 9'Z N N N M (2 M Equation (2.76) y i e l d s R(P. 0') = n(n + T m) (2.84) Applying Theorem 2.3.1 given by Linnik [17] i t i s readly seen that the vector c i s a random vector of n(n + m) degrees of freedom. To compute the degrees of freedom of the noise error vector V, (2.63) i s expressed i n the form V = (W M where W - I)N (2.85) i s given by (2.77). Let 'N The matrix (2.86) i s idempotent as can be seen from the following = W N - 2W = I - N + I N (2.87) Accordingly, there exists an orthogonal matrix 0^ such that (see reference [17]) where the matrix i s a diagonal one. The matrix $^ has the property that $ The matrix N*N " V N = 1 ( 2 ' 8 9 ) i s an idempotent matrix as can be seen from the following N N N N N N N N N N =z $'w $ = n (2.90) Hence: d . = d.. 2 n ii i.e., d.. =0 or 1 n (2.91) The number of d ^ s equals to unity i s given by the rank of the matrix W^, that i s , by n(n + m). From (2.89) i t follows that R(I - W ) = R ( ^ ( I N = Rd - - W )$ ) N N y = n(n + m)(N - 1) (2.92) The f i n a l r e s u l t follows from the fact that the t o t a l number of diagonal elements of i s n(n + m)N and of these n(.n + m) are equal to unity (see (2.91)). Applying the above r e s u l t given by (2.92) to (2.85) and making use of Theorem 2.3.1 [17], i t follows that the noise error vector V i s a random vector of n(n + m)(N - 1) degrees of freedom. From the r e s u l t reached i n the previous statement i t i s seen that V'V can be expressed i n the form n(n+m)(N-1) V'V = I i=l (2.93) Let Q N = E[c c'] (2.94) cc With c i s computed according to (2.12). The following r e s u l t s given by K i s h i [5] and Linnik [17] are If £ and E,^ are s t a t i s t i c a l l y independent Gaussian random 2 2 — variables with d i s t r i b u t e d as % and having m degrees of freedom, required. then the t - d i s t r i b u t i o n i s formed by the following r a t i o t = (2.95) g m Let (see Appendix IV) c. -c. K = 1 / (Q where (Q matrix Q (2.96) 1 CC ).. IX ).. i s the i - t h element of the main diagonal of the covariance cc 11 cc and c. stands f o r the i - t h element of the vector c. 1 Also l e t (see Appendix IV) m = n(n + m)(N - 1) .2 £q = V'V (2.97) Then using (2.96), (2.97) into (2.95) y i e l d s c. - c. 1 "n(n + m)(N - 1) ^ ) > > ( v cc 11 (2.98) 1 , v / n ( n + m ) ( N _ 1 } 42 Using the t - d i s t r i b u t i o n i t i s possible to compute the length of the i n t e r v a l about the parameter c. which would include c. with a certain i probability. 1 For example, l e t the assigned probability be 0.9, then ^ n C n + DCN - 1)1 4 Y > = 0.9 where P{*} stands for the p r o b a b i l i t y the brackets. density of the expression between The number y can be found from tabulated standard tables using the number n(n + m)(N - 1) as an index. Ic. 1 If (2.99) - c.I = / ( Q l Y x 1 ).. , Hence , sY V M cc i i n(n + m)(N - 1 (2.100) TT 1) the augmented states and controls are not used the denominator under the square-root sign i n (2.100) w i l l be n(N - n - m). Thus using the augmented states and controls decreases the estimation error | c\ - c^| and r e s u l t s i n a better estimation accuracy. The range 2A of the confidence i n t e r v a l for a p r o b a b i l i t y 0.9 i s given by of 2A = [c. ± YAQ ).. , . VZ TV l ' c c i i n(n + m)(N - 1) (2.101) % 2.8 Discussion For the observable case equation the unknown parameters. Equation It computed without For the non-observable case the matrix (2.59) i s to be used instead of order to estimate the unknown parameters. 6^ i s required. i s used to estimate The matrix P^ i s recursively any need to invert a matrix. P^ does not e x i s t . (2.12) (2.12) The generalized matrix in inverse i s desirable to compute G^j i n a way s i m i l a r to that developed i n this Chapter for computing P^ in the observable case. This requires a computation algorithm to evaluate 0^ recursively without any need f o r a matrix inversion. In the next chapter an algorithm w i l l be developed i n order to compute 6* r e c u r s i v e l y and no matrix inversion w i l l be needed. This algorithm i s computationally simple and i s suitable f o r real-time computations . A s i m p l i f i e d version of the algorithm that w i l l be developed i n the next chapter can be used i n order to compute P^. with another method to compute P^. This presents us Due to the s i m p l i c i t y of this algorithm i t can be used f o r state estimation applications i n cases where the plant i s represented by a discrete model. Since the generalized inverse has numerous applications besides those which occur i n state and parameter estimation, the important topic of recursive computation of a generalized matrix inverse i s dealt with i n the next chapter i n a general way. 3. RECURSIVE COMPUTATION OF THE GENERALIZED INVERSE OF A MATRIX 3.1 Introduction Because of i t s value i n the solution of minimum norm problems i n f i n i t e dimensional spaces, the generalized inverse (or pseudo-inverse) finds important applications i n l i n e a r system theory. Applications of the generalized inverse to estimation and f i l t e r i n g problems are given i n [18], [6], [19]. [21]. Applications to control problems are given i n [20], A number of methods have been proposed for computing the generalized inverse A''" of a matrix A. The method proposed by Ben-Israel et. [22] a l . requires the formation of C = A*A, of A, and a series expansion of C. where A* i s the conjugate transpose The method proposed by Pyle [7] requires a two-stage Gram-Schmidt orthogonalization process and subsequent solution of sets of l i n e a r equations. posed by Rust et. a l . [8]. A simpler algorithm has been pro- The algorithm i s based on a two stage Gram- Schmidt orthogonalization process. Unlike the other two methods men- tioned above, i t does not require the generation of powers of matrices or the solution of sets of l i n e a r equations. However, by eliminating the Gram-Schmidt orthogonalization process, an even simpler recursive algorithm can be derived. 3.2 Properties of the Generalized Inverse Let A be an m x n complex matrix. The generalized inverse of A, as defined by Penrose [23], i s the unique n x m a t r i x A* which s a t i s f i e s the following r e l a t i o n s : AA A (3.1) I (3.2) = A (AA ) * = AA (3.3) (A A)* (3.4) 1 AA 1 X Following Rust e t . a l . , suppose that A can be partitioned into the form A = (R,S) (3,5) where R i s an m x n^ matrix with l i n e a r l y independent columns and S i s an m x n R. 2 matrix whose columns are l i n e a r l y dependent on the columns of The following relations h o l d : RR = I (3.6) I R 1 S = (R*R) R* (3.7) -1 (3.8) -RU U = RS (I + UU*) (3.9) 1 R I U*(I + u u * ) " ^ (In (3.6) and (3.10) (3.10) 1 I i s a unit n^ x n^ matrix). The generalized inverse i s given i n p a r t i t i o n e d form by (3.10) where ^ by (3.7) and ( 3 . 9 ) , respectively. and U are given Rust et. a l . show that the inverses i n (3.10) can be found by a two stage Gram-Schmidt orthogonalization process (the f i r s t stage generates R*)• The p o s s i b i l i t y of replacing both orthogonalization stages by a recursive algorithm i s considered i n the next two sections. 3.3 Recursive Evaluation of R The following relationships are required to develop the recursive algorithm. Taking the conjugate transpose of (3.1) and using (3.3) y i e l d s A*AA = A* (3.11) Let A = (V,W) (3.12) be a p a r t i t i o n i n g of A where V i s composed of l i n e a r l y independent columns. Equations (3.7) and (3.11) can be applied to V: v 1 (3.13) = (V*V) V* -1 v*vv = (3.14) V* From (3.1) and the p a r t i t i o n i n g (3.12) i t i s seen that AA V = V (3.15) Taking the conjugate transpose of (3.15) and using (3.3) (3.13) and (3.14) y i e l d s V ^ Let r_., (j = 1, columns of R. 1 = (V*V) V* = V -1 (3.16) 1 n^) denote the l i n e a r l y independent Consider the p a r t i t i o n i n g k+1 (3.17) k+1 where R, = r , , R 1 1 n i t follows that 1 = R, and where the b.'s are row vectors. J From (3.17) (3.18) " W+i k+ik+i +r Applying (3.6) to b yields (3.19) = I Applying (3.11) to 8^+1 the p a r t i t i o n e d form (3.17) y i e l d s 1 1 1 (3.20) r* k+1 k+1 P a r t i t i o n i n g of (3.20) and solving f o r b^-fi gives \ + i -1 <WW = Replacing A by r k i + ( I " Vk+i> (3.21) and V by R^ i n (3.16) gives ^W^+l (3.22) " Substituting (3.18) into (3.22) and using (3.19) y i e l d s w =V 1 WW' (3.23) - \v (3.24) Substituting (3.23) into (3.21) gives W = wW 1 where °k+l = ( r k l k+l r + " r k l\< k l r + ) _ 1 (3.25) + is a scalar. The recursive algorithm (3.17) (3.23) (3.24) and (3.25). ,1 •i the matrix product R^^k. n f o r the computation of R i s based on I t i s desirable, however, to compute (3.24) from a recursive formula. Substituting (3.23) and (3.24) into (3.18) y i e l d s = Vk" v i + ( i - \^ k it i ) r + + ( i - v£ (3 - 26) Since (see (3.3)) (3.27) i t i s seen that the second matrix on the right-hand side of (3.26) i s the outer product of two vectors Summary of the Algorithm for Finding R"*" Initialize. \+i" b B k+l k+1 ( r ][ r \ l k l = r + ( " ^ ( I " ) R* -1 - 1 \J " I + = r , R* = (R*R t+A^W k+i k+i ~ vAi < Set k = 1, R R ) V l W = \*i v i + (I -A^wW 1 k+1 +1 k+1 3.. k. = k+1 4. If k < n 1 go to 2, The above algorithm replaces the f i r s t Gram-Schmidt orthogonalization stage of Rust's algorithm, which, i n the interest of accuracy, reorthogonalizes each column a f t e r i t i s f i r s t orthogonalized. An addi- t i o n a l advantage of the above algorithm i s the manner i n which A i s p a r t i t i o n e d into the form (3.5). In the Gram-Schmidt process, the de- 49 t e c t i o n of a zero vector a f t e r orthogonalization (but before normalization) indicates that the vector being orthogonalized i s l i n e a r l y dependent on the previous vectors. This vector i s consequently moved i n t o S and the next column vector of A i s chosen. In the algorithm presented here, l i n e a r dependency i s detected when CL ^ (see (3.25)) i s K."t"-L zero. Since a ^ ^ i s required i n (3.24) no a d d i t i o n a l computations are + required to test f o r l i n e a r dependency. I f 0^+1 ^ s z e r o » ^ fc e associated vector i s moved into S and the next column vector of A i s chosen. 3.4 Recursive Evaluation of A Consider the p a r t i t i o n i n g k+1 \+l. " » k+l?' ( R S (3.28) \+l " k+1 where S. = s i s the f i r s t column vector of S and S = S. The generalized 1 1 n I inverse R can be obtained by the recursive algorithm given i n Section n ? z (3.3). Applying (3.10) to A^, i t i s seen that J (3.29) 1 C. = U*B. J 11 (3.30) U. = ( , u , U]L 2 , u.) = R S. .J J (3.31) where D. = (I + U.U*) 1 The matrix inversion lemma [9] 11 (3.32) - \ ' \li is - D k \ + i ( I "kVA'W^A + used to develop a recursive formula for B^. (3.29) and using (3.33) B 1 (3.33) Setting j = k, k + 1 i n gives (3.34) k i= " ViVWk+i^k ( I + where e k + 1 = .(i + u g B R u +1 k k + 1 (3.35) )" The recursive algorithm for the computation of A* given based on (3.28) (3.30) (3.31) is (3.34) and (3.35). • Summary of the Algorithm for finding A* 1. Initialize. Set k = 1, B . = D ' V = [I - ] (1 + u*u )~ u u*]R 1 1 1 I -1 V i • - k i\ \ it i k (I e R + 3. k = k + 1 4. If k ^ n^, go to 2, 5. C n = (R S I 2 n )*B 2 n + 2 A 3.5 ) B + 1 = A nv. 1 Examples The f i r s t two examples are taken from the l i t e r a t u r e . This allows a comparison to be made showing the r e l a t i v e s i m p l i c i t y of the proposed algorithm. Example 1. This example i s given i n [22]. I t i s required to compute the generalized inverse of A = 1 -1 0 0 0 1 -1 0 0 0 1 -1 -1 0 0 1 (3.36) It i s seen that (3.36) has the form (3.5) where 0 1 -1 0 1 -1 0 0 R = 0 0 1 -1 (3.37) s = The algorithm of Section 3.3 i s used to f i n d R Step 1. Initialization. R l 1 -1 0 0 " l r R* = | ( 1 , -1, 0, 0) Step 2, a 2 b 0 B, R 2 2 R = 3 = ±(1, 1, -2, 0) |(2, -1, -1, 0) 1 3 2 -1 -1 0 -1 2 -1 0 -1 -1 2 0 0 0 0 0 Step 2. (repeated) b„ = f (1, 1, 1, -3) 1 "4 f3 -1 -1 - l ) [2 2 - 2 2 - J = R 3 -1 -1 -1 2 2-2-2 1 1 1 - 3 = R = Having found R*, the algorithm of Section (3.4) i s used to evaluate A*. Since S consists of a single column vector s^, only the i n i t i a l i z a t i o n stage i s required. u i " R I S I = -1 -1 -1 Step 1. •i-i 3-3-1 1 1 3 - 3 - 1 - 1 1 3-3 Step 5 C =U A l 1 A Example 2. - | ( 3, 1, -1, -3) l Pl = 1 8 3 -3 -1 1 1 3 -3 -1 -1 1 3 -3 -3 -1 1 3 This example i s given i n [ 7 ] . (3.38) It i s required to compute the generalized inverse of 1 0 1 1 0 1 - 1 0 1 1 0 1 (3.39) 53 It i s seen that (3.39) has the form (3.5) where 1 0 1 R = 0 1 1 , l" 0 1 1 -1 0 s = (3.40) The algorithm of Section (3.3) i s used to compute R*. Step 1. Initialization R l " , RJ = l r \ (1, 0, 1) Step 2. 0U = 2 3 f (-1, 2, 1) B 2 R 1 = = f (2, -1, 1) 1 2 = 1 , 2 3 -l -1 1 2 V K Having found R^, the algorithm of Section (3.4) i s used to evaluate A*, Step 1. Initialization U = D" 1 = I - u (up l 1 Step 2. ( l+ u j u ^ u * 3 [0 1 IJ u ) 2 3 v l 2 J B 2 = - i -v ( ° 3 15 - l 5 3 )} h Step 5. C 2 = * 2 - T5 3 "o U B A Example 3. ( 3 -1 4 3 2 " 15 0 3 5 4 -5 -1 0 3 (3.41) This example i l l u s t r a t e s a case where A does not have the form (3.5) i n i t i a l l y . I t i s required to compute the generalized inverse of A = 0 1 1 -1 1 0 0 -1 -1 0 0 0 0 0 0 1 (3.42) Note that the t h i r d column of A i s the sum of the f i r s t two columns. Generally, the detection of l i n e a r dependency by simple observation i s not possible. The algorithm of Section (3.3) i s applied to (3.42). It i s seen from (3.36) that there i s no change i n the computation given i n Example 1 u n t i l k = 2. At this stage the algorithm gives hence the t h i r d column of A i s shifted i s chosen. This gives = 0, into S and the l a s t column of A = 1 and b = (0, 0, 0, 1) 3 1 2 3 1 ( 3 -1 -1 0 1 -2 cr 2 - 1 - 1 0 1 1-2 0 0 0 0 3 The interchange of the l a s t two columns of (3.42) i s accomplished by a permutation matrix P given by 55 1 0 0 0 1 0 0 0 0 0 0 1 0 0 1 0 (3.43) The algorithm of Section (3.4) i s now used to f i n d A*' where A = AP, A = AP (3.44) The generalized inverse A''" i s then given by (see [8]) A Step 1. 1 = PA (3.45) 1 Initialization u^ = R s^ = 1 0 0 -1 1 0 0 -1 0 0 0 3 Step 5. U A 1 1 B 1 3 1 = I » °» (1 -1 1 0 0 _ 1 0 -1 0 -1 ' 0 ) 0 0 3 0 (3.46) The generalized- inverse A-^ is- found- from (-3-. 45) and- is- A = -1 1 0 0 0 -1 1 0 0 0 0 3 (3.47) 4. TWO-STAGE ESTIMATION OF TIME-INVARIANT AND TIME-VARYING PARAMETERS IN LINEAR SYSTEMS 4.1 Introduction Because of i t s importance i n system engineering, considerable e f f o r t has been devoted to the problem of estimating the parameters of systems whose dynamic behaviour can be described by d i f f e r e n t i a l Two equations. d i s t i n c t methods f o r estimating time-invariant parameters have been brought to a high degree of refinement. i s based on a steepest-descent The parameter tracking method parameter adjustment of a quadratic error function or a Liapunov function. The error can be taken to be a response error or a generalized equation error. Detailed descriptions of para- meter tracking systems are given by Lion [2], Pazera and Pottinger [3], Rose and Lona [4]. The second method for parameter estimation i s based on optimum f i l t e r i n g . In i t s simplest form an optimum f i l t e r minimizes a mean-square error cost function subject to a prescribed dynamical constraint. Detailed descriptions of optimum f i l t e r s and an extensive bibliography are given by Sage and Melsa [9]. The estimation of time-varying more d i f f i c u l t problem. often give acceptable parameters i s generally a much Parameter tracking systems and optimum f i l t e r s estimates provided that the parameters vary c i e n t l y slowly with respect to time. suffi- Lessing and Crane [10] have made use of i n t e g r a l transforms to extend the c a p a b i l i t y of parameter tracking systems to track time-varying the i n t e g r a l transform parameters. P r a c t i c a l implementation of approach requires the assumption that the parameters can be represented by polynomials of known order. The estimation of time-varying parameters by optimum f i l t e r s requires e s s e n t i a l l y the same assumption. differential meters. I t i s assumed that the parameters are solutions of equations of known form with unknown time-invariant para- The problem i s then reduced to the estimation of time-invariant parameters of a higher order augmented system. I t i s the purpose of this Chapter to develop a computationally simple and accurate method for estimating time-varying parameters which i s not r e s t r i c t e d by a p r i o r i assumptions about the dynamical behaviour of the parameters. 4.2 Problem Formulation Consider the continuous-time x p = Ax + Bu, p l i n e a r system x (t ) = x p o po (4.1) where x i s an n x 1 dimensional state v a r i a b l e vector and u i s an m x 1 P dimensional control vector which i s a known function of time. The system matrix A i s considered to be decomposable into a time-invariant component A^ and a time-varying component A^: A =• A f + A (-4.2) v The most general estimation problem associated with (4.1) i s to estimate x , A and B given an observed P continuous-time over an observation i n t e r v a l t measurement vector defined 4 t <_ t ^ , and given complete s t a t i s t i c a l information about measurement noise and system noise. However, i f emphasis i s placed on computationally simple methods then a less ambitious problem must be considered. Suppose that B i s known and that the state x^ i s 58 measured during the observation i n t e r v a l . The problem to be considered i s the estimation of A given x^ and given no information about noise statistics. The estimation i s performed sequentially i n two stages. During the f i r s t stage an estimate of the time-invariant A^ i s made given an estimate of A^. During the second stage a new estimate of A^ i s made using the estimate of A^ found during the f i r s t stage. This section deals with the second stage of i d e n t i f i c a t i o n . Let z. = x - x, 1 p 1 x, = A.x. + Bu, 1 f 1 It follows from equations z., = 1 (4.3) x.(t ) = x 1 o po (4.4) (4.1), (4.3) and (4.4) that A.Z.. f 1 + A x vp = A.z, .+ 6(x )a, t l p — z,(t ) = 0 1 o (4.5) In equation (4.5), a i s a vector whose components are the time-varying elements of A . v The equation A x = 0(x )a v p p - (4.6) defines a matrix 0 whose elements are l i n e a r functions of x . Let P g = x p - x x h A g(t) - g ( t - T ) T ( 4 > y ) In (4.7), x^ i s defined by (4.4), x^ i s assumed known (through measurements) and h. i s an approximation to g, with the accuracy of the approximation depending on the fixed-time increment T. Consider the following quadratic cost index J where z^ and = jfl f [||g-^ l| Q + 2 o 1 1 llh- (4.8) ^||" Q ]dT 2 2 are defined by (4.5) where the unknown vector a. be r e - placed by an estimate a_. The weighting matrices be p o s i t i v e d e f i n i t e . and are taken to I t follows from (4.3)-(4.8) that l i m J = 0 proT->0 vided that the time-varying parameter vector c i i s known exactly. Con- sequently, (4.5) and (4.8) can be used to formulate a problem i n optimum estimation. In the problem formulation the unknown time-varying para- meter vector i s considered to be a control vector f o r the system described by (4.5). The optimum control minimizes J where g and h are known functions of time. 4.3 Estimation of A^. Algorithm I. The minimization of J (see (4.8)) subject to the dynamical cons t r a i n t given by (4.5) i s a standard problem i n optimal control. Hamiltonian The f o r t h i s problem i s (prime denotes transposition). H = p ' C A ^ + 6(x ) a ) - - | ( g - z ) ' Q ( g - z ) - j ( h - A z -Q (x ) a) ' Q (h-A z.j-6 (x ) a) 1 1 1 f ± 2 f (4.9) The costate equation i s p = -H z = -A^p + Q ( g - z ) - A^Q (h-A 1 and the gradient condition i s 1 2 fZ;L - 6(x )a), p ( t ) = 0 p f (4.10) 60 H = 9' (x )[p+ Q,(h - A,z a p L r 1 G(x )a)] p — = 0 (4.11) Solving (4.11) for the optimum control a y i e l d s a = (9»(x )Q9(x )) V(x )[p + Q(h 2 p 2 - A z^ ] (4.12) f Existence of an optimal control for the formulated problem can be viewed as a c o n t r o l l a b i l i t y condition i n a function space S of elements (g, h ) , where h i s either equal to g or an approximation to g. Taking (4.8) to define the norm, then the system described by (4.5) defined to be controllable i n S i f element (z^ is a control ji exists such that the - g, z^ - h) has minimum norm. S a t i s f a c t i o n of this con- t r o l l a b i l i t y condition requires the existence of the matrix inverse in (4.12). Equations ( 4 . 5 ) , boundary-value problem. (4.10) and (4.12) constitutes a two-point Due to the l i n e a r i t y of the equations the solution i s e a s i l y found using standard methods. Let L(x ) ^ (6'(x )Q6(x )) V(x ) p 2 p p Substituting (4.12) into (4.5) I C. (4.13) p and (4.10) i t i s seen that 9(x )L + r (4.14 -r' 4 where C ± = (I C 2 = - Q 6(x ) l LQ ) 2 A f + A'Q^ 9(x )LQ h p Q 1 § 2 - A^Q (I 2 0(x ) LQJh (4.15) 61 and where I i s the unit matrix. The general solution of (4.14) has the form - *(t , t) f Q p(t ) P(t ) c f + (4.16) y 2 where A y = A 1 iKt, (4.17) T) r(x)dx o L2 J fc y and where i|;(t, x) f t = / i s the state t r a n s i t i o n matrix for (4.14). I t i s seen from (4.16) that y can be found by solving (4.14) taking z ^ ( t ) = 0 = p ( t ) Q for the i n i t i a l values. (4.14) taking -Kt Q The state t r a n s i t i o n matrix i s found by solving , t ) = I, r = 0. Substituting the boundary conditions into (4.16) and parti- tioning ^(t^,, t ) i n the form *(t , f t ) = "*11 |*12 *21 j*22 Q (4.18) i t i s seen that p(t ) Q .-1 -*22 y 2 (4.19) The i n i t i a l costate given-by (4.19) i n conjunction with (4.5), (4.10)' and t Q (4.12) defines the optimal estimate of a. over the observation i n t e r v a l <^ t ^ t^. The computational procedure w i l l be referred to as algorithm I. Minimization of a cost index of the form given by (4.8) subject to the dynamical constraint of (4.4) has found extensive applications i n state estimation (Sage [11]). in the problem formulations. There i s , however, a fundamental difference In the state estimation formulation, (4.5) and (4.8) are replaced by x = Ax + Bu, x ( t ) = x P P o po P ~ = J 1 t / t f [ ll*p - * M 1 2 q P + ll*p - % \ \ (4.20) 2 ^ <- > 4 21 The minimization of J subject to the dynamical constraint of (4.20) i s a problem i n optimum f i l t e r i n g . finding an approximation x P The problem can be considered as one of to z so that the error (z - x ) has minimum P P P norm subject to a cost penalty associated with the norm of the control effort. The s o l u t i o n of the two-point boundary-value problem defined by (4.20) and (4.21) can be obtained by invariant imbedding or by use of a R i c c a t i transformation. E i t h e r approach results i n the Kalman f i l t e r algorithm, f o r generating, the. optimal, estimate, x^ of the state. x the measurements z over an observation i n t e r v a l t < T < t. p o = = p given. Comparing the second terms i n the integrands of (4.8) and (4.21) i t i s seen that there i s a fundamental difference between the state and parameter e s t i mation problems. This difference manifests i t s e l f a n a l y t i c a l l y i n the c o n t r o l l a b i l i t y condition (the existence of L ( x ) , see (4.13)) which p does not a r i s e i n the s o l u t i o n of the state estimation problem. 4.4 Estimation of a P e r i o d i c A . , v Algorithm I I . S Several questions a r i s e concerning algorithm I. the estimation accuracy of The f i r s t question concerns the e f f e c t of measurement noise on the optimal estimate. I t i s seen from (4.7) that h, since i t approxi- mates g, may lose i t s s i g n i f i c a n c e i f the data representing g i s noisy. To some extent this loss of s i g n i f i c a n c e can be taken into account by decreasing the absolute value of the elements of the weighting matrix Q (see (4.8)). 2 A second question arises concerning the s a t i s f a c t i o n of the c o n t r o l l a b i l i t y condition. Algorithm I must be modified i f p r a c t i c a l answers to the above questions are to be found. When no information i s available about noise s t a t i s t i c s there i s no unique best method f o r f o r mulating an estimation problem, since the choice of a cost index, such as the quadratic cost index given by (4.8), i s an a r b i t r a r y one. Con- sequently, estimation algorithms can only be judged on the basis of trade-off between computational complexity and estimation accuracy. the above discussion i n mind, an algorithm i s developed With for estimating time-varying p e r i o d i c parameters. Let x 3 = z - x x 2 =-A x- ••+ e(-x^);a., f (4.22) 2 x- (t ) = 0* 2 2 Q * (4.2-3) It i s seen from (4.5) and (4.20) that X Equation 3 = A f 3 X + 9 ( x 2 + X 3 -' ) X 3 ( t o ) = (4.24) 0 (4.23) i s taken to represent the dynamical constraint for the new problem formulation. the dynamics are noise-free. I t i s seen from (4.1) and (4.23) that A cost index of the form given by (4.8) must be s p e c i f i e d to complete the problem formulation which i s based on the following assumptions: (2) (1). L(x^) i s well-conditioned (see (4.13)). an estimate I of a i s a v a i l a b l e and an estimate x^ of x^ can be found by solving (4.24) using successive approximations. x • x 3 0 Z = A x f 3 + Q(x )a, 2 = A j L + 0(x.)a, J. I. 1 — x 3(t ) = That i s , x^ i s defined 0 0 x_(t ) = 0 I o (4.25) (3) The cost index i s (see (4.8) and (4.22)) - H X 3 o 2 q 2 1 + M H " X 3 " X 2 H\] D (4.26) T Successive approximations converge i f the norm of a i s s u f f i c i e n t l y small. This r e s t r i c t i o n i s imposed to achieve convergence i n the two- stage estimation of A., and A . f v (If A., i s known, then the f i r s t estimation f stage i s unnecessary; (2.4) could then be used to define x^ and no norm constraints need be imposed on a ) . I t i s always possible to avoid ill- conditioning of L ( p ) by modifying the values of i t s time-varying elements. X The modifications, however, cannot be a r b i t r a r y . They must enter into the problem formulations i n a mathematically consistent manner. In the present problem formulation i t i s seen from (4.22) and (4.26) that L(x^) i s replaced by- L ( X j - ) . . Equations- (4. 23) , (4..25). and. (4..26). define-a--two- point boundary-value problem f o r the optimum estimate. Let T be the length of a sampling i n t e r v a l and l e t t,. - t be a fixed estimation i n t e r v a l . = NT At the end t ^ of an estimation i n t e r v a l a new measurement sequence i s taken over an observation i n t e r v a l t , < t < t + HT. f That i s , the next estimation i n t e r v a l i s t + £T < t < o (N + £)T (see F i g . (4.1)) . •OBSERVATION INTERVAL FOR NEW MEASUREMENTS 0I T Ha NT(UN)T h =H ESTIMATION INTERVAL OPERATION PERIOD flG.(4.l) THE RECURSIVE INTERVALS. OBSERVATION - ESTIMATION In order to solve Eqn. to extrapolate a^(t) . i s to f i t a polynomial (4.25) over this i n t e r v a l i t i s necessary This can be done i n a v a r i e t y of ways. to sampled values j i ( t o One method + kT) , k = 0, N, and to use the polynomial approximation for extrapolation. ..Alternatively, , i f i t i s known that the length of the estimation i n t e r v a l NT i s greater than the smallest period of the p e r i o d i c components of a_, then the known p e r i o d i c i t y of a_ can be used f o r extrapolation. The following i s a summary of Algorithm I I : 1. A measurement sequence g ( t Q + kT), h ( t Q a_(t + kT) over an estimation i n t e r v a l k = 0, + kT) and an estimate N are assumed known. Q Take a new measurement sequence over an observation i n t e r v a l k = N, N + I and define a_ over this i n t e r v a l by extrapolation. 2. new estimation i n t e r v a l for a new (4.26)). 4.5 Solve the two-point boundary-value problem associated with the estimate (Equations (4.23)*, (4.25) and The two steps of the algorithm are repeated r e c u r s i v e l y . Estimation of a Periodic A^. Algorithm I I I . Algorithm II makes recursive use of measurement data and vaious estimates. pre- I t i s based on the solution of an optimization problem having.-the noise free, dynamics, of. (4.. 23). and a cost index given by (4.26) where x^ i s generated i n a smoothed form by integrating (4.25). Conse- quently, i t can be anticipated that algorithm II w i l l be superior to algorithm I i f measurement noise i s present. There i s no guarantee, however, that the c o n t r o l l a b i l i t y condition i s s a t i s f i e d . L(x^) may be i l l - c o n d i t i o n e d . As mentioned i n the previous section, one method for overcoming i l l - c o n d i t i o n i n g i s to modify the elements of L by modifying That i s , the problem formulation. time-varying Let z^ be an estimate of z^ and l e t X 2 = A f 2 X + 9 ( x l + V - ' 2 x ( t o = J ° (4.27) represent the dynamics associated with a t h i r d problem formulation. The cost index i s taken to have the quadratic form (see (4.8)). J = \ / o t f t [||g - * \ \ \ + ||h - x | | Q ] d T 2 2 2 2 (4.28) Comparing (4.27), (4.28) with (4.5), (4.8) i t i s seen that z^ should be chosen so that x approximates z^. 2 (Ideally, i n the noise-free case and when L(x ) i s well-conditioned, then the choice z, = z, can be made. P 1 This r e s u l t s i n x 2 1 = z^ and Algorithm III degenerates into Algorithm I ) . I f L(x^) i s i l l - c o n d i t i o n e d Algorithm II w i l l f a i l . This i l l - c o n d i t i o n i n g can be overcome i n Algorithm I I I provided z^ i s chosen so that L(x^ + z^) i s well-conditioned. A further constraint i s imposed on z^ by the require- ment that the two-stage estimation procedure converges. This requirement, which was discussed i n the previous section, l e d to the use of successive approximations and the assumption that the norm of a_ must be s u f f i c i e n t l y small. I t i s evident that z^ i s constrained by several c o n f l i c t i n g r e - quirements. Some sort. of. mathematically tractable, approach must be taken, which accounts f o r the various trade-offs associated with i l l - c o n d i t i o n i n g , dynamical approximations and two-stage convergence. Solving (4.5) by successive approximations using z^ = 0 as an i n i t i a l i z i n g function and an estimate a. gives x 2 - A x f 2 + 6( )I, x (t ) = 0 XL 2 x„ = A.x„ + 9(x. + x ) a , 0 J 1 J 1 Z — Q x (t ) = 0 0 J O (4.29) 6; for the next two approximating functions. z 1 = ax 2 Let z^ be defined by + 3x (4.30) 3 The general structure given by (4.30) allows the trade-offs mentioned above to be investigated as a function of two scalar parameters a and $. The choice a = 0, B = 1 r e s u l t s i n z^ being defined by successive approximations given by (4.29). Numerous algorithms f o r the solution of l i n e a r algebraic equations are given i n [24]. convergence r e s u l t s when a = 2, 3 = -1. Abramov's method for accelerated (The j u s t i f i c a t i o n f o r using these algorithms i s e a s i l y seen by the d i s c r e t i z a t i o n and subsequent rearrangement of (4.5) into a l i n e a r algebraic equation). Algorithm I I I i s based on the solution of the two-point boundaryvalue problem defined by (4.27), (4.28), (4.29) and (4.30). parameters a. and p can. be. determined by numerical. experimentation., or, they can be chosen on the basis of successive approximations Abramov's method (a 4.6 The scalar (a = 0, 3 ='. 1) or = 2, 3 = -1). TwO-Stage Estimation of A^ and A^. During the f i r s t stage an estimate of A^ i s found given an estimate of A^. Three d i f f e r e n t algorithms have been developed i n the previous sections f o r the second stage during which an estimate of A^ i s found given an estimate of A^. Together The two stages are performed sequentially. they comprise a composite estimator f o r A. For reasons of s i m p l i c i t y , the discussion i s l i m i t e d to the development of a f i r s t algorithm which i s associated with second stage Algorithm I. stage The modifi- cations required i f second stage Algorithm II or I I I i s used are s t r a i g h t forward. Equation (4.4) defines the dynamics associated with the f i r s t stage of estimation. I f x^ were available through measurements, (4.4) could be d i s c r e t i z e d and estimated by use of standard minimum mean- square error estimation procedures. However, x^ i s not d i r e c t l y available. Furthermore, a two-stage estimation i s proposed and the i n t e r a c t i o n between stages must be accounted for i f a convergent sequence of e s t i mates i s to be obtained. Chapter 2 has suggested augmenting (4.4) by additional states and controls to reduce the e f f e c t of quantization noise which i s i n t r o duced by d i s c r e t i z a t i o n . It i s shown that the use of augmented states and controls results i n rapid convergence of the estimates. The augmented states and controls are defined by successive integrations of (4.4): y. = A y . x y i ( t o ) = 1 A y i + B u. i fx J = i I 0 i-l d T ft o , i > 1 (4.31) A y > u. = •/" u. ^dx, l o i - l t y l = X l u, = u 1 (4.32) For' composite estimation" of" A, successive integrations has' the b e n e f i c i a l - e f f e c t of reducing the influence of the time-varying components of A^ on the augmented states. The augmented system of equations (4.31) are d i s - cretized and rearranged into a set of l i n e a r algebraic equations. The optimum estimate f o r A^ i f taken to be the solution with minimum meansquare error (For d e t a i l s see Chapter 2). In the case of composite estimation the above procedure must be s l i g h t l y modified since x^ i s not d i r e c t l y available. Consequently, x^ must be estimated from e x i s t i n g data. From (4.3) and (4.5) i t is seen that XJL = x - (4.33) p where z-1= A fz, l+ 9(x p)a,—' z.(t 1 o) = 0 (4.34) c In (4.34), and a. are estimates of A^ and a_ respectively. The f i r s t stage of estimation consists of a minimum mean-square error solution of the d i s c r e t i z e d equations (4.31) using (4.33) and (4.34) to define x^. The new estimate of A^ i s used i n place of the old estimate provided that i t r e s u l t s i n improved prediction of the state of the system. The following cost index can be used as a basis f o r making this decision: J In (4.35), change of l = I t/ 'l*l' p f x " l X " S iM Q d" 2 T ( 4 3 i s a p o s i t i v e d e f i n i t e weighting matrix. ' 3 5 ) The rate of- between estimation i n t e r v a l s i s useful i n making decisions concerning the displacement length (£T) of the observation i n t e r v a l . displacement can be increased when the rate of decrease of The goes below a preset threshold. 4.7 Examples In a l l examples the following computational procedures were adopted. The two-point boundary-value problem was solved using a fourth- order Runge-Kutta method with an integration step-size of 0.025. The time-varying estimates a.(t) were generated i n the computer by storing the values at grid points and taking a piecewise constant approximation 70 a(t)=a(t +kT), kT < t<(k + 1)T). Q In the case of algorithms II and I I I , extrapolation of _a was performed using a very simple constant approximations, _a(t^ + kT) = a_(t^) , k = 0, , and the i n i t i a l i z i n g functions f o r successive approximations were taken to be i d e n t i c a l l y zero. 4.7.1 Example 1. Consider a f i r s t order system x = -(1 + b sino)t)x + u (4.36) where the control u i s taken to be a step function of magnitude equal to 3 and x(0) = 0 with T = 0.025, I = 4, N = 80 (4.37) The length of the whole operation period for a l l the examples i s taken to be 10 seconds. The length of the estimation i n t e r v a l , as seen from (4.37) i s NT = 2 seconds and the same has been considered for a l l the subsequent examples. of a(t) = b (4.1). The maximum percentage error of the estimate s i n wt for d i f f e r e n t b and w values are given i n Table Table (4.2) shows the r e s u l t s when a(t) i s a square wave: a(t) f 0.5,. =\ L-0.5, 0 < L< 0..5. 0.5 4 t < 1.0 (4.38) In the case of algorithms II and I I I , two estimation i n t e r v a l s were required to obtain the results given i n the tables. with successive approximations (a = 2, g = -1). (a = 0, 6 = 1 ) Algorithm I I I was and Abramov's method tried 4.7.2 Example 2. Consider the second order system x^ = - ( 1 + b^sin io^t)x + u^ 2 x where 2 = -(2 + t> s i n o) t)x 2 2 1 - x 2 + u (4.39) 2 i s taken to be a step function of magnitude equal to 2 and { t , 0 < t <_ 3 1 , t > 3 (4.40) In addition to (4.37) the following values are used: b^ = 0.4, b 2 = 0.6, co^ = 3, = 2. x^(0) = 0 = x ( 0 ) , 2 Table (4.3) gives the maximum per- centage estimation errors for the parameters a^(t) = b^ s i n co^t and a (t) = b 2 2 s i n o) t, respectively. In the case of algorithms II and I I I , 2 two estimation i n t e r v a l s were used. 4.7.3 Example 3. Consider the t h i r d order system x^ = (1 + b^ s i n to^t) x x 2 = -(2 + b x 3 = -(3 + a ( t ) ) x where the controls u^ and u 2 3 3 s i n ui^t) x 1 + u^ 2 - 2x 3 - x 2 3 + u 3 4- u 2 (4.41) are taken to be step functions of magnitudes 0.1 and 0.5, respectively and where ft , t< 3 " 2 3 l l . « > = (4.42) The parameter ^ ( t ) i s periodic and given by 1 , 0 < t < 0.25 a (t) U l , 0.25 < t <_0.5 3 (4.43) In addition to (4.37) the following values are used: x ( 0 ) = 0, b 3 ± = 0.3, b = 0.5, o> = 3., OJ = 2. 2 1 and 4.7.4 Table (4.4) gives the 2 maximum percentage estimation e r r o r s . x^(0) = ^ ( 0 ) = In the case of algorithms I I I I I , four estimation i n t e r v a l s were used. Example 4. The next two- examples i l l u s t r a t e the composite two-stage estimator. In both examples the f i r s t - s t a g e of estimation was i n i t i a l i z e d by taking a. = 0 as the i n i t i a l estimate f o r c i and only one set of augmented states and controls was used ( i = 1, 2, see (4.31) and (4.32)). Consider the second order system x^ = (1 + b^ s i n to^t)x + u^ 2 x = (b-- s i n o) t.)x - + 2 2 a 2 2 where the controls u^ and u 2 2 the following values are used: 2 3 = 1. + U 2 (4.-44)- are taken to be step functions of magnitude equal to 2. and 3., respectively. co^ = 1, c o = 2, a x 1 In (4.37), £ i s taken equal to two and x^(0) = = ^» ^1 = 0.3, t> = -0.5, 2 After two estimation i n t e r v a l s the maximum per- centage estimation error f o r a^ was 1.2% and a f t e r ten estimation i n t e r v a l s , the maximum percentage estimation errors f o r a 2 and a 3 were 1.8% 72 and 1.0%, r e s p e c t i v e l y . 4.7.5 Example 5. Consider the second order system s i n ui^t)x^ + x^ = (1 + X 2 = a 2 l X + a 3 2 X + U 2 (4-45) where b^ = 0.4, a^ = -2, a^ = - 1 , co^ = 2 and where a l l remaining data i s the same as i n Example 4. A f t e r the second estimation i n t e r v a l , the maximum percentage estimation errors f o r a^, a^, and a^ were 0.9%, 0.4% and 0.3%, r e s p e c t i v e l y . 4.7.6 Example 6. In the previous examples the measurements were noise-free. Example 1 has been t r i e d with zero mean gaussian measurement noise with unity standard deviation and the random noise generator has been randomly initialized. A peak r.m.s. s i g n a l to noise r a t i o of 10 was taken. following values were chosen: The b = 0.4, co = 1, a = 2, (3 = -1. The maximum percentage estimation errors f o r algorithms. I.,.. I I and. III. were 6%., 3.2% and 2%, respectively. In the case of algorithms II and I I I , seven and four estimation i n t e r v a l s were used, respectively. The examples show that the proposed algorithms and the two- stage estimation method can give accurate estimates of l i n e a r system parameters. I f the measurements are noise free, Algorithm I gives the most accurate estimates. In the case of measurement noise, Algorithms II and I I I give better estimates. In the above examples, no problem 11 was encountered with i l l - c o n d i t i o n i n g i n any of the above cases. 75 b 1 W / 0.5 1 1.0 0.6 % 3-0 10 % 2.0 0.7% 3.0 0.9 3.0 0.93% 0.75% 15 % JT 777 °(=2, 3=-1 777 TABLE TABLE (4.3) 0.85% % TABLE 1.0% 0.95% (4.1) TABLE (4.4) (4.2) 5. OPTIMUM CONSTANT FEEDBACK GAINS FOR LINEAR SYSTEMS AND THE PROBLEM OF TRAJECTORY SENSITIVITY REDUCTION 5.1 Introduction The optimum s o l u t i o n of the l i n e a r state regulator problem with quadratic cost has assumed a central p o s i t i o n i n the development of modern control theory. The optimum gain, however, i s i n general time-varying and this has severely r e s t r i c t e d p r a c t i c a l applications. P r a c t i c a l considerations generally specify constant or piece-wise stant feedback gains. con- Considerable e f f o r t has gone into the development of design procedures f o r constant feedback gains and i n the design of s p e c i f i c optimal c o n t r o l l e r s [25-29]. There are, however, computational as w e l l as p r a c t i c a l d i f f i - c u l t i e s associated with optimum constant gains. Specifying a constant gain r e s u l t s i n a constrained optimization problem and the s o l u t i o n of the associated two point boundary value problem i s generally more d i f f i cult. Furthermore, the optimum constant gain depends on the i n i t i a l conditions. One method that has been proposed to overcome these d i f f i - c u l t i e s i s to introduce the trace of a cost matrix which then allows ini- t i a l conditions to be eliminated i n the problem formulation by a r b i t r a r i l y assuming a uniform d i s t r i b u t i o n of i n i t i a l states over a sphere [25-27]. From a design standpoint, however, i t i s important to i n v e s t i - gate the e f f e c t of i n i t i a l conditions on optimum feedback gains before computing average suboptimum gains. Furthermore, i n many p r a c t i c a l s i t u - ations, state t r a j e c t o r i e s are controlled, and the effect of i n i t i a l conitions on optimum gains can become important. This occurs, f o r example, i n c e r t a i n schemes f o r combined estimation and c o n t r o l . During a f i n i t e observation i n t e r v a l , measurements of the state are taken and used to determine parameter estimates of a l i n e a r model. The l i n e a r model i s used to compute a control strategy which i s applied i n a subsequent control i n t e r v a l [11], [ 5 ] . systems. This scheme seems p a r t i c u l a r l y a t t r a c t i v e f o r adaptive I t s p r a c t i c a l implementation, however, requires an on-line computational algorithm f o r determining optimal or suboptimal Trajectory s e n s i t i v i t y minimization i s an important problem that has received considerable attention [30-37], gains. practical A widely adopted method i s based on augmenting the state vector with the s e n s i t i v i t y vector and formulating a l i n e a r state regulator problem f o r the augmented state. However, i f time-varying gains are permitted, the problem, i s i l l - p o s e d [35, 38]. A well-posed, traj.ectory. s e n s i t i v i t y minimization problem requires that the gain be constrained to be time-invariant. In general, the s e n s i t i v i t y vector i s affected by i n i t i a l conditions. Consequently, i f feedback i s to minimize trajectory s e n s i - t i v i t y , the e f f e c t of i n i t i a l conditions on the optimum constant must be evaluated. I t may be possible to choose a suboptimal gain average gain. A u n i f i e d treatment of several problem formulations dealing with time-varying gains and constant gains i s presented. Extensive" use i s made of the known r e s u l t s f o r the l i n e a r state regulator problem to develop a simple algorithm which results i n rapid convergence. An approximate method i s developed to replace the nonlinear two point boundary value problem by an algebraic optimization problem. This method appears to be suitable f o r on-line computation of suboptimum constant feedback gains. 5.2 Problem Formulation Consider a l i n e a r time-varying x = Ax + Bu, system x(0) = X (5.1) q where the state x and the control u are n x 1 and m x 1 vectors, respectively. The optimal control which minimizes the quadratic cost index (prime denotes transposition) 1 f = ± / (x'Q-x + u'Ru)dt, o 2 o 1 t J where (5.2) i s p o s i t i v e semidefinite and R p o s i t i v e d e f i n i t e , i s given i n feedback form by u = -R The time-varying gain matrix •''B'K.jX. (5.3) i s found by solving a matrix R i c c a t i d i f f e r e n t i a l equation and has the property that = K|. To ease the notational burden, argument values of functions are omitted when they appear unessential and only one parameter q i n the matrix A i s considered to have a s i g n i f i c a n t e f f e c t on s e n s i t i v i t y . Let z = represent the dq closed-loop t r a j e c t o r y sensitivity vector. A feedback structure of the form u = -R B'(K x + K z) (5.4) -1 2 i s postulated, where = K| and = K^. This s p e c i f i c structure i s postulated because i t resembles the feedback structure of the optimal control r e s u l t i n g from solving the state l i n e a r regulator problem (see (5.3)). Substituting (5.4) into (5.1) and neglecting the seconds ' order s e n s i t i v i t y function —j 9q yields z = 1^ x + (A - F K j z , 3q J. z(0) = 0 (5.5) where F = BR "*"B'. By defining an augmented vector y' = (x', z'), (5.1) and (5.5) can be combined into a composite vector d i f f e r e n t i a l y = (A - FK)y = Cy, y(0) = y equation (5.6) r where F 0 A = - A K , K = 3q Equation 1 0 J (5.6) represents 0 l K x 2 ' o y L 2 K K (5.7) = 3 (5.1) and (5.5) when = K^. This constraint i s relaxed i n (5.6) i n order to obtain a u n i f i e d treatment of the l i n e a r and nonlinear formulations. The design objective i s to choose gain matrices and so that (5.2) i s minimized subject to the a d d i t i o n a l constraint that the feedback (5.4) reduces, s e n s i t i v i t y vector z. i n some sense, the closed-loop t r a j e c t o r y A mathematically tractable approach that has been widely adopted i n s e n s i t i v i t y studies i s to neglect second-order s e n s i t i v i t y terms so that (5.6) i s v a l i d and to introduce a cost index J(K r K, t ) = j / f f (5.8) y'(Q + KFK)y dt where Q i s an augmented p o s i t i v e semidefinite weighting matrix [32-38], The optimum gain matrix minimizes (5.8) subject to the d i f f e r e n t i a l constraint (5.6). I t w i l l be seen that the s e n s i t i v i t y problem formu- lated by (5.6) and (5.8) i s well-posed only i f the gain matrices are constant. The argument l i s t i n the cost index (5.8) i s used to d i s t i n g u i s h between d i f f e r e n t problem formulations. where Let I(K^, K, t^) and I(K^, K, t ^ ) , i s a r b i t r a r y , represent the l i n e a r and nonlinear formulations, respectively ( I f (5.6) i s rewritten i n the form y = Ay + Fu, then for the l i n e a r formulation the matrices A and F are completely known and independent of K^. This r e s u l t s i n l i n e a r dynamic constraints. For the nonlinear formulation the parametric matrix A depends on the unknown gain matrix K^. This r e s u l t s i n nonlinear dynamic c o n s t r a i n t s ) . The Hamiltonian for e i t h e r formulation i s H = p'Cy - |- y'(Q + KFK)y (5.9) where the costate vector p i s defined by p = -C»p + (Q + KFK)y, Gradient matrices p(t ) = p f f = 0 (5.10) (see Appendix III) are used to derive the gradient condition for a minimum. In order to use gradient matrices a l l matrix elements must be independent and cannot be constrained by the symmetry condition K = K'. This i s e a s i l y accomplished by taking K = K + K', where K i s a matrix with a r b i t r a r y elements. In order to introduce a compact notation which proves useful i n subsequent discussions, i t i s assumed that F i s time-invariant. The gradient condition f o r the l i n e a r formulation I(K^, K, t^) can then be expressed 8J(K K,t ) i — — = 3K i n the form (see Appendix V) t / -_-5-dt=$ KF + FK$ J o 3K yy +F$ yy F = 0 +<S> py yp (5.11) where $ yv k o y Jv » dt (5.12) 81 (Equation (5.11) remains v a l i d i n the general case of time-varying F provided that KF and F are moved under the respective i n t e g r a l signs. Keeping this i n mind, i t i s seen that the form (5.11) i s a c t u a l l y not a r e s t r i c t i o n but rather a notational convenience). The gradient condition for the nonlinear formulation I(K^, K, t^) d i f f e r s from (5.11) i n an a d d i t i o n a l term which arises from the dependency of A on K^. Let p' = (p|, p p be a p a r t i t i o n i n g of the costate vector which corresponds to the p a r t i t i o n i n g of y and l e t K^ = L + L', where L i s a r b i t r a r y . The a d d i t i o n a l term which i s associated with K^ i s a t ~ f 3L o f Tr(-FK p!)dt 1 2 lZ r (5.13) = F$ + $ F p z zp„ Introducing the augmented matrix F$ E = + $ P z Z p 2 2 F ' 0 ! 0 1 (5.14) 0 allows the gradient condition for the nonlinear formulation to be expressed i n the form 8K " = E + $ KF + FKo + F'f + $ F = 0 yy yy py yp Minimizing J(K^, K, t^) i s a constrained optimization problem. (5.15) Conse- quently J(K*, K*, t ) = J ( K f 1 } K, t ) f where K* i s the optimum gain matrix f o r the nonlinear formulation. (5.16) If K^ = K*, i t i s seen from (5.16) that the minimum for the l i n e a r formulation 82 i s achieved when K = K*. The condition E* = 0, which i s s a t i s f i e d when F$ p z + $ 2 zp F = 0 (5.17) 2 i s an optimality condition. I t i s i n t e r e s t i n g to note that i f time- varying gains are permitted, the optimality condition (5.17) becomes a s i n g u l a r i t y condition which i s impossible to s a t i s f y [32, 37]. indicates that the s e n s i t i v i t y problem formulated with This time-varying gains i s i l l - p o s e d and cannot, i n general, meet the design objectives. Equations (5.11) and introducing a matrix S defined FS$ + $ yy (5.15) can be s i m p l i f i e d and u n i f i e d by by SF = -F$ yy - $ py F - E (5.18) yp (Several methods are available for solving an equation of the form (5.18) for S. See Appendix I I ) . Substituting (5.18) into (5.15) y i e l d s 9J(K K,t ) \ = $ (K - S)F + F(K - S)$ =0 9K yy yy (5.19) v Equations (5.6), (5.10) and (5.19) define a nonlinear two- point boundary value problem which must be solved i n order to determine the optimum gain matrix. 5.3 U n i f i e d Treatment of D i f f e r e n t Formulations Because i t i s a constrained optimization problem, the computation of optimum constant gain matrices computation of unconstrained i s generally more d i f f i c u l t than the time-varying gains. It i s only i n the s p e c i a l case of time-invariant systems and i n f i n i t e t f that the problems are computationally equivalent. From an engineering point of view, however, constant gain matrices are easier to implement and therefore preferable, despite the possible computational d i f f i c u l t i e s . evidently desirable to r e t a i n the comparative computational It i s simplicity of the unconstrained time-varying gain case, as f a r as this i s possible. This requires developing a u n i f i e d treatment f o r the various formulations and then making maximum use of the known results f o r the l i n e a r state regulator problem with time-varying gain and quadratic cost. Substituting the R i c c a t i transformation p = -Sy (5.20) into (5.10) and using (5.6) y i e l d s ' S* + SC + C'S = -Q - KFK, S = 0 (5.21) Consider the general case of a time-varying system with time-varying gains. In the l i n e a r formulation I(K^, K, t ^ ) , E = 0. Substituting (5.20) into (5.18) y i e l d s t / f (F(S - S)yy' + yy'(S - S)F)dt = 0 (5.22) Equation (5.22) i s s a t i s f i e d when S = S and (5.19) i s s a t i s f i e d when K - S =• S. Substituting this- constraint into (5.21)- gives- the- standard matrix R i c c a t i d i f f e r e n t i a l equation f o r the optimum time-varying gain. In the case of a constant gain matrix the above procedure must be modified. H e u r i s t i c considerations indicate that some sort of time average of S should r e s u l t i n a good suboptimum constant gain (For a time-invariant system, the steady-state solution of (5.21), taking K = S, often gives a good suboptimum result). The correct time average to use i s given by (5.22). Replacing Sy by p, (5.22) can be solved f o r S by use of Kronecker products (See Appendix I I I ) . The vector representation (S) s a t i s f i e s the l i n e a r of the matrix S equation (5.23) A method of successive approximations could be considered for solving the two-point boundary value problem. K = K^, Using a nominal gain (5.6) i s integrated i n the forward d i r e c t i o n and backward d i r e c t i o n . Equation (5.16) i n the (5.23) can then be solved for S and (5.19) used to determine the gain increment. The updated gain i s (5.24) and the procedure i s repeated i t e r a t i v e l y . While successive approximations' seems straightforward, they generally f a i l to converge i n the case of nonlinear two-point boundary value problems. The s p e c i a l structure of the formulated problem, however, allows some insight to be concerning 5.4 the convergence of successive obtained approximations. Gradient Descent and Successive Approximations Let (6K) v and 3J (—) be the vector representations for a matrix 3K v gain increment and the gradient matrix, respectively. The incremental cost associated with 6K i s given by (5.25) where (5.26) i s p o s i t i v e semidefinite (see Appendix I I I ) . Consequently, i f 6K = w(S - K), 0 < w = 1 (5.27) then 5J = -w(S - K)'W(S - K) (5.28) = 0 The gain i s updated by \+l " \ +W ( (5.29) \"V where the subscript k i s used to indicate an i t e r a t i o n stage. For the problem formulation under discussion, successive approximations (5.24) i s a s p e c i a l case (w = 1) of. gradient descent (5.29) with a weighting matrix (5.26). Consequently, the proposed successive approximations method i s equivalent to a deflected gradient method. In the case of quadratic cost indices i t is.known that some deflected gradient methods (such as conjugate gradients) r e s u l t i n s i g n i f i c a n t l y improved rates. convergence Having established the property that (5.24) i s a deflected gradient descent, i t i s of i n t e r e s t to see whether i t has the properties of improved convergence, such as quadratic convergence. I t i s known that successive approximations r e s u l t s i n quadratic convergence i n the case of unconstrained time-varying gains. Quadratic convergence can be proven by considering (5.21) f o r two successive i t e r a t i o n s . k+1 Introducing S, , the l i n e a r matrix d i f f e r e n t i a l equation f o r AS, can be solved using standard procedures, giving -AS (0) o / f c k *k+l >° C ( T f f K + ± ( T , 0 \+l ~ ) C ( ) - V ^ C t , ^ - W ~ + ( \ " \+l ) F 0)dx + S " \+l " ( k W '° ) ] ( T ) d T (5.30) ' In (5.30) ^ ( t , x) i s the s t a t e - t r a n s i t i o n matrix of (5.6) at the k-th i t e r a t i o n stage. I t i s e a s i l y shown that the cost index at stage k i s given by J , = — y'S.y . k 2 o k o The incremental cost A J , = J , , , - J , i s found k k+1 k from (5.30) and can be expressed i n the form -AJ k - T r t l y y C k K K ^ - \)H\ ~ VI +± + T r K ^ ) ^ - K^) ] (5.31) where - C Wk+i \+i - V f [ ( F + F( \+i - V W k + i ^ (5 - 32) and where the notation $ (k) i s used to indicate that (5.12) i s evaluated yy at stage k. The gradient condition (5.19) can be obtained from (5.31) by noting that Tr[(A + A')B] = 0 f o r a r b i t r a r y symmetric matrix B only when A + A' = 0 . Equating (5.32) to zero y i e l d s (5.19). The remarkable properties of the case with unconstrained timevarying gain are seen from (5.32). By setting K ^ = S k+ k > (5.32) vanishes and the decrease i n cost (5.31) i s quadratic [38]. This i s accomplished by solving (5.21) .using successive approximations. The same quadratic convergence would r e s u l t i n the case of constant gain i f chosen to make (5.32) vanish. However, y^-^ depends on K^-^* condition ( ^ ) ^ = could be In p r i n c i p l e this may be possible. Consequently, any attempt to s a t i s f y the 0 would r e s u l t i n poor convergence. It i s reasonable to look f o r an approximation by replacing 87 y ^ ^ i n (5.32) by y^. Let ( ^ r ) ^ represent the matrix a r i s i n g from this + substitution. I t i s seen that <M>k • V k ) ( V-l " V F+ f ( \ + l " V yy $ where S^ i s defined by (5.22) evaluated at stage k. (5.33) that (^O^ = 0 i f (set w = 1 i n (5.29)). ( k > ( 5 i 3 3 ) I t follows from i s determined by successive approximations Since ( 7 ^ ) . i s an approximation to ( 7 ^ ) , i t is reasonable to anticipate that successive approximations r e s u l t s i n rapid convergence. A rigorous proof of quadratic convergence would require shoxtfing that the second term of (5.31) i s either never negative or, i f i t i s negative for a certain gain sequence, that i t s absolute value i s less than a fixed f r a c t i o n of the f i r s t term, which for Kj -^ f 1 c+ positive. i s always From a p r a c t i c a l point of view such a proof i s not e s s e n t i a l . If the choice w = 1 should happen to r e s u l t i n divergence, then a choice w < 1 can be made which w i l l always give a convergent sequence (see (5.28)). The proposed computational algorithm i s the following. a nominal K = K^. Integrate (5.6) i n the forward d i r e c t i o n and (5.10) i n the backward d i r e c t i o n . use of (5.29). Choose Solve (5.18) for S^ and update the gain by The rate of convergence i s controlled by w. Computationally, the proposed method i s s i m i l a r to the conventional steepest descent procedure. The. essential, difference i s i n the use. of a variable weighting matrix (5.26) to give a deflected gradient. The associated change i n the gain matrix i s actually given by the simple successive approximations formula (5.29). Limited computational experience to date indicates that this variable gradient d e f l e c t i o n allows a large step (w = 1) to be taken r e s u l t i n g i n rapid convergence. Note that the l i n e a r and nonlinear formulations of the s e n s i t i v i t y problem are treated i n a u n i f i e d manner. This arises from the fact that i n both formulations a nominal gain K i s chosen to perform the computations. The above algorithm i s s u i t a b l e for o f f - l i n e design studies. Since the optimum constant gain depends on i n i t i a l conditions, there i s a problem of p r a c t i c a l on-line implementation. i n i t i a l conditions and deriving suboptimal given i n . [ 3 9 ] . A method for eliminating feedback control laws i s A method which i s s u i t a b l e f o r combined on-line e s t i - mation and control i s discussed i n Section (5.6). 5.5 Optimum Constant Gain: Trajectory S e n s i t i v i t y Reduction Because of i t s engineering importance, open-loop and closed-loop system s e n s i t i v i t i e s have been extensively investigated (see Reference [32] for a bibliography). However, for the closed-loop s e n s i t i v i t y formulation given by (5.6) and unconstrained (5.8), a paradox arises i f time-varying gains are permitted. This paradox has been discussed i n the l i t e r a t u r e but i t s e f f e c t on s e n s i t i v i t y reduction does not seem to have been f u l l y appreciated [33, 37]. example can be- used to i l l u s t r a t e - the paradox. x = -ax + u, A simple Consider- x(0) = x (5.34) o z = -x - az - aKjZ, z(0) = 0 u = -Hi K - K z = 0 (5.36) 2 J = -| /°°(x + Q z 2 2 2 (5.35) + Ru )dt 2 (5.37) and suppose that and R are both chosen large compared to unity. This would seem a reasonable choice to make i n order to keep both s e n s i t i v i t y and control e f f o r t small. The optimum time-varying gains, however, s a t i s f y the constraint (5.36) and are i n f i n i t e . taking K -> + » i n (5.35). Equations 1 aK K 2 = K This can be seen by (5.34)-(5.37) then y i e l d . x = 0 (5.38) (5.39) 2 (5.40) The optimum ( i n f i n i t e ) gain can r e s u l t i n a very s i g n i f i c a n t smaller cost than that associated with other suboptimal gains (see [37]). (but p h y s i c a l l y realizable) I t i s evident, however, that minimization of (5.37) has not succeeded i n the design objective, which i s to decrease system sensitivity. open-loop. The system (5.34) with the optimum control (5.36) operates The same d i f f i c u l t y i s experienced with a general n-th order system (see (5.1), (5.4) and (5.5)). The s e n s i t i v i t y problem defined by (5.6) and (5.8) i s consequently i l l - p o s e d i f unconstrained time-varying gains are permitted. From both a computational and an engineering point of view the simplest constraint i s that of time-invariant gains. It i s i n t e r e s t i n g to note that this choice was made by Kreindler i n h i s design studies which dealt exclusively with l i n e a r time-invariant systems and t ^ = + <» [3]. I t seems natural to choose constant gains i n this case. However, i f t ^ i s f i n i t e and/or the system i s time-varying, then the choice of constant gains i s not obvious. 90 5.6 On-Line Computation of Gains From an engineering point of view a very a t t r a c t i v e f e a s i b l e on-line adaptive c o n t r o l scheme i s to i d e n t i f y a process as a l i n e a r time-invariant system of the form (5.1) over a f i n i t e observation i n t e r v a l and then apply a c o n t r o l of the form (5.3) over a subsequent f i n i t e control i n t e r v a l . In such a combined estimation and control scheme i t i s often important to reduce the system s e n s i t i v i t y to incorrect parameter estimates. I t i s , however, not generally computationally f e a s i b l e to solve a nonlinear two-point boundary value problem on-line. Some form of suboptimal control must be introduced. It i s known that the steady-state solution of the matrix R i c c a t i equation (take K = S i n (5.21)) often results i n a good suboptimal control. With this i n mind l e t P = -My where M and N are constant matrices. - Ny + e (5.41) Substituting (5.41) into (5.10) and choosing M to be defined by MC + C M = -Q - K F K , (5.42) which i s the steady - state form of (5.21) gives e + Ce = -(NC + CN)y, e f = (M + N)y (5.43) f A For reasons of n o t a t i o n a l convenience the l i n e a r formulation I(K^, K, i s discussed. t^) Substituting (5.41) into (5.11) y i e l d s F(S - M - N)$ yy + $ yy (S - M - N)F = -F$ ey - $ ye F (5.44) The decomposition (5.41) does not uniquely define N and e. Consider the p o s s i b i l i t y of defining N i n terms of observational data so that the r i g h t hand side of (5.44) has a n e g l i g i b l e e f f e c t on S. The updated gain, as given by (5.24), i s then S = K 2 = M + N (5.45) It i s seen from (5.43) that e i s "small" when N and e^ are "small". Besides the obvious trade-offs involved i n keeping N and e^ simultaneously "small", there i s also the requirement f o r on-line computation using observational data. The following cost index accounts f o r the various factors which enter into an "optimal" decomposition (5.41): Je = ^ / 2 o In (5.46) W^ and f (e + C'e)*W. (e + C'e)dt + ^ e/W.e. 1 2 r 2 f are p o s i t i v e d e f i n i t e weighting matrices. (5.46) Sustituting (5.43) into (5.46) i t i s seen that Je i s an algebraic function of N, $ and y^y ^' Setting N = L + L', where L i s a r b i t r a r y , allows matrices 1 to be used to determine the minimum. 1^oL gradient This y i e l d s = TC + CT + W (M + N)i|) + iK(M + N)W_ = 0 2 t r 2 1 0 (5.47) where T = W S$ + $ SW, 1 yy yy 1 1 S = NC + C'N <J> = y y ^ f f Equation (5.47) can be solved f o r N by use of Kronecker products. seen that (5.48) It is (S) (T) v v = G. (N) 1 v = G (S) 2 (TC* + CT) v [ G 1 2 3 G G + 4 G ] ( N ) v v = G.(T) 3 v ~ 4 = G ( M ) v ( 5 * 4 9 ) where G 1 = C 2 G 4 , © I + I © C ' =G yy 1 = <|» ® W f 2 1 + W O ^ 2 3 yy (5.50) f E f f i c i e n t algorithms are available f o r evaluating (5.50) [40]. The solution of (5.49) f o r N involves one matrix inversion. The system parameters... in. the. C. matrix- are- determined, by the i d e n t i f i e r from, ob.servational data. The matrices $ and i i , can be found by a forward yy f integration of (5.6). A l t e r n a t i v e l y , i t may be possible to use measured states or state estimates to evaluate these matrices. The gain compu- tation (5.45) i s then an algebraic problem and replaces the two-point boundary value problem. On-line use of (5.42), (5.49) and (5.45) to update the gain i s computationally equivalent to solving these equations by successive approximations. optimum gain. The f i r s t i t e r a t i o n must result i n an improved subI t i s desirable, but no e s s e n t i a l , that successive i t e r a - tions r e s u l t i n a convergent Since and W 2 suboptimum sequence. are a r b i t r a r y p o s i t i v e d e f i n i t e weighting matrices, i t i s not possible to reach general conclusions concerning the convergence of successive approximations. I t i s easy to show, how- ever, that given a nominal K, weighting matrices e x i s t which r e s u l t Suppose that the optimum gain K* i s determined i n one-step convergence. by o f f - l i n e computations. Equation (5.42) can be solved for M. (see (5.45)) N = K* - M and take W^ = I. Equation Let (5.47) i s a Liapunov type equation which can be solved for a p o s i t i v e d e f i n i t e l ^ . With the weighting matrices so defined there i s one-step convergence to the o optimum gain. It should therefore be possible to choose weighting, matrices so that the f i r s t i t e r a t i o n r e s u l t s i n an improved suboptimum gain for small v a r i a t i o n s about nominal parameter values. I f the system i s time-invariant and i f N i s small i n the sense that terms involving products of the elements of N can be neglected, then i t i s seen by subs t i t u t i n g (5.45) into (5.42) that M i s the steady-state gain. In this case (5.42) can be solved for the steady-state gain and N evaluated by (5.44). 5.7 Examples To give a n o n - t r i v i a l i l l u s t r a t i o n of the successive approximations as w e l l as the approximation method discussed i n Section (5.6) t ^ should be chosen smaller than the s e t t l i n g time of the open-loop response for the nominal gain, and the nominal gain should d i f f e r s i g n i f i c a n t l y from the optimum gain. These requirements are s a t i s f i e d i n the case of the second order system C x l 2 = " 2 qX = x + U ' x , x l(°) = (5.51) 1 , 0 x (0) = 2 -1.0 t J o = \ f I o f (x 2 1 + 3x 2 + u )dt 2 2. 94 by taking q = 1, t f = 2.0, and = (-10, 10). (The system i s an o s c i l l a t o r for u = 0 and i s unstable f o r the chosen nominal gain). Successive approximations ((5.6), (5.10), (5.23) and (5.24)) required eight i t e r a t i o n s to converge. No attempt was made to determine weighting matrices which would r e s u l t i n rapid convergence for the on-line approximation method ((5.42), (5.49) and (5.45)). = W 2 = I was made. An a r b i t r a r y choice Six i t e r a t i o n s were required. The r e s u l t s are the following: 1. Successive approximations K^** = (-1.06, - .215), 2. J ** = 1.518 On-line approximation method. K ** = (-1.713, - .952), J * * = 1.565 q (5.52) The minimum cost for unconstrained time-varying gain i s J * = 1.503. It i s seen that both methods give very good suboptimal constant gains. The increases are only 1% and 4% for case (1) and (2), respectively. To i l l u s t r a t e the s e n s i t i v i t y problem consider (5.51) to be augmented by (5.5) and take J = 5 f o tf ( z + z )dt + J 1 2 o 2 (5.53) 2 A comparison i s made with Lee's et a l . method which uses a time-varying gain and the choice t ^ = 10 i s made. The following results are obtained (a) Lee's method [35] (augmented state feedback; t = t ) u = -2x - 1.5x + .25z - .5z , J * = 1.94 ± 2 1 2 (5.54) (b) Constant gain (augmented state feedback) u = -6.8x (c) i - 6.6x + 4.66 2 Z;L + 20z 2> J * = 1.70 (5.55) Constant gain (Non-augmented state feedback) u = -1.812X - 1.22x , J = 2.04 (5.56) 2 Figures (5.1)-(5.3) i l l u s t r a t e the improvements that result byuse of constant gain and augmented state feedback. The magnitude of the control i s smaller, the cost index i s smaller and there i s a s i g n i f i c a n t reduction i n trajectory s e n s i t i v i t y . That the design objective of t r a - jectory s e n s i t i v i t y reduction i s achieved can be seen by comparing the s e n s i t i v i t y functions of case (b) and case (c) (case (c) i s based on (5.51) taking t f = 10). The optimality condition (5.17) was .0004). checked and gave (-.0002, I t should be noted that Lee's method uses a time-varying gain which i s a suboptimum solution of an i l l - p o s e d s e n s i t i v i t y problem. (The gain values given i n (5.54) are for t = t ). It i s therefore not surprising that an optimum constant gain can give s i g n i f i c a n t l y better results. Note that increasing t ^ from 2 to 10 and using s e n s i t i v i t y augmentation does not result i n an appreciable increase i n the cost index (compare (5.55) and (5.52). The convergence of successive approximation (5.29) was gated by choosing K ± as the i n i t i a l gain. optimum investi- = (-1.5, -1.0), K 2 = (-0.1, -0.5) This choice i s s i g n i f i c a n t l y different gain (see (5.55)). (5.57) from the Only f i v e i t e r a t i o n s were required to converge to the optimum. The following example i s used to i l l u s t r a t e the o n - l i n e approximation method when i t i s applied to a r e a l i s t i c system of moderate complexity. Figure (5.4) i l l u s t r a t e s a block diagram for the feedback control of a power generator connected to an i n f i n i t e bus [41]. The matrices i n (5.1) are 0.0 -0.676 0.0 0.0 0.0 25.0 0.0 0.0 1.0 -.25 A= , 0.0 0.0 0.0 -0.06 (5.58) B= 2.0 0.0 0.0 -2.0 2.0 -2.0 The components of the state vector represent the following v a r i a b l e s ; x^ = power angle deviation, x 2 = angular frequency deviation, x^ = mechanical power deviation, x^ = governor p o s i t i o n v a r i a t i o n , u = governor p o s i t i o n c o n t r o l . 3.0, The i n i t i a l state i s taken to be x^ = ( 0 . , 2 . 0 , 0.0) and i n (5.2) and (5.46) the choice Q = W = W = I and ± R = 1 i s made. 2 The choice t^ = 1 sec. i s made to accentuate convergence problems, i n case they should e x i s t . The minimum cost using the optimum time-varying gain S(t) (5.21)) i s J * = 33.215. (see Successive approximations for the optimum con- stant gain i s i n i t i a l i z e d using K = S(0). The results for the cases (a) K = M, (b) K = M + N, are as follows. k = 1 2 3 4 (a) J = 45.404 39.167 36.636 36.121 (b) J = 38.035 36.267 35.601 34.682 (a) u* = -.66x (b) u* = -.069x ± (5.59) - 1.22x -8.44x -3.35x 2 3 4 -1.47x -10.67x -5.47x 2 3 4 (5.60) The suboptimal choice K = M i n case (a) amounts to taking the solution of the i n f i n i t e time case f o r the constant gain. I t i s known that (5.42) converges quadratically to the optimum gain f o r the i n f i n i t e time case. The suboptimum choice K = M + N i n case (b) i s based on the approximation technique discussed i n Section (5.6). It i s seen to have a quadratic-type convergence and r e s u l t s i n a smaller cost index than f o r case (a). In fact, the f i r s t i t e r a t i o n results i n a s i g n i f i c a n t decrease. A one-step improvement of this kind could be of p r a c t i c a l s i g n i f i c a n c e i n an on-line a p p l i c a t i o n . a FIG.(5.1) EXAMPLE I CONTROL AS A FUNCTION OF TIME FOR: (a) LEE'S METHOD, (b) CONSTANT GAIN AND AUGMENTED STATES, (c) CONSTANT GAIN AND NON- AUGMENTED STATES. 0.24 0.161 ^ 0.081 0.0 to co-0.08 ^~ o UJ £ -0.24 FIG. (5.2) EXAMPLE 1. TRAJECTORY SENSITIVITY FUNCTION 3 FOR(a) LEE'S METHOD- (b) CONSTANT GAIN AND AUGMENTED STATES, (c) CONSTANT GAIN AND NON - AUGMENTED STATES. 7 VO VO FIG. (5.3) EXAMPLE I TRAJECTORY SENSITIVITY FUNCTION Z (a) LEE'S METHOD, (b) CONSTANT GAIN AND AUGMENTED STATES, (c) CONSTANT GAIN AND NON- AUGMENTED STATES. 2 F 0 R : SPEED GOVERNOR FIG. (5.4) TURBINE EXAMPLE 2. BLOCK DIAGRAM VOLTAGE CONTROL SYSTEM. GENERATOR OF A SYNCHRONOUS GENERATOR 6. 6.1 OPTIMAL ADAPTIVE CONTROL Introduction In this chapter we w i l l discuss one possible strategy for optimal adaptive control which i s p a r t i c u l a r l y a t t r a c t i v e point of view. from an engineering The strategy i s to i d e n t i f y system parameters at the end of an observation i n t e r v a l and then to use the parameters to derive an "optimal" control for a subsequent control i n t e r v a l . This chapter develops one version of this strategy which i s based on the i d e n t i f i e r s and optimum controls developed i n Chapters 2, and 6.2 5. Optimum Control During the F i r s t and Subsequent Optimization Intervals The p a r t i c u l a r optimization technique to be used during the f i r s t optimization i n t e r v a l i s of importance since the i n i t i a l information about the system i s usually i n s u f f i c i e n t to derive a good c o n t r o l . parameter estimates are poor then the optimum constant gains may s u f f i c i e n t l y close to the correct values. If not be The problem of convergence can become serious during the f i r s t optimization i n t e r v a l . Hence, i t i s re- commended to use the successive approximations method with a step size constant w < 1 (see equation (5.27)). Figures (6.1-a, b, c, d) show the time intervals associated with a recursive observation-optimization strategy which would be suitable for the optimal adaptive systems. operating at time t (6.1-a, The system i s considered to s t a r t and to terminate i t s operation at time t ^ . d) (T) i s the length of the sampling i n t e r v a l . The In Figures length of each observation i n t e r v a l i s (NT) seconds and the s t a r t of the i - t h INITIALIZATION INTERVAL OPERATION R PERIOD -sa-f«a J t -NT FIG. (6.1-a) 0 L to THE INITIALIZATION INTERVAL AND THE OPERATION PERIOD. rT rT t -NT 0 FIG. (6.1-b) FlG.(6.1-c) t tj° 0 NT t°+NT THE SUCCESSIVE INTERVALS. OBSERVATION THE SUCCESSIVE OPTIMIZATION INTERVALS (l = r). IT FIG. (6.1-d) t , T c THE SUCCESSIVE OPTIMIZATION INTERVALS (l = 2r). observation i n t e r v a l i s denoted by the instant ( t ^ ) . The displacement between two successive observation i n t e r v a l s i s (rT) seconds (see F i g . (6.1-b)) where r i s an integer number. The length of each optimization i n t e r v a l i s (T ) seconds and the s t a r t of the i - t h optimization i n t e r v a l i s denoted by the instant (t^) (see Fig's (6.1-c, d ) ) . The displacement between two successive optimization i n t e r v a l s i s (&T) seconds. meter I i s a multiple of r , i . e . , I = mr. The para-, The value of m need not be the same f o r every two successive optimization i n t e r v a l s . The choice of the parameter (m) i s based on the rate of v a r i a t i o n of certain perfor- mance indices to be subsequently discussed (see subsection(6.3.4)). During the i n i t i a l i z a t i o n i n t e r v a l (t -NT < t < t ) s p e c i a l input test Q signals can be used to drive the system. Measurements of the input and output at successive sampling instants are used to determine i n i t i a l estimates of the unknown parameters of the system. The choice of the i n i t i a l conditions at the s t a r t of the i n i t i a l i z a t i o n i n t e r v a l (t = t -NT) o i s completely a r b i t r a r y . At the end of each observation i n t e r v a l new estimates of the plant parameters are computed. At the beginning of each optimization i n t e r v a l new estimates of the "optimum" constant gains are made using the l a t e s t estimates of the plant unknown parameters. These new gains are used f o r the coming control i n t e r v a l of IT seconds. 6.3 Examples In this section a l i n e a r , continuous system with time-invariant parameters i s considered. Some of the system parameters are unknown and an optimal adaptive control i s computed f o r the system. The recursive observation-optimization strategy described i n the previous section i s followed. The i d e n t i f i c a t i o n of the unknown parameters i s performed using the d i g i t a l i d e n t i f i e r developed i n Chapter 2. The number of augmented sets of states for Examples 1 and 2 are one and two, respectively. The test inputs used during the i n i t i a l i z a t i o n i n t e r v a l are taken to be step functions with unity magnitude, i n both examples. at The i n i t i a l conditions the s t a r t of the i n i t i a l i z a t i o n i n t e r v a l are taken to be i d e n t i c a l l y zero i n both Example 1 and 2. It should be noted that use of a d i g i t a l i d e n t i f i e r r e s u l t s i n quantization noise. Consequently, even though the examples treat deter- m i n i s t i c cases, the system models are perturbed by quantization noise. This noise also a f f e c t s the value of the optimum constant which are computed using the i d e n t i f i e d parameters. noted that the optimum constant It should be further gains depend on the system state at the s t a r t of an optimization i n t e r v a l . of" the proposed" adaptive"strategy 6.3.2 feedback gains The examples i l l u s t r a t e convergence i n the presence of quantization noise. Example 1 The system considered i n this example i s represented by the dynamic model l X x 2 = a l x 2 + U =- a^x^ (6-. 1) where the nominal values for a^ and a^ are -1.0 The i n i t i a l conditions are x^CO) = ^(0) of the recursive observation-optimization = 1.0. and 1.0, respectively. The e s s e n t i a l parameters strategy (see Figs. (6.1-a, d)) are taken to be: T = 0.2, T t , = 10., f = NT = 2.0 t = 0., r = I = 2, N = 10 o (6.2) It i s required to compute the optimum constant feedback gains and such that the performance index 2+t? J. = \ f l (x + x + u )dt 1 L 1 2 • C L 2 (6.3) 2 I i s minimized using the new estimates of a^ and a^ which are obtained at the end of the i - th observation i n t e r v a l . The optimum constant gains are computed using the successive approximations method developed i n Section (5.2) with w = 1. The optimum solution i s assumed to be obtained when the inequality ' i " i J J _ 1 | = 1 0 " ' 5 i s s a t i s f i e d , where J . i s the value of J . at the end of the k-th l ( 6 , 4 ) iteration. l The time t ^ (see Fig.'s (6.1-c, d)) denotes the beginning of the i - t h , optimization i n t e r v a l . During the f i r s t optimization i n t e r v a l the i n i t i a l guess f o r the feedback gains and i s taken to be = Y.^ = -10. For the subsequent optimization i n t e r v a l s , the i n i t i a l guess f o r the i - t h optimization i n t e r v a l i s taken to be the optimum values obtained f o r the ( i - l ) - t h optimization i n t e r v a l . after s i x i t e r a t i o n s The stopping rule (6.4) was s a t i s f i e d f o r the f i r s t optimization i n t e r v a l . No more than two i t e r a t i o n s were required to s a t i s f y the stopping rule given by (6^4) for the subsequent optimization i n t e r v a l s . The performance index of the optimal adaptive system i s taken to be J* = • a 25 E AJ* i=l (6.5) 1 whe re the optimum value of J ^ (see (6.3)) i s denoted by J * and AJ* i s i i the value of J * f o r an integration i n t e r v a l t ^ < t < t ^ + 0.4. Since the optimum constant gains computed at the s t a r t of the i - t h optimization i n t e r v a l are only used f o r a period of (JIT) seconds s t a r t i n g from the c instant t . , then J * as given by (6.5) i s the t o t a l operation cost. X For cL the sake of comparison the optimum value of the performance index J o = \ f I o (x . + xl + u ) d t 1 l 2 (6.6) 2 i s computed using (6.1) with a^ and a^ being equal to t h e i r respective nominal values and the control feedback gains are taken to be the optimal time-varying gains. The optimum value of J q i s found to be 2.046 while the value of J * (see (6.5)) i s equal to 2.202, i . e . , the increase i n the performance index value does not exceed 8%. Figures (6.2) and (6.3) show that the i d e n t i f i c a t i o n errors do not exceed 5%. 6.3.3 Example 2 The system considered i n this example i s represented by the dynamical model x^ = a^c^ + 2u x 2 = a x 3 + 2u x 3 = a x 2 - u 2 3 where the nominal values f o r a^, a , and a 2 respectively. 3 (6.7) are -1.0, 2.0, and -2.0, The i n i t i a l conditions are x^(0) = x ( 0 ) = ^ ( 0 ) = 1.0. x 2 The e s s e n t i a l parameters of the recursive observation-optimization strategy are as given by (6.2). gainsK , K 2 and K 3 It i s required to compute the optimum constant such that the performance index 108 2+t° J.=|J 1 (x 2 + x 2 + 4 u )dt (6.8) 2 + i i s minimized during the i - t h optimization i n t e r v a l . The optimum solution i s assumed to be obtained when the inequality (6.4) i s s a t i s f i e d . In this example, using the successive approximations method with w = 1 (see equation (5.27)) during the f i r s t optimization i n t e r v a l r e s u l t s i n a diverging sequence of i t e r a t i o n s . Thus w i s taken to be and the i n i t i a l guess f o r the feedback gains are taken to be = -2, and = -3. rule (see (6.4)). 0.1 = -1, Twenty i t e r a t i o n s are required to s a t i s f y the stopping For the subsequent optimization i n t e r v a l s , the i n i t i a l guess f o r the i - t h optimization i n t e r v a l i s taken to be the optimum values obtained f o r the ( i - l ) - t h optimization i n t e r v a l . No more than two iterations, are required to satisfy. (6.4) for the subsequent optimization i n t e r v a l s and w i s taken equal to unity. For the sake of comparison the optimum value of the performance index J 1 10 2 2 2 2 = i / (xf + x^ + x^ + u )dt o 2 o 1 2 3 i s computed using (6.7) with a^, a^, (see (6.5)) was The found to be 0.655 while the value of J * a equal to 0.705, that i s , the increase i n the performance index value does not exceed 7.5%. Figures (6.6), (6.7), and (6.8) show the i d e n t i f i c a t i o n results for the unknown parameters a^, a^, and a^, respectively. « and a^ being equal to their respective nominal values and the control feedback gains are time-varying. optimum value, J * . of J was ^ ' o o (6.9) 1 U The maximum i d e n t i f i c a t i o n errors do not exceed 5%. The rapid rate of i d e n t i f i c a t i o n i s evident i n both this example and the previous one since a l l the unknown parameters have been i d e n t i f i e d 109 to a good degree of accuracy a f t e r only one observation i n t e r v a l . In Figures (6.4), and (6.5), the piece-wise constant curves which are drawn as s o l i d l i n e s represent the optimum constant gain K* and K*, r e s p e c t i v e l y . The parameters a^, and (see (6.1)) are updated at the end of each observation i n t e r v a l and a new estimate of the feedback gains i s made at the s t a r t of the following optimization i n t e r v a l . Since the optimum values of the constant feedback gains depend on the plant state at the s t a r t of the optimization i n t e r v a l , these optimum values are not constant during the whole operation period (see Figure (6.1-a)). The curves drawn as dashed l i n e s i n Figures (6.4), represent the optimum time-varying feedback gains K^, and and (6.5) respectively. These gains minimize the performance index given by equation (6.6) subject to the dynamical constraints given by (6.1) with a^ and a^ being equal to t h e i r respective nominal values. Figures (6.9), K*, K* and K*, (6.10), and (6.11) show the optimum feedback gains respectively. the optimum constant gains. The curves drawn as s o l i d l i n e s represent The parameters a^, a^, and a^ (see (6.7)) are updated at the end of each observation i n t e r v a l and a new estimate of the feedback gains i s made at the s t a r t of the following optimization i n t e r v a l . The optimum constant gains are not constant along the whole operation period but unlike the s i t u a t i o n i n Example 1 the feedback gains for Example 2 (see Figures (6.9), (6.10), (6.11)) change less The curves draxm as dashed l i n e s i n Figures (6.9), frequently. (6.10), and (6.11) represent the optimum time-varying feedback gains K^, K^, and These gains minimize the performance index given by equation (6.9) respectively. subject to the dynamical constraints given by (6.7) with a^, a^, and a^ being equal to t h e i r respective nominal values. 110 The recursive observation-optimization Section (6.2) strategy described i n i s characterized by i t s f l e x i b i l i t y . due to the fact that i t This f l e x i b i l i t y i s possible to adjust the displacement distance ( £ T ) between the successive optimization intervals c, d)) is (see Figures (6.1-b, according to the behaviour of the system under consideration. The shortest time i n t e r v a l during which the gains remain constant i s equal to (2rT) i n Example 1 (see Figures Example 2 (see Figures ( 6 . 9 ) , (6.4), (6.10), 4r for Examples 1 and 2, respectively, with less computational e f f o r t . (6.5)) and equal to (4rT) (6.11)). in Hence, using I = 2r and results i n a s a t i s f a c t o r y performance The value of % used i n both Example 1 and 2 was % = r. In the following discussion two performance indices w i l l be discussed. not i t These performance indices can be used to decide whether or i s necessary to update the i d e n t i f i e d parameters and/or the con- stant feedback gains at the end of each observation 6.3.4 interval. Discussion The following performance index can be used to reach decisions concerning the necessity for updating model parameters: NT 1 2 o p m p m In Eqn. (6.10) Q i s a p o s i t i v e d e f i n i t e weighting matrix, x p i s the measured plant state during the l a t e s t observation i n t e r v a l , and x m i s the corresponding model state during the same i n t e r v a l using the l a t e s t i d e n t i f i e d values of the plant unknown parameters. The following performance index can be used together with given by (6.10) to reach decisions concerning the necessity for updating feedback gains: x'(t?)x (t?) - x * ( t ..)x (t? C j In Eqn. = P 1 P P 1 (6.11), x (t^) I s t n e i - t h optimization i n t e r v a l (the present X - X P ^ ) (6.11) plant state at the s t a r t of the i n t e r v a l ) , and x (t^- -j) i s the plant state at the s t a r t of the ( i - l ) - t h optimization i n t e r v a l . The assignment of values for the upper bounds of and J^y respectively, would depend upon the degree of f a m i l i a r i t y with the p a r t i cular system under i n v e s t i g a t i o n (as can be seen from Figures (6.4), (6.5), (6.9), (6.10), and (6.11)). choice for the system considered Taking IT = 2rT would be a reasonable i n Example 1. £T = 4rT would be a reasonable choice. be performed whenever or For Example 2 taking Parameter or gain updating would exceeded t h e i r respective bounds. 0-4. 0.2 0-0. 1 2 3 4 5 6 7 8 9 W t 0.2 'I 0.4 . -0.5. -0.8 . -10. FIG. (6.2) EXAMPLE 1. IDENTIFICATION OF PARAMETER QUANTIZATION NOISE. o, WITH 1.0. 0.8. 0.6 . 0-4V 2 0.2 0.0 10 0.2 0.4 FIG.(6.3) EXAMPLE I. IDENTIFICATION OF PARAMETER QUANTIZATION NOISE. a WITH 2 t OPTIMAL SOLUTION ADAPTIVE SOLUTION 0.2 0.0 -0.2\ i i 1 1 1 1 2 3 4 5 6 — 1— I 1 7 8 r9 -0-4 «1 -0.6 -0-0 -1.0 -1.2 -1.4 FIG. (6.4) EXAMPLE 1. THE OPTIMUM PIECEWISE CONSTANT FEEDBACK GAIN -0.61 FIG.(5.5) EXAMPLE J. THE OPTIMUM GAIN K . 2 PIECEWISE CONSTANT FEEDBACK -0.2 I -0.4 -OS -as -t.o FIO.(SJS) EXAMPLE 2. IDENTIFICATION OF PARAMETER a, WITH QUANTIZATION NOISE 2.0 IS 1.2 08 0-4 10 -0-4 -0-8 FI61S-7) EXAMPLE 2. IDENTIFICATION QUANTIZATION NOISE. OF PARAMETER aj WITH as. 0-4. ao. -0.4 "3 -as. -1-2 . -IS -2X1. FIGIS-S) EXAMPLE 2. IDENTIFICATION QUANTIZATION NOISE. OF PARAMETER o, WITH I OPTIMAL ADAPTIVE SOLTN . SOLTN 02 0.0 -0.1 -0.4 FIG. (5.9) EXAMPLE CAIN K 2. THE OPTIMUM PIECEWISE CONSTANT FEEDBACK t OPTIMAL ADAPTIVE 0.4 SOLTN — SOLTN — — 0.2 0.0 I 2 3 4 S S 7 S 9 -0.2 -0-4 -OS -08 -t-0 FIG.t6.IO) EXAMPLE GAIN Kj 2. THE OPTIMUM PIECEWISE CONSTANT FEEDBACK 1.0. OS OS 0.4 3 0.2. OO -0.2 10 t OPTIMAL ADAPTIVE SOLTN SOLTN - - -0.4 FIG.(S-U) EXAMPLE 2. THE OPTIMUM GAIN Kj. PIECEWISE CONSTANT FEEDBACK 7. CONCLUSIONS In Chapter 2 the generalized equation error has been shown to be applicable to a mean-square method of rapid d i g i t a l estimation of l i n e a r system parameters. Due to the imposed structure of the estimator the manipulation of high order matrices i s avoided. by the estimator have a block diagonal structure. The matrices used Consequently only the lower order matrix within a block need be considered. It has been shown that the matrix can be recursively evaluated and that matrix inversions are not required. Examples i l l u s t r a t e the effectiveness of the estimator f o r a v a r i e t y of cases dealing with quantization noise as w e l l as measurements noise. In Chapter 3 a simple algorithm f o r computing the generalized inverse of a matrix has been developed. The advantage of this algorithm over other methods i s that i t eliminates the need for Gram-Schmidt orthogonalization and i t s associated ( i n the interest of accuracy) a l i z a t i o n as new vectors are introduced. reorthogon- Inherent i n the method i s a very simple check on the l i n e a r dependence of the d i f f e r e n t matrix columns. In Chapter 4 a two-stage method f o r estimating, time-invariant and time-varying parameters i n l i n e a r systems has been developed. During the second stage, the time-varying parameters are considered as unknown control inputs to a l i n e a r subsystem of known dynamics. This method f o r estimating time-varying parameters i s computationally simple i n that the estimates are obtained by the solution of a l i n e a r two-point boundaryvalue problem. Standard methods of solution r e s u l t i n one-step convergence A method of stage and control augmentation i n conjunction with d i s c r e t i - zation and mean-square error minimization has been shown to be e f f e c t i v e for the f i r s t stage of estimation. This approach reduces undesirable coupling between the two-stages which could r e s u l t i n convergence problems. Numerous examples i l l u s t r a t e the effectiveness of the method. In Chapter 5 method of successive approximations f o r solving the two point boundary-value problem f o r optimum constant gain matrices has been developed. The method i s shown to be computationally equivalent to a def l e c t e d gradient method. Convergence can always be achieved by choice of a scalar step-size parameter. Rapid convergence i s important i n o f f - l i n e design studies when the e f f e c t of various weighting matrices on optimal control performance i s investigated. A close relationship be- tween successive approximations and the solution technique for the l i n e a r time invariant case with i n f i n i t e t ^ i s exploited so that i n many cases rapid convergence i s achieved. The problem of trajectory s e n s i t i v i t y reduction by augmented state feedback i s shown to be well-posed i f gain matrices are constrained to be constant. The problem of on-line implementation of a constant gain matrix whose elements depend on i n i t i a l conditions i s discussed. The simplest approach i s to take a weighted average of the computed optimum gain matrices. feedback control law. back control This r e s u l t s i n an average suboptimal l i n e a r A second approach i s to introduce nonlinear feed- [39]. An on-line approximate method i s developed which appears suitable for systems whose parameters must be i d e n t i f i e d . The two point boundary-value problem i s replaced by an algebraic problem, the solution of which gives a suboptimal constant gain. I t does not appear f e a s i b l e to give a formal proof of convergence or a proof of improved response. By proper choice of weighting matrices, however, i t does seem possible to achieve an improvement with the f i r s t i t e r a t i o n . Examples and para- meter values were chosen to d e l i b e r a t e l y accentuate situations where convergence problems could a r i s e . countered and rapid convergence was No d i f f i c u l t i e s , however, were enachieved i n a l l cases investigated. A natural development of the present research i s to augment the s e n s i t i v i t y vector with the s e n s i t i v i t y functions with respect to the i n i t i a l conditions. This w i l l reduce the s e n s i t i v i t y of the optimum constant gains with respect to the state of the plant at the s t a r t of the optimization i n t e r v a l . Accordingly, the optimum constant gains w i l l be updated less frequently. Further research i s needed concerning the "best" choice of the constrained structure of the feedback gains. One p o s s i b i l i t y has been investigated i n this thesis where the gains are taken to be invariant. Let time- A more general choice i s to use s p e c i f i c time-varying gains. the feedback gains be constrained by the s p e c i f i c structure K(t) = FG(t), where the time-varying matrix G(t) i s to be s p e c i f i e d before s t a r t i n g the optimization process. The time-invariant gain matrix F i s to be computed according to the algorithm developed problem i s how i n Chapter 5. to choose G(t) systematically such that the r e s u l t i n g time-varying gain matrix K(t) i s "optimum". A s p e c i f i c time-varying gain could result i n better r e s u l t s concerning s e n s i t i v i t y function reduction. The APPENDIX I INITIAL CONDITIONS LINEAR LEAST-SQUARES ESTIMATOR The values of the optimum constant feedback gains as computed by the algorithm developed i n Chapter 5 depend on the plant state at the s t a r t of the optimization i n t e r v a l . Accordingly, a good estimate of the plant i n i t i a l conditions i s required. In a noise free s i t u a t i o n , i t may be possible to obtain the i n i t i a l conditions by measuring the state of the plant. However, i f the measurements of the plant state are conta- minated by noise then an i n i t i a l conditions estimator i s required. possible choice of such an estimator i s the l i n e a r least-squares One one. This estimator does not require a p r i o r i knowledge about the noise s t a t i s t i c s and i s characterized by i t s simplicity. Considering-the- i-ttv optimization i n t e r v a l (see Figures- 6vl-a, b, c, d), the instant (t^) indicates the s t a r t of this i n t e r v a l as well as the completion of the i - t h observation i n t e r v a l . Let the plant be represented c c (for t ^ - NT < t < t^) by x = Ax + B ± where A^, and (1.1) u ± are the i d e n t i f i e d values of the plant parameteric matrices A, and B, r e s p e c t i v e l y . These values are obtained according to the measurements taken during the i ~ t h observation i n t e r v a l . term on the right hand side of equation (1.1) B.u = -F.K 1 1 i jX can be expressed The second as (1.2) where the d e f i n i t i o n of the matrix F. i s as given i n Chapter 5 (F = BR "*"B') and represent the piecewise constant optimum feedback gain matrix being used during the i - t h observation i n t e r v a l . Substituting (1.2) into (1.1) y i e l d s x = (A. - F.K.)x i 1 1 (1.3) The t r a n s i t i o n matrix i^(t^, x) of the dynamic system represented by (1.3) can be obtained by solving the following matrix d i f f e r e n t i a l equation <L = -(A. - F .K. ) ' ii; r \ l l i K t ? , t?) = I T 1 T x (1.4) 1 Let X j , represent the noisy measurements of the plant-state at /s C t = t ^ - KT, and l e t x^ be an estimate of x^ defined by *K *K = x ( t (1.5) i> where x(t^) i s the estimate of x(t^) which minimizes the quadratic error function e = (x -r 2 a - x ) ' (x a a - x ) a (1.6) where *K-*<^ The state vector x ' i " 1 ^ X a = (1.7) i s given by a b J x a = ij, x(t?) a l (1.8) where (1.9) N The value of x(t^) which minimizes the error function e given by (1.6) i s known to be x(t ) = C 1 W a faJ ' V xa ,a (1.10) Equation (I.10) gives the best estimate of the plant state at t = t ^ which i s the i n i t i a l condition for the i - t h optimization (the present interval). to compute new interval These i n i t i a l conditions are required i n order values for the optimum constant feedback gains. APPENDIX I I ANALYTICAL SOLUTION FOR THE LIAPUNOV MATRIX EQUATION [42] This algorithm i s characterized by i t s s i m p l i c i t y since no matrix inversion i s needed and the solution i s obtained i n a closed form and not i n a numerical or i t e r a t i v e way. Consider the matrix equation S'X + XS = -Q where S, X and Q are a l l of dimension (II.1) (n x n) Let A S, = k k-i Z X.S i=0 k (II.2) 1 where the parameters X ^ ( i = 0, n-l) are the c o e f f i c i e n t s of the c h r a c t e r i s t i c equation of the matrix S. These parameters can be computed e a s i l y using Faddeev's method. X. = - ~ Tr[SS. . ] , S. = SS. . + X.I, i = 1, 1 x i - l' l i - l x ' ' X Q = 1 S Q = I, S ' n = 0 (II.3) The l a s t r e l a t i o n i s used as a check. We may write down the solution to equation (II.1) as X = 1 2 n " 2 j j - i E (-l) R(S, -S) k=0 n K k-l (S, -S) E s ; ( - l ) *S £=0 k A k Q 1 % k % (II.4) 123 where: R(S, -S) = ( - l ) 2 n H (S) = n n T°- H (S) n A,1 A o A3 0 A2 A (II.5) 2 Ao Ai 1 0 0 (II.6) n The determinant H (S) i s known as the Hurwiz determinant, r A (S, - S ) k odd o K+1 _ l» 2 where G^, H (S) —7 G, (S) An k n f o r k even, i s the (- + l ) t h . k even (II.7) Cofactor of the f i r s t column of H (S). n Hence equation (II,. 4) can be written i n a s i m p l i f i e d form as 1 X = 2*^ n n (S) _ 1 E i=0 2 1 G (S)£,=0 E 2± I 2i-£ (-iyS\QS (II.8) n4 From (II.8) i t i s clear that the condition for a unique solution i s H (S) ^ 0 which- i s equivalent, to. the condition, that, the dynamical system. (X = SX) i s asymptotically stable. For the n o n - t r i v i a l case where the matrix Q i s symmetrical, the term 2 1 E (-D A=O I S;QS (II.9) 124 w i l l be symmetrical and since G„.(S) and H (S) are scalers, X i s also 2 i n J symmetrical and we need only to solve n(n + l)/2 entries of X instead 2 of n . Thus f o r this case we have X U 1 1 where = - r T - i V G (S)Z 2 H (S) 1=0 N L <-!>*<SjQS Si=0 % i s the upper t r i a n g l e of X and ( pQ 2i-X?]J s1 S l i s ) X % t r i e ' (11.10) U u PP e r triangle f o r (SlQS„. „•) f o r the case i = £, while f o r i ^ I i t i s a triangular matrix with i t s main diagonal equal to twice the main diagonal of (SlQS ) and every element ( i , j ) i s equal to the sum of the corresponding element ( i , j ) i n the upper t r i a n g l e of ( £ Q S element i n the lower t r i a n g l e . s 2 i _ ) P £ l u s the ( i + 1, j - 1) 125 APPENDIX I I I The trace properties and gradient matrices used i n deriving the gradient conditions are the following [43], Tr[AB] = Tr[BA], Tr[AB ] = Tr[BA ] f 1 Tr[AKB] = A'B' Tr[AK'B] = BA Tr[AKBK] = A'K'B' + B'K'A' dK Tr[AK'BK'] = BK'A + AK'B dK. ^ TrfAKBK'] = A'KB' + AKB (III.l) dK The following properties of Kronecker matrix products have been u t i l i z e d i n the algebraic manipulation of matrix equations [44], (AB) v = ( 1 0 A) ( B ) = (B' © I) ( A ) v v (A ® B) (C © D) = (AC) © (BD) (A + B) ® (C + D) = A ® C + A ® D + B ® C + B ® D If Aj,, i = 1, .. ., n, and u , j = 1, of A and B, respectively, then n, are the eigenvalues u^. are the eigenvalues of A @ B. Because of the p a r t i t i o n i n g (5.7) of F, can be found from the eigenvalues of F. equations: (III. 2) the eigenvalues of F Consider the following eigenvalue Fz. = X.z. x $ V. yy J x x = (111,3) U.V. i i It i s seen from the defining equations f o r F and $ that (see (5.5) and (5.12)) (B'z.)'R / o 1 (B'z.) = X z ' z. = 0 i (III.4) i 2 (v.'y) dt = u.v.'v. = 0 1 111 The matrix W i s therefore p o s i t i v e semi-definite. Since ("j^) = W(S - K)^ i t follows that the equality sign i n (5.28) holds only when the gradient condition i s s a t i s f i e d . 127 APPENDIX IV The confidence i n t e r v a l length deduced i n Subsection(2.7.3) (see (2.101)) i s based on the assumptions that the random vector (c - c) i s normally d i s t r i b u t e d and that the random vector V i s s t a t i s t i c a l l y r. independent of (c - c) and i s normally d i s t r i b u t e d . Also the random 2 variable V'V i s assumed to be distributed as % . The above mentioned assumptions can only be j u s t i f i e d i f the matrix 6^ (see (2.11)) i s deterministic. In the following we w i l l prove the above assumptions f o r the case where 6^ i s a deterministic quantity: Using (2.12) and (2.38) y i e l d s c - c = P 0'Z - c N NN M M = N N N P 6 ( 6 C + N ) " C = P 6^N (IV.1) N where, f o r convenience, the random vector N i s taken to be normally d i s t r i b u t e d with E[N] = 0 and E[N N'] = 1 . I t follows d i r e c t l y from (IV.1) that (c - c) i s a normally distributed random vector. Using (2.63) and(EV.l) y i e l d s . E[(c - c ) V ] = E[P 0^NN'(6 P 6^ - I ) ] N = P N N N N N N ~ =0 6 (9 P e N I } (IV.2) Thus the random vectors (c - c) and V are s t a t i s t i c a l l y independent. The vector V i s normally distributed as can d i r e c t l y be seen from (2.63) where 128 V = (0 N N N " P I ) S = (W N " ^ ( 2 ' 6 3 ) Since 8JJ and P ^ are deterministic matrices then V i s normally distributed as N . Using (2.63), (2.64), and (2.65) y i e l d s W - D*(W - I)N = N'(W N N = N'(W - I) N 2 N = N ( I - W )N (IV.3) r N Considering the orthogonal matrix <J>^ (see (2.66)), equation (IV.3) can be written as V'V = N$ $'(I - W )$ $'N V N V N N V V V N'$ (I - D ) ^ N N N = EyJ (IV.4) where y. i s a l i n e a r combination of a l l the elements of the random vector (the matrix (I - D^) i s a diagonal one with the main diagonal elements being equal to 1 or 0 (see (2.69))). Thus the random vector y^ i s d i s t r i b u t e d as the vector $' N. The vector $'N i s distributed as N N the random vector N which i s a d i r e c t result of Fisher's Theorem [17]. Accordingly, the random vector 2 2 Ey^ i s d i s t r i b u t e d as % . 0^. i s normally distributed l i k e N and The above proofs are based on a deterministic I f 0 ^ i s a random v a r i a b l e , i t i s no longer possible to prove that the random vectors (c - c) and V are normally d i s t r i b u t e d along with the assumption that V V i s d i s t r i b u t e d i s % . r There appears to be no mathe- matically tractable approach to derive confidence i n t e r v a l s i n the case 129 where 0^ contains quantization noise. Equation To r e t a i n the basic s i m p l i c i t y of (2.101) the following procedure i s suggested. In equations (IV.1 - 4 ) , i t i s possible to replace Z ( k ) by Z ( k ) where ^ ( k ) i s N obtained by use of the best estimate (see N N c at stage k, that i s by use of (2.10)). z N =? ; dv.5) N The "smoothed" values 0^ of 0^ are used i n equations (IV.1-4). Since the randomness i n 6^ tends to be eliminated, i t seems j u s t i f i a b l e to consider 0^ as a deterministic quantity. APPENDIX V The gradient condition as given by equation (5.11) can be derived as follows: Let H denotes the Hamiltonian, then using (5.6) and (5.9) yields H = p'(A - F(K + K'))y - | y » ( Q + (K + K')F(K + K'))y (V.l) = p'Ay - p'FKy - p'FK'y - | y'Qy - | y'KFKy - j y'KFK'y - - | y'K'FKy - j y'K'FK'y Using ( I I I . l ) the Hamiltonian can be written as H = p'Ay - TrfFKyp'] - Tr[FK'yp'] - ~ y'Qy - |- Tr[yy'KFK] - | Tr[yy'KFK'] - | Tr[FKyy'K'] - |- T r [ y y ' K ' F K ' ] (V.2) The following gradient matrix operations are required g | Tr[FKyp'] = F(py') ^ Tr[FK'yp'] = (yp')F " \ ~k T r tyy' K f K ] = " f ( y y ' ) K ' F - \ FK' (yy') - \ g | Tr[yy'KFK'] = - ( y y ' ) K F - T r [ F K y y ' K ' ] = -FK(yy') " \ a! [ y y ' ' T r K f K '] ( = - \ ^ ' ( y y ' ) - f(yy')K'F The above operations allow ~ to be written as v . 3 ) 13 3H — 1 1 = - (py») - (yp')p - | ( « ) ' F - -| FK»(yy*) - (yy')KF F yy K - FK(yy') - | FK'(yy') - | (yy')K'F = - F(py') - (yp')F - (yy*)K'F - FK'(yy') - (yy')KF - FK(yy«) = - F(py') - (yp*)F - (yy')(K + K')F - F(K + K )(yy') , = - F(py') - (yp')F - (yy')KF - FK(yy*) (V.4) where K = K + K ? From (V.4) and f o r time-invariant F and K equation (V.4) gives t a H _ t £ - 9K = F / 0 t t t py'dt + / yp'dt F + / yy'dt K F + F K / yy'dt = 0 0 0 o (V.5) Using the d e f i n i t i o n A & = / yv t J yv'dt O Equation (V.5) can be written as t gxr _ _ _ / - -r^r = F<£> + $ F + $ KF + FK $ =0 9K py yp yy yy q This completes the proof of equation (5.11). (V.6) 132 REFERENCES 1. P. V. Kokotovic', J . B. Cruz, J r . , J . E. H e l l e r and P. Sannuti, "Synthesis of Optimally Sensitive Systems", Proc. I.E.E.E., Vol. 56, No. 8, August, 1968. 2. Lion, P., "Rapid I d e n t i f i c a t i o n of Linear and Nonlinear Systems", Joint Automatic Control Conference, August, 1966, pp. 605-615. 3. Pazde ra, J . , and Pottinger, H., "Linear System I d e n t i f i c a t i o n Via Liapunov Design Techniques", Joint Automatic Control Conference, August,1969, pp. 795-801. 4. Rose, R., and Lona, G., "An Approximate Steepest Descent Method f o r Parameter I d e n t i f i c a t i o n " , Joint Automatic Control Conference, August, 1969, pp. 483-487. 5. K i s h i , F. H., "On-Line Computer Control Techniques and Their A p p l i cation to Re-entry Aerospace Vehicle Control", C. T. Leondes, ed. Advances i n Control Systems Theory arid Applications, V o l . 1, Academic Press, New York, 1964. 6. C. P r i c e , "The Matrix Pseudo-Inverse and Minimal Variance SIAM Rev. 6, 1964, pp. 115-120. 7. L. D. Pyle, "Generalized Inverse Computations Using the Gradient Projection Method", Journal of ACM, 11, 1964, pp. 422-428. 8. B. Rust, W. R. Burrus, C. Schneeberger, "A Simple Algorithm for Computing the Generalized Inverse of a Matrix", Journal of ACM, 9, 1966, pp. 381-387. 9. Sage, A., and Melsa, J . , Systems I d e n t i f i c a t i o n , Academic Press, 1971. Estimates',' 10. H. C. Lessing and D. F. Crane, "The Use of Integral Transforms i n the Estimation of Time Variable Parameters", NASA Techn. Note, NASA TN D-5194, May 1969. 11. Sage, A. P., Optimum Systems Control, P r e n t i c e - H a l l , Inc., 1968, pp. 507. 12. Mayne, D. 0., "Optimal Non-Stationary Estimation of the Parameters of a Linear System with Gaussian Inputs", Journal of E l e c t r o n i c s and Control, V o l . 14, 1963, pp. 101--112. 13. Lee, R. C. K., Optimal Estimation and Control, M.I.T. Press, Cambridge, Mass., 1964. 14. Lee, R. C. K., and Ho, Y. C , " I d e n t i f i c a t i o n of Linear Dynamic Systems", Information arid Control, Vol. 8, No. 1, February, 1965. 15. Leathrum, J . , "On-Line Recursive Estimation", Joint Automatic Control Conference, August, 1969, pp. 173-178. 16. Bakke, R., "Case Study i n the Paper-Making i n System Control, 1968, pp. 113-192. 17. Y. V. Linnik, Method of Least Squares and the P r i n c i p l e s of the Theory of Observation, Pergamon, New York, 1961. 18. M. Aoki, Optimization of Stochastic Systems, Academic Press, 1967, pp. 318-324. 19. B. J . Prochaska, "Applications of Generalized Inverses to E s t i mation i n a Dynamic System", JACC, 1971, pp. 820-829. 20. D. Luenberger, Optimization by Vector Space Methods and Sons, 1969. 21. J . S. Meditch, "A Class of Suboptimal Linear Controls", Trans, on Automatic Control, July 1966, pp. 433-439. 22. A. Ben-Israel, A. Charnes, "Contributions to the Theory of General i z e d Inverses", SIAM Rev. 11, 1963, pp. 667-699. 23. R. Penrose, "A Generalized Inverse for Matrices", Proc. Cambridge Philosophical Society, 51, 1955, pp. 406-413. 24. D. K. Faddeev and V. N. Faddeeva, Computational Methods of Linear Algebra, W. H. Freeman and Co., 1963, 25. D. L. Kleinman and M. Athans, "The Design of Suboptimal Linear Timevarying Systems", I.E.E.E. Trans., Automatic Control, V o l . AC-13, pp. 150-159, A p r i l 1968. 26. D. L. Kleinman, T. Fortmann and M. Athans, "On the Design of Linear Systems with Piecewise Constant Feedback Gains", I.E.E.E. Trans., Automatic Control, V o l . AC-13, pp. 354-361, August 1968. 27. W. Levine and M. Athans, "On the Determination of the Optimal Constant. Output. Feedback..Gains for Linear Multivar.iab.le.. Systems", I.E.E.E. Trans., Automatic Control, V o l . AC-15, pp. 44-48, February 1970. 28. M. M. Newmann, "Specific. Optimal Control of the Linear Regulator Using a Dynamical Controller Based on the Minimal-Order Luenberger Observer", Intern. J . Control, V o l . 12, pp. 33-48, No. 1, 1970. 29. R. L. Kosut, "Suboptimal Control of Linear Time-Invariant Systems Subject to Control Structure Constraints", Proc. JACC, pp. 820-828, 1970. 30. T. Hendericks and H. D'angelo, "An Optimal Fixed Control Structure Design with Minimal S e n s i t i v i t y f o r a Large E l e c t r i c Booster", Proc. of the 5th Annual A l l e r t o n Conference, pp. 142-148, 1967. Industry", Case Studies 9 John Wiley 134 31. E. Kreindler, "On Minimization of the Trajectory S e n s i t i v i t y " , Intern. J . Control, Vol. 8, No. 1, pp. 89-96, 1968. 32. E. Kreindler, "Synthesis of F l i g h t Control Systems Subject to Vehicle Parameters Variations", Grumman A i r c r a f t Engg. Corp., Bethpage, N. Y., Techn. Rept. AFFDL-TR-66-209, A p r i l 1967. 33. E. Kreindler, "Formulation of the Minimum Trajectory S e n s i t i v i t y Problems", I.E.E.E. Trans., Automatic Control, Vol. AC-14, pp. 206-207, A p r i l 1969. 34. J . F. Cassidy J r . and I. Lee, "On the Optimal Feedback Control of a Large Launch Vehicle to Reduce Trajectory S e n s i t i v i t y " , Proc. of 1967 JACC, pp. 587-595. 35. H. J . Dougherty, I. Lee and P. M. DeRusso, "Synthesis of Optimal Feedback Control Systems Subject to Parameter Variations", 1967 Proc. JACC, pp. .125-133. 36. P. Sannuti and J . B. Cruz J r . , "A Note on Trajectory S e n s i t i v i t y of Optimal Control Systems", I.E.E.E. Trans., Automatic Control, Vol. AC-13, pp. 111-112, February 1968. 37. A. J . Bradt, " S e n s i t i v i t y Functions i n the Design of Optimal C o n t r o l l e r s " , I.E.E.E. Trans., Automatic Control, V o l . AC-13, pp. 110-111, February 1968. 38. D. L. Kleinman, "On an I t e r a t i v e Technique for R i c c a t i Equation Computations", I.E.E.E. Trans., Automatic Control, Vol. AC-13, . pp. 114-115, February 1968. 39. A. G. Longmuir and E. V. Bohn, "The Synthesis of Suboptimal Feedback Control Laws", I.E.E.E. Trans., Automatic Control, Vol. AC-12, pp. 755-758, December 1967. 40. C. F. Chen and L. S. Shiek, "A Note on Expanding PA + A'P = -Q", I.E.E.E. Trans., Automatic Control, V o l . AC-13, pp. 122-123, February 1968. 41. I i . K. Kirchmayer, Economic Control of Interconnected Systems , New York, John Wiley and Sons, 1959. 42. P. Mueller, "Solution of the Matrix Equations AX + XB = -Q and S'X + XS = -Q*", SIAM J . Appl. Math., Vol. 18, No. 3 pp. 682-687, May 1970. 43. M. Athans and F. C. Schweppe, "Gradient Matrices and Matrix Calculations", M.I.T. L i n c o l i n Lab., Lexington Mass. Techn. Note 1965-53, November 1965. 44. H. Neudecker, "A Note on Kronecker Matrix Products and Matrix Equation Systems", SIAM J . Appl. Math., V o l . 17, No. 3, pp. 603-606, May 1969.
- Library Home /
- Search Collections /
- Open Collections /
- Browse Collections /
- UBC Theses and Dissertations /
- Linear systems identification and optimization with...
Open Collections
UBC Theses and Dissertations
Featured Collection
UBC Theses and Dissertations
Linear systems identification and optimization with application to adaptive control Hanafy, Adel Abdel Raouf 1972
pdf
Page Metadata
Item Metadata
Title | Linear systems identification and optimization with application to adaptive control |
Creator |
Hanafy, Adel Abdel Raouf |
Publisher | University of British Columbia |
Date | 1972 |
Date Issued | 2011-03-21T18:30:58Z |
Description | This thesis is concerned with the problem of identifying and controlling linear continuous systems. Algorithms which are feasible for real-time digital computations are developed for the design of both the controller and the identifier. The generalized equation error is shown to be applicable to a mean-square method of rapid digital estimation of linear system parameters. Due to the imposed structure of the estimator the manipulation of high order matrices is avoided. Examples illustrate the effectiveness of the estimator for a variety of cases dealing with quantization noise as well as measurements noise. In some cases, this digital identifier requires the computation of the generalized inverse of a matrix. A simple algorithm for computing the generalized inverse of a matrix is developed. This algorithm eliminates the need for Gram-Schmidt orthogonalization and its associated (in the interest of accuracy) reorthogonalization as new vectors are introduced. A two-stage estimator is developed for estimating time-invariant and time-varying parameters in linear systems. During the second stage, the time-varying parameters are considered as unknown control inputs to a linear subsystem of known dynamics. For the first stage, the digital identifier is shown to be effective in identifying the time-invariant parameters. Numerous examples illustrate the effectiveness of the method. To design a feedback controller a method of successive approximations for solving the two point boundary value problem for optimum constant gain matrices is developed. The method is shown to be comptationally equivalent to a deflected gradient method. Convergence can always be achieved by choice of a scalar step-size parameter. An online approximate method is developed which appears suitable for systems whose parameters must be identified. The two point boundary value problem is replaced by an algebraic problem, the solution of which gives a sub-optimal constant gains. The problem of trajectory sensitivity reduction by augmented state feedback is shown to be well-posed if constrained structure of the gain matrices is taken. The simple structure of constant gain matrices is considered. Based on the developed identifier and controller, one possible strategy for optimal adaptive control which is particularly attractive from an engineering point of view is studied. The strategy is to identify system parameters at the end of an observation interval and then to use the parameters to derive an "optimal" control for a subsequent control interval. Two examples are considered in order to illustrate the effectiveness of this strategy. |
Genre |
Thesis/Dissertation |
Type |
Text |
Language | eng |
Collection |
Retrospective Theses and Dissertations, 1919-2007 |
Series | UBC Retrospective Theses Digitization Project |
Date Available | 2011-03-21 |
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.0101279 |
URI | http://hdl.handle.net/2429/32645 |
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 |
Campus |
UBCV |
Scholarly Level | Graduate |
Aggregated Source Repository | DSpace |
Download
- Media
- [if-you-see-this-DO-NOT-CLICK]
- UBC_1972_A1 H35.pdf [ 5.83MB ]
- Metadata
- JSON: 1.0101279.json
- JSON-LD: 1.0101279+ld.json
- RDF/XML (Pretty): 1.0101279.xml
- RDF/JSON: 1.0101279+rdf.json
- Turtle: 1.0101279+rdf-turtle.txt
- N-Triples: 1.0101279+rdf-ntriples.txt
- Original Record: 1.0101279 +original-record.json
- Full Text
- 1.0101279.txt
- Citation
- 1.0101279.ris
Full Text
Cite
Citation Scheme:
Usage Statistics
Country | Views | Downloads |
---|---|---|
China | 11 | 25 |
United States | 11 | 0 |
France | 4 | 0 |
Sweden | 3 | 0 |
Mexico | 2 | 0 |
Russia | 2 | 0 |
Canada | 1 | 0 |
Brazil | 1 | 0 |
Czech Republic | 1 | 0 |
City | Views | Downloads |
---|---|---|
Unknown | 10 | 10 |
Ashburn | 5 | 0 |
Beijing | 5 | 0 |
Harbin | 3 | 0 |
Shenzhen | 3 | 25 |
Stockholm | 3 | 0 |
Mountain View | 3 | 0 |
Salvador | 1 | 0 |
London | 1 | 0 |
Sunnyvale | 1 | 0 |
Mladá Boleslav | 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-0101279/manifest