Model Validation for Robust Control

Robert Andrew Davis

Churchill College Cambridge

A dissertation submitted to the University of Cambridge for the degree of Doctor of Philosophy

January 1995

To Helen and my parents

Preface

This work has been funded by the Engineering and Physical Sciences Research Council, to whom I am grateful. I would also like to thank my parents, the University of Cambridge Department of Engineering and Churchill College for additional financial assistance.

I am indebted to Professor Keith Glover for motivating the work in this thesis and providing me with constant inspiration. I have learnt a great deal from him, especially perseverance with the more difficult problems.

I am grateful to Giles Wood for supplying me with his model of the flexible beam, experimental data and much insight into the “real” system. For the linear models of the Harrier I must thank Rick Hyde, and for the PC implementation of the nonlinear model, Sanjay Lall. Rick Hyde was also instrumental, together with the Defence Research Agency at Bedford, in obtaining the flight test data. I am also grateful to the Ministry of Defence for sponsoring the Harrier VAAC Programme, which made the Harrier project possible.

I would like to thank the members of the control group, both past and present, for providing company during my PhD, especially Glen Vinnicombe, Markus Brandstetter and Phil Goddard. Glen has provided much input into section 5.3 and been an example of an excellent researcher. Markus has generously spent many hours reading preliminary versions of the thesis, and provided much helpful criticism. With Phil I have shared an office for the duration of my studies. I would like to thank him for his tolerance, friendship and never letting me forget what a wonderful place England is.

On a more personal level I would like to thank my parents, who have been a continuous source of love and support, and have provided the security of a dependable home. I would also like to thank Mark Gostick and Giles Pateman who have been, and continue to be, great friends.

Finally I must thank Helen, who has never lost faith in me and has given me endless patience, love and support.

As required by University Statute, I hereby declare that this dissertation is not substantially the same as any that I have submitted for a degree or diploma or other qualification at any other university. I further state that no part of this dissertation has already been, or is being currently submitted for any such degree, diploma or other qualification. This dissertation is the result of my own work and includes nothing which is the outcome of work done in collaboration.

Abstract

This thesis is concerned with validation of the models used for robust control. Theory is developed for validating the models used for robust control design and applied to data from a flexible beam and Harrier VSTOL aircraft.

Results from Carathéodory-Fejér interpolation theory are used to extend published results to a general class of models used for robust control. The computational tractability of these results is then examined using techniques from computational complexity. A large class of model validation problems are proved to be NP-hard, which means they are computationally at least as hard as a class of problems that are recognised to be computationally demanding. Conditions are also obtained for when the validity of a model can be tested by solving a convex feasibility problem.

Two sets of models, that are equivalent for robust control, are shown to be different for the purposes of validation. This is shown to be a consequence of the difference between certain balls defined in the gap metric. A study of the gap metric provides an interesting interpretation of model validation and motivates the validation of models defined as balls in the ν-gap metric. Recent results in the ν-gap metric are used to prove validation results for the largest set of models that are stabilizable by a certain set of robust controllers. These may be considered as the best possible validation results for a certain set of robust controllers.

The theoretical results are used to test the validity of models of a flexible beam and Harrier aircraft. It is shown that a model of the flexible beam cannot account for the observed data, but that a modified model can. Linear models of the Harrier are validated using data simulated by a nonlinear model. The same linear models are also validated using flight test data from the aircraft.

Keywords: model validation, robust control, H-infinity, nu-gap metric, flight test data, flexible beam

Contents

1 Introduction 1 1.1 Invalidation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3 1.2 Central themes . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4

1.2.1 Interpolation . . . . . . . . . . . . . . . . . . . . . . . . . . . 4 1.2.2 Convexity . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5

1.3 Summary of contents . . . . . . . . . . . . . . . . . . . . . . . . . . . 6

2 Foundations 8 2.1 Notation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8 2.2 Signals and systems . . . . . . . . . . . . . . . . . . . . . . . . . . . 9

2.2.1 Signal spaces and function spaces . . . . . . . . . . . . . . . . 10 2.2.2 Systems over finite time . . . . . . . . . . . . . . . . . . . . . 11

2.3 Linear fractional transformations . . . . . . . . . . . . . . . . . . . . 12 2.3.1 Structured uncertainty . . . . . . . . . . . . . . . . . . . . . . 14 2.3.2 Model-sets . . . . . . . . . . . . . . . . . . . . . . . . . . . . 15

2.4 Model validation problems . . . . . . . . . . . . . . . . . . . . . . . . 17

3 Literature Survey 19 3.1 Identification . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19 3.2 Frequency domain . . . . . . . . . . . . . . . . . . . . . . . . . . . . 20 3.3 Time domain . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22 3.4 Comparison . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 24

4 General Model Validation 25

4.1 LFT model-sets . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25 4.1.1 Unstructured uncertainty . . . . . . . . . . . . . . . . . . . . 26 4.1.2 Including noise . . . . . . . . . . . . . . . . . . . . . . . . . . 28 4.1.3 Structured uncertainty . . . . . . . . . . . . . . . . . . . . . . 30

4.2 Computational complexity . . . . . . . . . . . . . . . . . . . . . . . . 31 4.2.1 Central concepts . . . . . . . . . . . . . . . . . . . . . . . . . 31 4.2.2 Model validation is NP-hard . . . . . . . . . . . . . . . . . . 32 4.2.3 Comments . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34

4.3 Convexity . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34 4.3.1 Unstructured uncertainty . . . . . . . . . . . . . . . . . . . . 35 4.3.2 Structured uncertainty . . . . . . . . . . . . . . . . . . . . . . 38

CONTENTS

4.4 Computational Issues . . . . . . . . . . . . . . . . . . . . . . . . . . 39 4.4.1 Time-varying uncertainty . . . . . . . . . . . . . . . . . . . . 40 4.4.2 Time-invariant uncertainty . . . . . . . . . . . . . . . . . . . 40 5 Gap Metric and ν-Gap Metric 42

5.1 Normalized coprime factor models . . . . . . . . . . . . . . . . . . . 42 5.1.1 Left and right factorizations . . . . . . . . . . . . . . . . . . . 43 5.1.2 Including noise . . . . . . . . . . . . . . . . . . . . . . . . . . 45 5.1.3 Computational issues . . . . . . . . . . . . . . . . . . . . . . 47

5.2 The gap metric . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 48 5.2.1 Implications . . . . . . . . . . . . . . . . . . . . . . . . . . . . 49 5.2.2 Limitations . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51

5.3 The ν-gap metric . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 52 5.3.1 Implications . . . . . . . . . . . . . . . . . . . . . . . . . . . . 53 5.3.2 Including noise . . . . . . . . . . . . . . . . . . . . . . . . . . 56

6 Flexible Beam 58

6.1 Model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 59

6.1.1 Model-set . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 59 6.2 Experiment design . . . . . . . . . . . . . . . . . . . . . . . . . . . . 61

6.2.1 Data length . . . . . . . . . . . . . . . . . . . . . . . . . . . . 62 6.2.2 Computational issues . . . . . . . . . . . . . . . . . . . . . . 62

6.3 Experimental data . . . . . . . . . . . . . . . . . . . . . . . . . . . . 63 6.3.1 Measuring noise . . . . . . . . . . . . . . . . . . . . . . . . . 64 6.3.2 Invalidation . . . . . . . . . . . . . . . . . . . . . . . . . . . . 64 6.3.3 Modifying the model . . . . . . . . . . . . . . . . . . . . . . . 66 6.3.4 Validation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 66

6.4 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 68

7 Harrier VSTOL Aircraft 71 7.1 Model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 71 7.2 Experiment design . . . . . . . . . . . . . . . . . . . . . . . . . . . . 72

7.2.1 Computational issues . . . . . . . . . . . . . . . . . . . . . . 74 7.3 Simulated data . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 75

7.3.1 Measuring noise . . . . . . . . . . . . . . . . . . . . . . . . . 75 7.3.2 Validation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 77

7.4 Flight-test data . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 82 7.4.1 Measuring noise . . . . . . . . . . . . . . . . . . . . . . . . . 82 7.4.2 Validation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 82

8 Conclusions 86 8.1 Contributions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 86 8.2 Further Research . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 87

Bibliography 88

Chapter 1

Introduction

This thesis is concerned with the relationship between physical systems and mathematical models^{1}. The relationship is examined by testing whether or not a mathematical model could have produced data from the physical system. Testing mathematical models in this way is called model validation.

A physical system can never be modelled exactly by a mathematical model; there is always uncertainty. Uncertainty can arise from two sources; unknown inputs (eg noise) and unknown dynamics (eg high frequency dynamics). Both of these sources of uncertainty can be included in the mathematical model. For example, consider a discrete-time, time-domain model described by a linear operator P, with additive uncertainty and output noise, as shown in figure 1.1. A typical assumption about the uncertainty Δ is that Δ is a linear operator with induced 2-norm, which we denote by ∥ · ∥_{i2}, bounded by some positive number γ. A typical assumption about the noise signal n is that it is a discrete-time signal with ∥n∥_{∞} ≤ δ, for some positive real number δ. The output y, depends not only on the known input u, but also the

+ ^{+} ^{+} Δ  y u  n P  ^{+} v

Figure 1.1: Additive uncertainty model with sensor noise.

noise signal n, and the uncertainty signal, v. Hence, for a single input there are a set of possible outputs, so it is more accurate to talk about a set of models. We will use the term model-set to describe the set of possible models. In this example we denote the model-set by A_{1}(P, Δ, γ), where Δ is the set of stable, rational, transfer functions and

A_{1}(P, Δ, γ) := {[ I P + Δ ] : Δ ∈ Δ, ∥Δ∥_{i2} ≤ γ} .

Hence an element of A_{1}(P, Δ, γ) maps the inputs [ \begin{matrix}n \\ u \end{matrix} ], to the measured output y.

1We will use the terms “mathematical model” and “physical system” with the same meaning as discussed in [DFT92, ch.1]

Using model-sets, rather than a single model, enables a controller to be designed that is robust to modelling errors. A typical objective of robust controller design is to minimize the induced 2-norm of a closed loop transfer function. This is done to reduce the effect of unmeasured disturbances, or uncertainty, in the system dynamics. Denote the closed loop transfer function by T_{yu}, then a controller will make ∥T_{yu}∥_{i2} < γ for the physical system if:

1. ∥T_{yu}∥_{i2} < γ for every element of the model-set.

2. The physical system can be modelled exactly by an element of the model-set.

By saying a physical system can be modelled exactly by an element of the modelset, we mean that for every possible input, the corresponding output of the model is identical to the “actual” output of the system. If noise is present then “actual” inputs and outputs cannot be measured exactly, so we can only require the output of the model to be close to the measured output. For example, suppose a physical system is modelled by A_{1}(P, Δ, γ), and (u, y) is a set of observed data. Then we say the physical system can be modelled exactly by an element of the model-set if there exists P_{1} ∈ A_{1}(P, Δ, γ), and n with ∥n∥_{∞} < δ, such that

y = P_{1} [ \begin{matrix}n \\ u \end{matrix} ] .

Clearly it is not possible to check that a physical system can be modelled exactly by an element of the model-set, since it is not possible to measure every possible input and output of the physical system. However, it is possible to check whether some data from a physical system can be modelled exactly by an element of the model-set. If an element of the model-set models exactly measured data, then our confidence in the model-set is increased. If no element of the model-set can account for the data, then the model-set must be modified.

The design of controllers that are guaranteed to satisfy some performance objective for every element of a model-set, such as ∥T_{yu}∥_{i2} < γ, has received considerable attention from the control community, and controllers can now be designed for many different model-sets, for example see [DGKF89]. However, testing whether or not a physical system can be modelled exactly by an element of the model-set has received less attention from the control community. This is the subject of the thesis, and the main contribution is:

To develop techniques for testing the validity of model-sets used for robust control, and to apply these techniques to data obtained from a flexible beam, and Harrier VSTOL^{2} aircraft.

Controllers have been designed, and tested successfully, for both the flexible beam [WWa] and the Harrier [HG93, Hyd93]. It is therefore tempting to conclude that the model-sets used for design contained an element that modelled exactly the physical system^{3}. However, this is not necessarily true, as a simple example illustrates.

Suppose we must design a feedback controller to stabilize an unstable physical system, whose input-output behaviour is modelled by the transfer function P (s) =

2Vertical and Short Take-Off and Landing.

3If the controller did not satisfy the property that it was designed to achieve for the model-set, then the model-set is invalidated. An approach to model invalidation based on this observation is described in [LDF94].

1.1. INVALIDATION

\frac{1}{s−1}. Suppose the uncertainty in the model is thought to enter additively, so the model-set is A(P, Δ, γ), where

A(P, Δ, γ) := {P + Δ : Δ ∈ Δ, ∥Δ∥_{∞} ≤ γ} .

The set of controllers that robustly stabilizes^{4} A(P, Δ, γ) can be characterized using a well known result:

Theorem 1.1 ([DFT92]) A controller C(s) robustly stabilizes A(P, Δ, γ) if and only if

∥\frac{C(s)}{1 + P (s)C(s)}∥_{∞} < γ^{−1}.  (1.1)

Hence the simple controller C(s) = 2 will stabilize every element in A(Δ, \frac{1}{2}).

Now suppose that the controller C(s) is implemented on the physical system, and the resulting closed loop transfer function is stable. Can one conclude that the model-set contains an element that exactly models the physical system? The answer is no, which is demonstrated by taking, for a small positive number ϵ,

Δ(s) = \frac{1}{s − 1 − ϵ} − \frac{1}{s − 1}.

The closed loop transfer function from u to y is \frac{1}{s+1−ϵ}, which for ϵ < \frac{1}{2} is clearly stable. A simple calculation also shows that the closed loop transfer function \frac{1}{s+1−ϵ} implies that Δ(s) = \frac{1}{s−1−ϵ} − \frac{1}{s−1}. Hence, although the closed loop transfer function is stable, there is no element of the model-set that describes exactly the “true” system.

The example motivates a more general question: “What is the largest set of linear, time-invariant systems that can, a priori, be guaranteed to be stabilized by any controller which stabilizes a nominal model and satisfies a specified H_{∞}-norm bound on some closed loop transfer function?” It turns out that in some specific cases, any H_{∞} controller that is designed to tolerate a class of stable Δ will also stabilize some unstable Δ [Glo86, Vin92]. A characterization of the largest set of linear, time-invariant systems with this property is derived in [Vin92, ch.7]. These results are generalized in [Vina] to include certain sets of nonlinear and time-varying systems that are guaranteed to be stabilized. Validation of model-sets of this type will be a central theme in this thesis.

1.1 Invalidation

Model validation is a misleading phrase. A model can never be validated as data may be collected in the future which will invalidate it. This is exactly the situation in the physical sciences. A model of a system is conjectured, which remains consistent with experimental observations until new experiments reveal discrepancies between the predicted, and actual, observations. Newtonian dynamics is a classic example of this. This model accurately described the dynamics of particles, until particles travelling at very high speeds were observed. Einstein’s theory of relativity accurately explained these new observations, so Newtonian dynamics was replaced.

4Robustly stabilize in the sense that every element of A(P, Δ, γ) is internally stabilized by C.

See [DFT92] for a definition of internal stability.

1.2. CENTRAL THEMES

A more correct term than model validation would be “not invalidated”, as a set of data can never validate a model, only invalidate it. Therefore a model is either invalidated or not invalidated. However, the phrase model validation is widely used throughout science and engineering. So to conform with the commonly accepted term “model validation” will be used, but should be interpreted as meaning not invalidated.

1.2 Central themes

Two central themes in this thesis are interpolation and convexity. We will use results from interpolation theory to prove virtually all model validation results, and will use the property of convexity to distinguish between problems that are computationally easy and difficult.

1.2.1 Interpolation

To see why interpolation theory is important for model validation consider a simple model validation problem: Suppose the model-set used to describe a SISO^{5} physical system is G, where

G := {G(λ) : G(λ) = h_{0} + h_{1}λ + h_{2}λ^{2} + …, λ ∈ ℂ, sup_{|λ|<1} |G(λ)| < 1} .

Note that for λ = z^{−1}, every element of G is now a transfer function, with z the variable in the standard Z-transform. Hence the coefficients h_{i} are the Markov parameters of the system with transfer function G(z) [Kai80].

Suppose some data has been recorded from the system being modelled by G. Let the input to the system be the sequence u = (u_{0}, u_{1}, …, u_{l−1}) and the output y = (y_{0}, y_{1}, …, y_{l−1}). Then the model is validated for this data if there exists an element of the model-set G, that produces the output y when the input is u.

The output of an element in the model-set is obtained by a convolution of the Markov parameters with the inputs [Kai80]. This can be expressed, using a Toeplitz matrix formed from the Markov parameters, by the equation

[ y_{0} ; y_{1} ; ⋮ ; y_{l−1} ] = [ \begin{matrix}h_{0} & 0 & … & 0 \\ h_{1} & h_{0} & ... & ⋮ \\ ⋮ & & ... & 0 \\ h_{l−1} & h_{l−2} & … & h_{0} \end{matrix} ] [ u_{0} ; u_{1} ; ⋮ ; u_{l−1} ] .  (1.2)

This equation shows that, for u_{0} ≠ 0, the inputs and outputs uniquely determine the first l Markov parameters of the transfer function G(z). The model is validated if the remaining Markov parameters can be chosen so that G(z) ∈ G.

Carathéodory-Fejér interpolation is the branch of interpolation theory that is applicable to the model validation problem. The tangential Carathéodory-Fejér interpolation problem, as stated in [FF90], is:

5Single Input and Single Output

Problem 1.1 (Tangential Carathéodory-Fejér) Let u and y be two matrices of the form,

u^{T} = [u_{1}^{T} , u_{2}^{T} , …, u_{n}^{T} ] and y^{T} = [y_{1}^{T} , y_{2}^{T} , …, y_{n}^{T} ] .

Find necessary and sufficient conditions for the existence of an infinite contractive Toeplitz matrix Â_{∞} satisfying y = Â_{n}u, where Â_{n} is the n by n analytic Toeplitz matrix in the upper left hand corner of Â_{∞}.

From the statement of this problem it is clear that the solution to the tangential Carathéodory-Fejér problem also solves the model validation problem just described.

The tangential Carathéodory-Fejér problem has been solved, for example see [FF90, Thm.1.6] or [Dym89, Thm.6.4], for MIMO^{6} systems, where the h_{i} are matrices. The elegant solution is:

Theorem 1.2 ([FF90]) There exists a solution to the tangential Carathéodory-Fejér interpolation problem if and only if Y^{T} Y ≤ U^{T} U, where

U := [ \begin{matrix}u_{0} & 0 & … & 0 \\ u_{1} & u_{2} & ... & ⋮ \\ ⋮ & & ... & 0 \\ u_{l−1} & u_{l−2} & … & u_{0} \end{matrix} ] and Y := [ \begin{matrix}y_{0} & 0 & … & 0 \\ y_{1} & y_{2} & ... & ⋮ \\ ⋮ & & ... & 0 \\ y_{l−1} & y_{l−2} & … & y_{0} \end{matrix} ] .

This result is extremely important and we use it throughout the thesis. Toeplitz matrices with the structure in equation 1.2 correspond to time-invariant linear systems, and we will also use a similar interpolation result for time-varying systems.

1.2.2 Convexity

In the example shown in figure 1.1, the output y depends on the known input u and the signals v and n. Given a set of measured data (u, y), and norm bounds γ

and δ, the model-set is validated if there exists a signal [ \begin{matrix}v \\ n \end{matrix} ] satisfying the equality

constraint

y = P u + v + n,  (1.3)

and certain inequality constraints, eg ∥n∥_{∞} < δ, V^{T} V ≤ γ^{2}U^{T} U, where V and U are block lower Toeplitz matrices formed from v and u respectively. In general there

will be a set of possible [ \begin{matrix}v \\ n \end{matrix} ] satisfying both equation 1.3 and the norm bounds

(which may be empty), call this set V. The model will be validated if there exists a

[ \begin{matrix}v \\ n \end{matrix} ] ∈ V. If V is a convex set then validating the model-set is equivalent to solving

a convex feasibility problem.

In the optimization community, convex feasibility problems are considered to be computationally tractable, in both theory and practice. The reasons are succinctly described in [Roc93, p.194]:

6Multiple Input and Multiple Output

1.3. SUMMARY OF CONTENTS

One distinguishing idea which dominates many issues in optimization theory is convexity … An important reason is the fact that when a convex function is minimized over a convex set every locally optimal solution is global. Also, first-order necessary conditions turn out to be sufficient. A variety of other properties conducive to computation and interpretation of solutions ride on convexity as well. In fact the great watershed in optimization isn’t between linearity and nonlinearity, but convexity and nonconvexity.

It is also stated in [BBFE93]:

In the past , a “solution to a problem” generally meant a “closed-form” or “analytic” solution. We believe that in the future, our concept of “solution” should be extended to include many forms of convex programming.

Experience has shown that although a problem can be solved by computing the solution to a convex programming problem, the reality of present day computer hardware and software deems the problem intractable. Consequently, the future is not yet here and convex problems can not be considered as “solved”.

1.3 Summary of contents

Chapter 2: Foundations Here we lay the mathematical foundations upon which the thesis is built. Once the notation has been stated we define a class of signals and systems. All signals and systems are in discrete time since the interpolation results required to prove the results are for discrete time. The model-sets used in robust control are defined and the model validation problem stated.

Chapter 3: Literature Survey Little has been published on model validation for the model-sets used in robust control, but there are two main approaches which may be distinguished by their setting; one is in the time domain and the other in the frequency domain. The frequency domain approach contains results for a general class of model validation problems and the results have implications for the computational complexity of model validation. The time domain approach is the one adopted in this thesis and is based on results in Carathéodory-Fejér interpolation theory.

Chapter 4: General Model Validation We extend the results from [PKT^{+}94] to a class of general model-sets used in robust control. The remainder of the chapter is devoted to testing the computational complexity of the various generalizations. We employ results from computational complexity theory to show that it is unlikely that a polynomial time algorithm will be found to solve a general model validation problem. Hence, important special cases of model validation problems are considered and their convexity analysed. Finally, we discuss issues relating to the computation of solutions.

Chapter 5: Gap Metric and ν-Gap Metric The results from chapter 4 have an interesting interpretation in the gap metric [GS90]. This arises when considering the class of model-sets obtained by perturbing the normalized coprime factors of the nominal model. This class of model is important in robust control and the chapter begins with results for this class of model-set. The gap metric is then described and results from [GS90] used to interpret model validation in the gap metric. Model validation is then considered in the ν-gap metric [Vin93], which allows validation of the largest set of models than can be guaranteed to be stabilized by a certain set of H_{∞} controllers.

Chapter 6: Flexible Beam A flexible beam is taken as a case study. The model is computationally simple, having one input and one output, but the beam provides a challenging robust control problem. Various practical issues are discussed, such as the effect on the results of the length of data and the effect of noise. The model is invalidated for a set of data before a modified model is validated.

Chapter 7: Harrier VSTOL Aircraft In this chapter we consider the validity of models of a Harrier. In contrast to the flexible beam, the models of the Harrier are more complicated, having three inputs and three outputs. Details of the validation experiments are given and the linear model is validated using data obtained from a nonlinear model. Finally, flight-test data is used to validate the models used for robust controller design.

Chapter 8: Conclusions We summarize the contributions and suggest possibilities for future work.

Foundations

In this chapter we lay the mathematical foundations of the thesis and state the model validation problems examined. This requires definitions of the fundamental concepts of signals, systems and model-sets. All signals and systems will be in the time domain and over discrete time, since in the time domain causality is more easily imposed, and discrete time is natural for sampled measurements of a physical system.

Toeplitz operators map the inputs to the outputs in the time domain [You88, ZDG]. It is therefore not surprising that Toeplitz operators play a fundamental role in this thesis. When these operate over finite time, which is the case in validation problems as data is finite, Toeplitz operators simplify to Toeplitz matrices. Hence properties of Toeplitz matrices, and more specifically block lower triangular Toeplitz matrices, are important.

The model-sets considered are those used for H_{∞} controller design. A general class of model-sets are those that can be expressed as a linear fractional transformation on a set of transfer functions with bounded norm [DGKF89]. An important subset of these are the models that can be expressed as norm bounded perturbations to the normalized coprime factors of the nominal model [MG90]. We will call these norm bounded perturbations the uncertainty in the model-set, which is consistent with its robust control interpretation.

Different types of uncertainty may be considered, and allowing extra structure leads naturally to the structured singular value, μ, introduced in [Doy82]. This function, together with its associated theory, has important implications for model validation. Firstly, results from μ-theory are used in the frequency domain approach to model validation described in the next chapter. This was the first approach to model validation from a robust control perspective. Secondly, results on the computational complexity of μ have consequences for the computational complexity of model validation, which we will describe in chapter 4.

2.1 Notation

We will use ℝ, ℂ, ℚ, ℕ and ℤ to denote the spaces of real numbers, complex numbers, rational numbers, natural numbers and integers respectively. We will use the subscript + or − to denote the restriction to non-negative, or negative numbers respectively, eg ℝ_{+} = {x : x ∈ ℝ, x ≥ 0}, ℝ_{−} = {x : x ∈ ℝ, x < 0}. ϕ will denote the empty set. Calligraphic letters, for example S and T I, denote sets and spaces.

2.2. SIGNALS AND SYSTEMS

Roman upper case letters, for example A, B, X, denote operators or matrices, and Roman lower case letters, such as u, y, v, w, denote vectors or sequences. Greek letters, such as γ, δ are reserved for scalars.

The symbol := should be read as “is defined to be”, ⊕ denotes the direct sum, A^{T} denotes the transpose of a vector or matrix, A, and A^{∗} the conjugate transpose. Also σ̄(A) will denote the largest singular value of the matrix A. For a system with transfer matrix G(z), G(z)^{∗} := G(1/z)^{T}. P^{p×m} will denote the set of all rational transfer function matrices with m inputs and p outputs.

Linear Matrix Inequality is abbreviated to LMI and is a matrix inequality of the form

F (x) := F_{0} + ∑_{i=1}^{m} x_{i}F_{i} > 0,

where x ∈ ℝ^{m} is a variable and F_{i} = F_{i}^{T} ∈ ℝ^{n×n}, i = 0, 1, …, m are given. The set of matrices {F_{i} : i = 0, 1, … m} is called the LMI basis, as the LMI is a linear combination of elements of the basis. A key property of LMI’s, that is used to prove convexity of model validation problems, is that the set {x : F (x) > 0} is convex [BEFB94]. Note that LMI’s should be more correctly called Affine Matrix Inequalities as expressions of the form F_{0} +∑_{i=1}^{m} x_{i}F_{i} are affine in x_{i} and not linear. However, LMI has become the commonly accepted term.

2.2 Signals and systems

As signals and systems will be in discrete time, it is natural to consider sequence spaces. Let S^{m} denote the set of sequences with elements is ℝ^{m}, that is

S^{m} := {(…, u_{−1}, u_{0}, u_{1}, …) : u_{i} ∈ ℝ^{m}, i ∈ ℤ} .

The one sided sequence spaces, denoted by Sm_{+} and Sm_{−}, consist of sequences of the form (u_{0}, u_{1}, …) and (…, u_{−2}, u_{−1}), respectively. We will deal almost exclusively with signals in Sm_{+}.

Two elementary operations on elements of Sm_{+} are given in the following definition:

Definition 2.1 The k-step truncation operator, denoted by π_{k}, is the map

π_{k} : Sm_{+} → Sm_{+},
(u_{0}, u_{1}, …) ↦ (u_{0}, u_{1}, …, u_{k−1}, 0, 0, …).

The right shift operator, denoted by S, is the map

S : Sm_{+} → Sm_{+},
(u_{0}, u_{1}, …) ↦ (0, u_{0}, u_{1}, …).

The notation π_{l}Sm_{+} will be used as a concise way of denoting the set of sequences with l elements. This will often be used for denoting the set that contains input and output data, eg x ∈ π_{l}Sm_{+} means x = (x_{0}, x_{1}, …, x_{l−1}), where x_{i} ∈ ℝ^{m} for i = 0, 1, …, l − 1.

A linear system, H, will be regarded as a linear operator from Sm_{+} to Sp_{+}, which is causal if, given any v and w in Sm_{+}, π_{k}v = π_{k}w implies π_{k}Hv = π_{k}Hw for all k ∈ ℕ. A linear system, H, is time-invariant if it commutes with the right shift operator, that is HS = SH.

2.2.1 Signal spaces and function spaces

One of the key properties of a signal is its’ size. This can be quantified using norms, which naturally arise when considering inner-products. The notion of norm on a signal can also be used to define the notion of norm on a system, producing the concept of induced norm.

An inner-product is defined as the map

< ·, · >: S^{m} × S^{m} → ℝ,
(u, v) ↦ ∑_{i=−∞}^{∞} u_{i}^{T} v_{i},

when the series converges. This can be used to define the 2-norm, ∥ · ∥_{2}, of a sequence u ∈ S^{m}; ∥u∥_{2} := \sqrt{< u, u >}. The infinity-norm, ∥ · ∥_{∞}, can also be defined for sequences in S^{m}, but this is not induced by an inner-product [You88, p.23]. Define the infinity-norm of an element in ℝ^{m} as the maximum absolute value of its elements. Then for any u ∈ S^{m}, ∥u∥_{∞} := sup_{i∈ℤ} ∥u_{i}∥_{∞}.

The Hilbert space l_{2}^{m} is defined as

l_{2}^{m} := {u ∈ S^{m} : ∥u∥_{2} < ∞} .

The Hilbert spaces l_{2+}^{m} and l_{2−}^{m} are defined similarly for signals in Sm_{+} and Sm_{−} . The induced 2-norm of a system, H, is given by:

∥H∥_{i2} := sup_{u∈l_{2+}^{m},u≠0} \frac{∥Hu∥_{2}}{∥u∥_{2}} ,

if the supremum is finite. If ∥H∥_{i2} is finite, we will say H is stable. A stable system maps l_{2+}^{m} into l_{2+}^{m}, which motivates the study of the set of systems that map l_{2+}^{m} into l_{2+}^{m}. This set of systems turns out to be a Hardy space [You88] which we denote by H_{∞}.

The Hardy space H_{∞} is usually defined in the frequency domain, and the time and frequency domains are related using the discrete Fourier transform [You88, ch.13]. The Fourier transform provides an isometric isomorphism, denoted by ≅, between time and frequency domain spaces,

l_{2} ≅ L_{2},
l_{2+} ≅ H_{2},
l_{2−} ≅ H_{2}^{⊥}.

The Banach space L_{∞} denotes the space of essentially bounded matrix functions on the unit circle, with the norm of an element H given by

∥H∥_{∞} := ess sup_{θ∈[0,2π]} σ̄(H(e^{iθ})).

The Hardy space H_{∞} is the subspace of L_{∞} that contains those elements that are analytic inside the unit disc.

The prefix R on any of these spaces will denote the subspace consisting of rational functions. For example if H(λ) ∈ RH_{∞}, H(λ) will consist of those elements of H_{∞} that can be expressed as rational functions of λ. Note that we will use the argument λ when H_{∞} corresponds to functions analytic inside the unit disc. This is in contrast to the standard Z-transform, here λ = z^{−1}. We will use the argument z to describe functions in numerical examples, eg P (z) = \frac{7z−9}{3z−1} is a stable transfer function analytic outside the unit disc.

2.2.2 Systems over finite time

One of the most frequently used objects in this thesis is the Toeplitz matrix. Each Toeplitz matrix is associated with a transfer function or vector:

Definition 2.2 Let P ∈ P^{p×m} and H = (H_{0}, H_{1}, …, H_{l−1}), where H_{i} ∈ ℝ^{p×m} (i = 0, 1, …, l − 1), be the finite impulse response of P. Then the block lower Toeplitz matrix of P, T_{P}, is defined as

T_{P} := [ \begin{matrix}H_{0} & 0 & … & 0 \\ H_{1} & H_{0} & ... & ⋮ \\ ⋮ & ⋮ & ... & 0 \\ H_{l−1} & H_{l−2} & … & H_{0} \end{matrix} ] .

Definition 2.3 Let v = (v_{0}, v_{1}, …, v_{l−1}) be a sequence of l vectors, where v_{i} ∈ ℝ^{m} for i = 0, 1, …, l − 1. Then the block lower Toeplitz matrix of v, T_{v}, is defined as

T_{v} := [ \begin{matrix}v_{0} & 0 & … & 0 \\ v_{1} & v_{0} & ... & ⋮ \\ ⋮ & ⋮ & ... & 0 \\ v_{l−1} & v_{l−2} & … & v_{0} \end{matrix} ] .

We will use vec(v) to denote the vector in ℝ^{lm} formed by stacking the elements of v,

vec(v) := [ v_{0} ; v_{1} ; ⋮ ; v_{l−1} ] .

Using these definitions we can neatly express the relationship between the input of a linear system and its output, in the time domain. Suppose the input to a system P ∈ P^{p×m} is the sequence u = (u_{0}, u_{1}, …, u_{l−1}), then the output sequence y = (y_{0}, y_{1}, …, y_{l−1}) is obtained from the equation

vec(y) = T_{P} vec(u),

assuming initial conditions are zero. For convenience vec will be omitted in most equations. This should not cause confusion as whenever an equation contains a sequence u it should be interpreted as being the vector vec(u).

To illustrate the above notation consider the following example:

Example 2.1 Consider a linear, time-invariant system described by the set of discrete time, state space equations,

x_{k+1} = Ax_{k} + Bu_{k},
y_{k} = Cx_{k}.

The finite impulse response of this system is (0, CB, CAB, CA^{2}B, …, CA^{l−2}B). Hence the output of the system, with x_{0} = 0, is given by the equation

[ y_{0} ; y_{1} ; ⋮ ; y_{l−1} ] = [ \begin{matrix}0 & 0 & … & … & 0 \\ CB & 0 & & & 0 \\ CAB & CB & ... & & ⋮ \\ ⋮ & & ... & ... & ⋮ \\ CA^{l−2}B & CA^{l−3}B & … & CB & 0 \end{matrix} ] [ u_{0} ; u_{1} ; ⋮ ; u_{l−1} ] .

2.3. LINEAR FRACTIONAL TRANSFORMATIONS

For ease of notation we will allow the k-step truncation operator, π_{k}, to operate on Toeplitz matrices, the result being

π_{k}T_{P} := [ \begin{matrix}H_{0} & 0 & … & 0 \\ H_{1} & H_{0} & ... & ⋮ \\ ⋮ & ⋮ & ... & 0 \\ H_{k−1} & H_{k−2} & … & H_{0} \end{matrix} ] .

This allows π_{k}T_{P} u to be expressed as π_{k}T_{P} π_{k}u. Note that this is only possible because T_{P} is block lower triangular. For a full matrix V, π_{k}V u can only be expressed as V_{k}u where V_{k} consists of the first k rows of V.

As the Toeplitz matrices in definition 2.2 correspond to linear time-invariant systems, the space of all matrices of this form will be denoted by T I_{l}^{p×m}. The space of all Toeplitz matrices corresponding to time-varying systems will be denoted by T V_{l}^{p×m}, which consists of matrices of the form

[ \begin{matrix}A_{1,1} & 0 & … & 0 \\ A_{2,1} & A_{2,2} & ... & ⋮ \\ ⋮ & ⋮ & ... & 0 \\ A_{l,l} & A_{l,2} & … & A_{l,l} \end{matrix} ] .

In both T I_{l}^{p×m} and T V_{l}^{p×m}, l will be allowed to be infinite where T I_{∞}^{p×m} and T V_{∞}^{p×m} denote the space of Toeplitz operators with the appropriate structure.

2.3 Linear fractional transformations

A general model-set used in robust control consists of models that can be represented as a linear fractional transformation (LFT) on an unspecified, but norm bounded, transfer function. Additive, multiplicative and coprime factor uncertainty descriptions can all be expressed as an LFT on the uncertainty, with a suitable choice of the coefficient matrix.

P_{12} P_{22} w  v y  u P_{21} P_{11} Δ

Figure 2.1: Upper LFT of P on Δ

An upper LFT of the coefficient matrix P on the uncertainty matrix Δ is represented symbolically in figure 2.1. The diagram represents the equations

w = P_{11}v + P_{12}u,  (2.1)
y = P_{21}v + P_{22}u,  (2.2)
v = Δw,  (2.3)

where the P_{ij} (i, j = 1, 2) are compatibly dimensioned, possibly complex, matrices and u, y, w and v vectors. Let Δ also be a matrix, then if det(I − P_{11}Δ) ≠ 0 the map from u to y, denoted by F_{u}(P, Δ), can be calculated to be

F_{u}(P, Δ) = P_{22} + P_{21}Δ(I − P_{11}Δ)^{−1}P_{12}.

This definition of an LFT is useful as it can be used to define LFT’s in the frequency domain or the time domain. In the frequency domain P and Δ are transfer functions and the LFT of P (λ) and Δ(λ) is defined at each frequency using the matrices P (e^{jθ}) and Δ(e^{jθ}). In the time domain P and Δ are Toeplitz matrices so equations (2.1),(2.2) and (2.3) should be interpreted as

w = T_{P_{11}}v + T_{P_{12}}u,
y = T_{P_{21}}v + T_{P_{22}}u,  (2.4)
v = T_{Δ}w.

Which meaning is relevant will be made clear from the context in which it is used.

Lower LFT’s can be defined similarly with the Δ matrix feeding back “underneath” P. If det(I − P_{22}Δ) ≠ 0 then the map from u to y will be denoted by F_{l}(P, Δ) and can be calculated to be

F_{l}(P, Δ) = P_{11} + P_{12}Δ(I − P_{22}Δ)^{−1}P_{21}.

The LFT of figure 2.1 can be modified to account for noise, or other unknown inputs. This is accomplished by considering an extra input, n, which is unknown but bounded in some norm. The block diagram of this model is shown in figure 2.2, which represents the equations

w = P_{11}v + P_{12}n + P_{13}u,
y = P_{21}v + P_{22}n + P_{23}u,  (2.5)
v = Δw.

Again these equations should be interpreted with their dual meaning.

w  Δ  v P_{23} P_{13} P_{22} P_{12} P_{21} P_{11} y  u n

Figure 2.2: Upper LFT with noise

The transfer function, or Toeplitz matrix, Δ is used to allow for the uncertainty in the model. If the only restriction placed on Δ is that it lies in some space, for example RH_{∞} or T V_{∞}^{p×m}, and has norm bounded by some number then it is said to be unstructured. However, in some situations it is natural to impose additional structure on Δ, and it is then said to be structured.

2.3.1 Structured uncertainty

The extra structure in the uncertainty can arise in a variety of ways and corresponds to uncertain parameters and unmodelled dynamics [Pac88]. In the frequency domain this gives rise to the repeated real, repeated complex and full complex blocks of the mixed μ setup in [You93].

In [You93] the block structure is described by the triple of integers (m_{r}, m_{c}, m_{C}), where m_{r} is the number of repeated real blocks, m_{c} the number of repeated complex blocks and m_{C} the number of full complex blocks. The dimensions of these blocks are given by the block structure K(m_{r}, m_{c}, m_{C}) where

K := (k_{1}, …, k_{m_{r}}, k_{m_{r}+1}, …, k_{m_{r}+m_{c}}, k_{m_{r}+m_{c}+1}, …, k_{m_{r}+m_{c}+m_{C}}).

Then the set of allowable perturbations, X_{K}, is given by

X_{K} := {Δ = blockdiag(δ_{1}^{r}I_{k_{1}}, … δ_{m_{r}}^{r}I_{k_{mr}} , δ_{1}^{c}I_{k_{mr+1}}, …, δ_{m_{c}}^{c}I_{k_{mr+mc}} ,

Δ_{1}^{C}, …, Δ_{m_{C}}^{C}) : δ_{i}^{r} ∈ ℝ, δ_{i}^{c} ∈ ℂ, Δ_{i}^{C} ∈ ℂ^{k_{mr+mc+i}×k_{mr+mc+i}}} .

Note that X_{K} assumes that the uncertainty blocks are square. Given a system with non-square blocks, zero rows or columns can be appended to the blocks, and zero inputs and outputs added, to make every block square. This simplifies the notation and we will always assume structured blocks are square.

Allowing structure in the time domain is similar but more difficult to write down precisely. When the uncertainty is time-invariant the Δ will be a block Toeplitz matrix in T I_{l}^{p×m}, but with each block possessing additional structure. This structure will change between the blocks on the diagonal of Δ and those off it; the memoryless elements corresponding to the repeated real blocks will only appear in the diagonal blocks.

Stating this formally let (m_{r}, m_{c}, m_{C}) denote the number of blocks with the different structures, as before. Also let K(m_{r}, m_{c}, m_{C}) denote the block structure and S_{K} the set of matrices corresponding to the repeated real blocks, where

S_{K} := {Δ = blockdiag(δ_{1}^{r}I_{k_{1}}, …, δ_{m_{r}}^{r}I_{k_{mr}} , 0_{k_{mr+1}}, …, 0_{k_{mr} +k_{mc}+k_{m_{C}}}) :
δ_{i}^{r} ∈ ℝ} ,

and F_{K} the set of matrices corresponding to the repeated complex and full blocks where

F_{K} := {Δ = blockdiag(0_{k_{1}}, …, 0_{k_{mr}} , δ_{1}^{c}I_{k_{mr+1}}, …, δ_{m_{c}}^{c}I_{k_{mr+mc}} , Δ_{1}, …, Δ_{m_{f}}) :

δ_{i}^{c} ∈ ℝ, Δ_{i} ∈ ℝ^{k_{mr+mc+i}×k_{mr+mc+i}}} .

Then we denote the set of allowable structured time-invariant perturbations (over l time steps) by ST I_{K}(l) where

ST I_{K}(l) := \{[ \begin{matrix}A_{1} & 0 & … & 0 \\ A_{2} & A_{1} & ... & ⋮ \\ ⋮ & ⋮ & ... & 0 \\ A_{l} & A_{l−1} & … & A_{1} \end{matrix} ] : A_{1} ∈ \{S_{K} + F_{K}\} , A_{i} ∈ F_{K}, i = 2, ..., l\} .

Note that as for T I_{l}^{p×m}, we will allow the l in ST I_{J} (l) to be infinite, to denote the space of operators of the same form. Also, the above notation can easily be generalized to allow the perturbations to be time-varying (or a combination of timevarying and time-invariant), but that will not be necessary for the thesis.

Example 2.2 Let m_{r} = 1, m_{c} = 1, m_{C} = 1 and K = (2, 1, 2) then

S_{K} = \{Δ = [ \begin{matrix}δ_{1} & 0 & 0 & 0 & 0 \\ 0 & δ_{1} & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 \end{matrix} ] : δ_{1} ∈ ℝ\} ,

F_{K} = \{Δ = [ \begin{matrix}0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & δ_{2} & 0 & 0 \\ 0 & 0 & 0 & Δ_{11} & Δ_{12} \\ 0 & 0 & 0 & Δ_{21} & Δ_{22} \end{matrix} ] : δ_{2} ∈ ℝ, Δ_{ij} ∈ ℝ\} .

The class of allowable structured perturbations, ST I_{K}(l), is the set of block lower Toeplitz matrices formed from both S_{K} and F_{K} on the diagonal and F_{K} below the diagonal.

A useful projection that is required to describe the results in [Smi90] is denoted by R_{i}. We will use it to split up the signals v and w compatibly with the block structure X_{K}, that is R_{i}w are the inputs to the i-th block and R_{i}v are the outputs from the i-th block. It is defined as

R_{i} := [ 0_{k_{1}} … 0_{k_{i−1}} I_{k_{i}} 0_{k_{i+1}} … 0_{k_{mr+mc+m_{C}}} ] .

With time domain LFT models R_{i} is more difficult to write down, but will have the same action; of partitioning v and w conformally with the block structure ST I_{J} (l).

2.3.2 Model-sets

Two model-sets often used for robust control design are a nominal model with additive uncertainty, A(P, Δ, γ), and a nominal model with multiplicative uncertainty, M(P, Δ, γ). Given P ∈ P^{p×m}

A(P, Δ, γ) := {P + Δ : Δ ∈ Δ, ∥Δ∥_{i2} < γ} ,
M(P, Δ, γ) := {P (I + Δ) : Δ ∈ Δ, ∥Δ∥_{i2} < γ} .

More general model-sets can be obtained using LFT’s. Given P ∈ P^{(n_{w}+n_{y})×(n_{v}+n_{u})} the LFT on the uncertainty Δ is defined by equations (2.4). If I −P_{11}Δ is invertible for all Δ ∈ Δ then the model-set, which we will denote by LFT (P, Δ, γ), is defined as

LFT (P, Δ, γ) := {F_{u}(P, Δ) : Δ ∈ Δ, ∥Δ∥_{i2} < γ} ,

with P partitioned appropriately. Similarly, to represent models with noise whose norm is bounded δ, given P ∈ P^{(n_{w}+n_{y})×(n_{v}+n_{n}+n_{u})} the LFT on the uncertainty Δ is defined by equations (2.5). Again, if I − P_{11}Δ is invertible for all Δ ∈ Δ then the model-set can be defined, which we now denote by LFT (P, Δ, γ, δ). Note that in these definitions all the P_{ij} are frequency domain transfer functions. In the time domain P, for an impulse response of length l, will be in T I_{l}^{(n_{w}+n_{y})×(n_{v}+n_{u})}, T V_{l}^{(n_{w}+n_{y})×(n_{v}+n_{u})}, T I_{l}^{(n_{w}+n_{y})×(n_{v}+n_{n}+n_{u})} or T V_{l}^{(n_{w}+n_{y})×(n_{v}+n_{n}+n_{u})}. These dimensions will be assumed for the remainder of the thesis.

An important subset of LFT (P, Δ, γ) is the set of all models than can be expressed as perturbations in the normalized coprime factors of the model, as in figure 2.3. Before defining this model-set, coprime factorizations must be defined.

u  y u  u w_{1}  v  w_{2}  v_{1}  w  v_{2}^{+} ^{+}^{+}  ^{+}^{+}^{+}^{−} (a)  (b)^{−}  N ΔÑ  ΔM̃ M̃^{−1}  Δ_{M}  Δ_{N} M^{−1} Ñ

Figure 2.3: Left (a) and right (b) normalized coprime factor models

Definition 2.4 Given M̃, Ñ ∈ H_{∞} with the same number of rows then M̃ and Ñ are left coprime if there exists X̃ ∈ H_{∞} and Ỹ ∈ H_{∞} satisfying

ÑỸ + M̃X̃ = I.

Definition 2.5 (Ñ, M̃) is a left coprime factorization of P ∈ P^{p×m} if

i) M̃ is invertible in P^{p×m}, ii) P = M̃^{−1}Ñ, iii) Ñ and M̃ are left coprime.

Definition 2.6 (Ñ, M̃) is a normalized left coprime factorization of P ∈ P^{p×m} if (Ñ, M̃) is a left coprime factorization of P and

ÑÑ^{∗} + M̃M̃^{∗} = I.

Right NCF models are defined similarly and the precise definitions can be found in [Vid87]. Results from [Vid87] also ensure that any P ∈ P^{p×m} has a normalized right coprime factorization and a normalized left coprime factorization. Using this result the set of left and right normalized coprime factor models can be defined. For any P ∈ P^{p×m} let (N, M) be a right normalized coprime factorization of P and (Ñ, M̃) be a left normalized coprime factorization.

We will use N CF(Ñ, M̃, Δ, γ) and N CF(N, M, Δ, γ) to denote left and right NCF model-sets respectively. Given P ∈ P^{p×m} with normalized left coprime factorization (Ñ, M̃), normalized right coprime factorization (N, M) and γ ≤ 1.

2.4. MODEL VALIDATION PROBLEMS

N CF(Ñ, M̃, Δ, γ) := {(M̃ + ΔM̃)^{−1}(Ñ + ΔÑ) : [ ΔÑ ΔM̃ ] ∈ Δ,

∥ [ ΔÑ ΔM̃ ] ∥_{i2} < γ} ,

N CF(N, M, Δ, γ) := \{(N + Δ_{N})(M + Δ_{M})^{−1} : [ \begin{matrix}Δ_{N} \\ Δ_{M} \end{matrix} ] ∈ Δ,

∥[ \begin{matrix}Δ_{N} \\ Δ_{M} \end{matrix} ]∥_{i2} < γ\} .

We will use ϵ_{max}, as in [MG90], to denote the largest value of γ, such that there exists a single controller that stabilizes^{1} every element of N CF(Ñ, M̃, Δ, γ). It turns out that ϵ_{max} has the same value for N CF(N, M, Δ, γ) [MG90].

The same as with the LFT model-sets, we distinguish sets with noise by the additional parameter, δ, which is the bound on the norm of the noise. The noise can enter at a variety of places in the model and the consequences for validation will be described in chapter 5. Considerations of convexity result in the NCF model-sets being restricted to a specific form. Given P ∈ P^{p×m} with normalized left coprime factorization (Ñ, M̃)

N CF(Ñ, M̃, Δ, γ, δ) := {[ M̃^{−1} (M̃ + ΔM̃)^{−1}(Ñ + ΔÑ) ] :
[ ΔÑ ΔM̃ ] ∈ Δ, ∥[ ΔÑ ΔM̃ ]∥_{i2} < γ} .

Any element of N CF(Ñ, M̃, Δ, γ, δ) will map the noise and input signal, [ \begin{matrix}n \\ u \end{matrix} ], to

the output, y.

The model-sets A(P, Δ, γ), M(P, Δ, γ), and N CF(Ñ, M̃, Δ, γ), where Δ = T I_{∞}^{p×m} or T V_{∞}^{p×m}, have been used for validation in [PKT^{+}94, ZK92] and one of the contributions of this thesis is the generalization of the results in [PKT^{+}94, ZK92] to the model-set LFT (P, Δ, γ).

2.4 Model validation problems

Now the notation has been defined we can state the validation problems studied in this thesis. The generic question is:

Given a model-set and a set of input-output data, does there exist a model in the model-set that could have produced the data?

We will consider different model-sets, a general one being formed from an LFT of P on Δ. The model validation problem for models of this form, which we will call the LFT model validation problem, is

Problem 2.1 Given a model P ∈ P^{(n_{w}+n_{y})×(n_{v}+n_{u})}, an uncertainty set Δ, a bound on the induced 2-norm of the uncertainty, γ, and a set of input-output data u ∈ π_{l}Sn_{u}_{+} and y ∈ π_{l}Sn_{y}_{+}, does there exist Δ ∈ Δ such that the following equations have a solution for w

(I − T_{P_{11}}π_{l}Δ)w = T_{P_{12}}u,
y = T_{P_{21}}π_{l}Δw + T_{P_{22}}u,
∥Δ∥_{i2} ≤ γ?

1For precise definitions of stability see [MG90] or chapter 5.

This problem does not consider the possible effects of noise. If this is included the question becomes

Problem 2.2 Given a model P ∈ P^{(n_{w}+n_{y})×(n_{v}+n_{n}+n_{u})}, an uncertainty set Δ, a bound on the induced 2-norm of the uncertainty, γ, a bound on the 2-norm of the noise, δ, and a set of input-output data u ∈ π_{l}Sn_{u}_{+} and y ∈ π_{l}Sn_{y}_{+}, does there exist Δ ∈ Δ such that the following equations have a solution for w and n

(I − T_{P_{11}}π_{l}Δ)w = T_{P_{12}}n + T_{P_{13}}u,
y = T_{P_{21}}π_{l}Δw + T_{P_{22}}n + T_{P_{23}}u,
∥Δ∥_{i2} ≤ γ,
∥n∥_{2} ≤ δ?

We will use two different norms to bound the noise, the 2-norm and the infinitynorm. All results will be stated with a 2-norm bound on the noise but, with the exception of when computational issues are considered, the 2-norm can be replaced by the infinity-norm without changing convexity properties of the problem. Which norm is chosen is usually a matter of preference, but choosing the infinity-norm has some computational advantages, described in section 4.4.

All of these problems are decision problems, as they have a yes or no answer. We will called them model validation decision problems (MVDP’s). An alternative approach is to ask the question;

Given a model-set (parameterized by γ) and a set of input-output data, what is the smallest γ such that there exists a model in the model-set that could have produced the data?

We will call this type of problem a model validation optimization problem (MVOP), since γ is optimized. For LFT model-sets including noise the MVOP is,

Problem 2.3 Given a model P ∈ P^{(n_{w}+n_{y})×(n_{v}+n_{n}+n_{u})}, an uncertainty set Δ, a bound on the 2-norm of the noise, δ, and a set of input-output data u ∈ π_{l}Sn_{u}_{+} and y ∈ π_{l}Sn_{y}_{+}, what is the smallest value of γ such that there exists Δ ∈ Δ such that the following equations have a solution for w and n

(I − T_{P_{11}}π_{l}Δ)w = T_{P_{12}}n + T_{P_{13}}u,
y = T_{P_{21}}π_{l}Δw + T_{P_{22}}n + T_{P_{23}}u,
∥Δ∥_{i2} ≤ γ,
∥n∥_{2} ≤ δ?

Clearly answering an MVOP provides an answer to the corresponding MVDP. Also by answering an MVDP, the corresponding MVOP can be solved by doing a bisection search over γ, and solving an MVDP at each iteration. Hence, from a computational viewpoint, MVDP’s are equivalent to MVOP’s, since being able to solve one enables the other to be solved. Consequently, we will use model validation problem to refer to both MVDP’s and MVOP’s, the distinction only being made when necessary.

Literature Survey

In this chapter we review the approaches to model validation in the literature and describe the results that are relevant our problem. The main approach to model validation in the literature is from an identification viewpoint and model validation is a step in the identification process. However, identification methods are not as suitable for validating the model-sets we consider as recent work in the frequency domain [Smi90, SD90, SD92] and time domain [PKT^{+}92, ZK92, ZK93, PKT^{+}94]. These two approaches provide results for validating special cases of A(P, Δ, γ), M(P, Δ, γ), N CF(Ñ, M̃, Δ, γ) and LFT (P, ST I_{K}(∞), γ).

Many papers have been published on model validation, in fields ranging from engineering [KGED94, MFM94] to medicine [Bai94] and agriculture [WAK94]. The approaches taken vary with the field, depending on the type of models used and the physical data available. For example, the approach taken in [BT86a, BT86b] is to distort parameters in the model until the simulated data matches the actual data. If the distortion required was small then the model was considered validated, and if not it was considered invalidated. This is similar to the approach we adopt, where the model-set is searched for an element that matches the data exactly. However, the model-sets in this thesis are complicated by the fact that they have infinite dimensions.

3.1 Identification

Identification, as described in [Lju87], is the method of constructing a mathematical model from a set of observations, or data. It is usually split into three phases:

1. Obtaining the data record.

2. Selecting a set of models from which the model of the system is to be chosen.

3. Identifying the best model in the model-set, based on the data.

The validation problems we consider clearly lie in phase 3 of this scheme. For example, when validating an LFT model-set, data is assumed to be given and the model-set is LFT (P, Δ, γ). Solving the MVOP for this model-set means finding the “best” model in the model-set, which is the one with the smallest uncertainty bound γ.

3.2. FREQUENCY DOMAIN

Many papers have been published on system identification and it is stated in [ÅE71] that “New methods” are suggested en masse, and, on the surface, the field appears to look more like a bag of tricks than a unified subject.

We will distinguish between different approaches based on the assumptions about the noise in the physical system. The typical assumption [Lju87, Hja93] on a (SISO) physical system is that it is a linear system with the output corrupted by additive noise,

y(t) = G_{0}(z)u(t) + v(t),

where u(t) is the input, y(t) the output, G_{0}(z) the nominal model and v(t) the noise. We will call “standard identification” those techniques with stochastic assumptions on v(t), such as v(t) is a realization of a zero mean, stationary stochastic process with a given spectrum, and call “H_{∞} identification” those techniques with a bound on v(t), such as |v(t)| ≤ δ. For an interesting discussion on the relationship between stochastic and non-stochastic assumptions see [Hja93, ch.1].

Standard identification techniques [vdB93, SVdH92, Cor89, WL92, GGN92] are not suitable for validating the model-sets we consider because of the stochastic assumptions on the noise. In LFT model-sets, the noise is assumed to be bounded in some norm, which is the assumption in “H_{∞} identification”. Many papers have appeared on “H_{∞} identification” [MR85, HJN91, Mäk91, Par91, RL92, GK92, Par93] with the main approach based on the problem introduced in [HJN91]. In this approach models are identified based on a number of points in the frequency response, that are corrupted by additive noise which is bounded in magnitude. With the model comes an H_{∞}-norm bound on the error between the identified model and “true” system, so “H_{∞} identification” may be considered as a method of identifying model-sets of the form of A(P, Δ, γ).

“H_{∞} identification” techniques can be used as a method of model validation. Given a model-set and a set of points in the frequency response, “H_{∞} identification” can be used to identify a model-set. If there is a non-empty intersection between an identified model-set and the assumed model-set then the assumed model-set cannot be invalidated. The problem with identifying model-sets in this way, however, is that the data must consist of several points in the frequency response of the system. In our approach to model validation we would like to consider arbitrary data and more general model-sets.

3.2 Frequency domain

The first study of model validation for LFT model-sets was carried out by Smith and Doyle [Smi90, SD90, SD92]. The analysis was in the frequency domain and the model-sets included block structured uncertainty and noise (LFT (P, X_{K}, γ, δ), with K(0, 0, m)). The model validation problem was broken down into a series of “constant matrix” validation problems at each frequency. Each constant matrix problem was formulated as an optimization by minimizing the norm bound on the noise, subject to the constraint that an element of the model-set interpolated the data. For a small number of uncertainty blocks, the optimization was solved using Lagrange multiplier techniques. For a larger number of uncertainty blocks, upper and lower bounds were obtained using a generalization of the structured singular value, μ.

The model validation problem analysed in [Smi90] has the same interconnection structure as figure 2.2, with the bounds on the uncertainty and noise are normalized to 1. This is achieved by scaling the signals and absorbing the scalings into the model.

Problem 3.1 ([Smi90]) Let P ∈ P^{(n_{w}+n_{y})×(n_{v}+n_{u})}, with μ(P_{11}) ≤ 1, then given a model F_{u}(P, Δ), and an input-output datum (u, y), does there exist (n, Δ), ∥n∥ ≤ 1, Δ ∈ X_{K}, ∥Δ∥_{∞} < 1, such that

y = F_{u}(P, Δ) [ \begin{matrix}n \\ u \end{matrix} ]?

Any pair (n, Δ) meeting the conditions of this problem is said to be admissible.

Smith expressed problem 3.1 as a feasibility problem by proving that an admissible (n, Δ) exists if and only if there exists a noise signal n, and outputs from Δ, such that the model interpolated the data. Define the projection matrix, R_{i}, compatibly with the block structure, as in sub-section 2.3.1, and let S_{i} := R_{i}^{T} R_{i}. Then the main validation result in [Smi90] is

Theorem 3.1 [Smi90] There exists an admissible (n, Δ) for the model validation problem:

y = F_{u}(P, Δ) [ \begin{matrix}n \\ u \end{matrix} ] , ∥n∥ ≤ 1, Δ ∈ X_{K}, ∥Δ∥_{∞} < 1,

if and only if there exists x such that:

1. x^{∗} [ \begin{matrix}S_{i} & 0 \\ 0 & 0 \end{matrix} ] x ≤ [ x^{∗} u^{∗}] P^{∗} [ \begin{matrix}S_{i} & 0 \\ 0 & 0 \end{matrix} ] P [ \begin{matrix}x \\ u \end{matrix} ] , i = 1, 2, …, m.

2. x^{∗} [ \begin{matrix}0 & 0 \\ 0 & I \end{matrix} ] x ≤ 1.

3. y = [ 0 I] P [ \begin{matrix}x \\ u \end{matrix} ] .

It is easy to see how to prove this result; condition 1 corresponds to the constraints implied by the uncertainty structure, condition 2 the constraint implied by the norm bound on the noise, and condition 3 the constraint that the inputs produce the output.

Any feasibility problem can be expressed as an optimization problem by treating a constraint as the objective function. In [Smi90] the norm bound on the noise was taken to be the objective function. The equality constraint (3), in theorem 3.1, was

removed by parameterizing all solutions for x, ie y = [ 0 I] P [ \begin{matrix}x \\ u \end{matrix} ] if and only if

x ∈ χ_{e}. Inequalities (1) and (2) were then restricted to the set of x ∈ χ_{e}. This technique is also used in the following chapter to prove results in the time domain. Theorem 3.1 can thus be written as the following optimization problem:

3.3. TIME DOMAIN

Problem 3.2 ([Smi90])

min_{x∈χ_{e}} f(x) subject to g_{i}(x) ≤ 0, i = 1, 2, …, m,

where

f(x) = x^{∗} [ \begin{matrix}S_{i} & 0 \\ 0 & 0 \end{matrix} ] x

and

g_{i}(x) = x^{∗} [ \begin{matrix}S_{i} & 0 \\ 0 & 0 \end{matrix} ] x − [ x^{∗} u^{∗}] P^{∗} [ \begin{matrix}S_{i} & 0 \\ 0 & 0 \end{matrix} ] P [ \begin{matrix}x \\ u \end{matrix} ] .

For a single uncertainty block (m = 1) the optimization problem was solved using Lagrange multiplier techniques. For m > 1, the optimization problem is not necessarily convex and Lagrange multiplier techniques are not guaranteed to find the optimum point. The conditions required to guarantee that an optimal point can be found are closely related to a generalization of μ, called Ψ_{s}(M, χ).

The structured singular value μ can be difficult to calculate unless there are a small number of uncertainty blocks. However, good upper and lower bounds exist that are readily computable, see [BDG^{+}]. This motivated Smith to try a similar method for the function Ψ_{s}(M, χ), a generalization of μ. To define Ψ_{s}(M, χ) let χ be a given set and let I_{s} and Ī_{s} be two indexing sets that partition the set of uncertainty blocks, that is I_{s} ∪ Ī_{s} = {1, 2, …, m} and I_{s} ∩ Ī_{s} = ϕ. Define Λ := {x : x ∈ χ, ∥R_{i}x∥ ≤ ∥R_{i}Mx∥, i ∈ I_{s}}, then

Ψ_{s}(M, χ) = \{ \begin{matrix}sup_{γ,∥x∥=1,x∈χ} \{γ : \begin{matrix}∥R_{i}x∥ & ≤ & ∥R_{i}Mx∥ & , i ∈ I_{s}, \\ ∥R_{i}x∥γ & ≤ & ∥R_{i}Mx∥ & , i ∈ Ī_{s} \end{matrix} \} , & if Λ ≠ ϕ, \\ 0, & if Λ = ϕ. \end{matrix}

Smith proved that computing Ψ_{s}(M, χ), for a specific choice of M and χ, solves the LFT model validation problem in theorem 3.1. Not surprisingly it is not possible to compute Ψ_{s}(M, χ) for most model validation problems, but upper and lower bounds can be computed, the upper bound being expressible as an LMI [NS91].

3.3 Time domain

A time domain approach to model validation has recently been carried out by Poolla et al [PKT^{+}92, PKT^{+}94], and independently by Zhou and Kimura [ZK92, ZK93]. In this approach, results from Carathéodory-Fejér interpolation theory are used to validate model-sets of the form A(P, Δ, γ), M(P, Δ, γ) and N CF(Ñ, M̃, Δ, γ). The results were proved by calculating the inputs and outputs to the uncertainty, from the input-output data, and applying interpolation results to show a suitable norm bounded uncertainty exists. In the noise free case, necessary and sufficient conditions for validation were expressed as a test on the positive definiteness of a matrix formed from the data. When noise is included the conditions involve the solution of a convex feasibility problem.

We discussed the importance of interpolation results in chapter 1 and stated the solution to the tangential Carathéodory-Fejér interpolation problem. This result is also fundamental in [PKT^{+}92, PKT^{+}94, ZK92, ZK93] and was proved in [PKT^{+}92] and , for SISO systems, [ZK92].

Theorem 3.2 [PKT^{+}92] Given sequences u ∈ π_{l}Sm_{+} and y ∈ π_{l}Sn_{+}, and a positive real number γ, there exists a stable, causal, linear, time-invariant operator Δ satisfying

∥Δ∥_{i2} ≤ γ,

π_{l}Δ [ u_{0} ; u_{1} ; ⋮ ; u_{l−1} ] = [ y_{0} ; y_{1} ; ⋮ ; y_{l−1} ] ,

if and only if

T_{y}^{T} T_{y} ≤ γ^{2}T_{u}^{T} T_{u}.

The analogous result for time-varying uncertainty, proved in [PKT^{+}92], is:

Theorem 3.3 [PKT^{+}92, PKT^{+}94] . Given sequences u ∈ π_{l}Sm_{+} and y ∈ π_{l}Sn_{+}, and a positive real number γ, there exists a stable, causal, linear, time-varying operator Δ satisfying

∥Δ∥_{i2} ≤ γ,

π_{l}Δ [ u_{0} ; u_{1} ; ⋮ ; u_{l−1} ] = [ y_{0} ; y_{1} ; ⋮ ; y_{l−1} ] ,

if and only if

∥π_{k}y∥_{2} ≤ γ∥π_{k}u∥_{2},

for all k = 1, 2, …, l.

The following remark was also stated in [PKT^{+}94].

Remark 3.1 The above result for linear, time-varying operators also holds for nonlinear operators.

This is an extremely important remark. It means that theorem 3.3 can also be used to validate model-sets that have nonlinear norm bounded uncertainty. As modelsets of this form can be guaranteed to be stabilized by certain H_{∞} controllers (see chapter 5), it is important that these sets can be validated.

Interpolation results are also described in [PKT^{+}92, PKT^{+}94] for time-invariant and time-varying operators with a bounded induced-infinity norm. These results are useful for model validation applied to l^{1} control problems [DP87]. However, we only consider the model-sets used in H_{∞}-control so these results will not be needed.

The interpolation results in theorems 3.2 and 3.3 are applied in [PKT^{+}92, PKT^{+}94, ZK92, ZK93, ZK94] to additive, multiplicative and normalized coprime factor model-sets. Additive noise was also considered, and model validation then becomes equivalent to a convex feasibility problem. The result for A(P, T I_{∞}^{p×m}, γ) is

3.4. COMPARISON

Theorem 3.4 [PKT^{+}92] Given the model-set A(P, T I_{∞}^{p×m}, γ) and a set of inputoutput data u ∈ π_{l}Sm_{+}, y ∈ π_{l}Sp_{+} with the output corrupted by an additive disturbance n, which comes from a convex set N, let

û = (û_{0}, û_{1}, …, û_{l−1}) := π_{l}u,
ŷ = (ŷ_{0}, ŷ_{1}, …, ŷ_{l−1}) := y − π_{l}P u.

Then there exists P_{1} ∈ A(P, T I_{∞}^{p×m}, γ) satisfying

y = π_{l}P_{1}u,

if and only if there exists n ∈ π_{l}N such that

σ̄ [(Tŷ − T_{n})(T T ; û Tû)^{\frac{1}{2}} ] ≤ γ.

3.4 Comparison

It is difficult to compare the time and frequency domain approaches as they study different problems. The frequency domain approach of Smith and Doyle studies LFT model-sets with structured uncertainty, whereas the time domain approaches only consider additive, multiplicative and coprime factor uncertainty. It is therefore not surprising that the results in the time domain are stronger and easier to calculate.

The fundamental difference between the time domain and frequency domain approaches is the causality of the uncertainty, Δ. In the frequency domain approach, an admissible Δ need not be causal. It is therefore possible for a model-set to be validated when no causal Δ is admissible. This cannot happen in the time domain approach as causality is guaranteed. In fact, unless additional assumptions are made on the inputs for future time, causality is necessary in the time domain approach. This is because, given a finite set of input-output data, it is easy to construct a noncausal operator of arbitrarily small norm that will interpolate the data, by assuming future inputs are very large.

A limitation of both the frequency domain, and time domain approaches, is that they only consider stable uncertainty. It was demonstrated in the example in chapter 1 that realistic problems may be considered where the uncertainty can be unstable. We will derive results that enable us to validate unstable uncertainty when considering the ν-gap metric in chapter 5.

General Model Validation

In this chapter we consider generalizations of the time-domain validation results described in the previous chapter. The results in [PKT^{+}94] are easily generalized to LFT model-sets with unstructured uncertainty and no noise, unstructured uncertainty with noise and structured uncertainty. However, as the results are generalized, convexity of the corresponding feasibility problem may be lost. We examine conditions under which convexity is preserved and prove that validating a general LFT model-set is NP-hard. This means that the problem of validating a general LFT model-set is computationally at least as hard as a class of problems that are recognized to be computationally demanding.

4.1 LFT model-sets

In this section we prove results for progressively more complicated model-sets. We begin with results for LFT model-sets with unstructured uncertainty, then obtain results for the same model-sets including noise, and finally state a result for LFT model-sets with structured uncertainty. The results are all proved in a similar fashion: For the model-sets considered in [PKT^{+}94, ZK92], the inputs and outputs from

P_{12} P_{22} w  v y  u P_{21} P_{11} Δ

Figure 4.1: LFT without noise

4.1. LFT MODEL-SETS

the uncertainty are uniquely determined (in the absence of noise) from the inputoutput data u and y. For LFT model-sets, shown in figure 4.1, signals v and w are only uniquely determined from u and y if T_{P_{21}} is square and invertible. If T_{P_{21}} is not square and invertible there are a set of v satisfying

T_{P_{21}}v = y − T_{P_{22}}u,  (4.1)

and for each v there is a corresponding w, given by

w = T_{P_{11}}v + T_{P_{12}}u.  (4.2)

Hence the model is validated if there exists a v in the solution space of equation (4.1), and a corresponding w, such that there exists a Δ, with ∥Δ∥_{i2} ≤ γ, that maps w to v.

4.1.1 Unstructured uncertainty

Proposition 4.1 Given a model P ∈ P^{(n_{w}+n_{y})×(n_{v}+n_{u})}, a bound on the induced 2norm of the uncertainty, γ, and a set of input-output data u ∈ π_{l}Sn_{u}_{+} and y ∈ π_{l}Sn_{y}_{+}, there exists Δ ∈ T I_{∞}^{n_{v}×n_{w}} such that the following equations have a solution for w

(I − T_{P_{11}}π_{l}Δ)w = T_{P_{12}}u,
y = T_{P_{21}}π_{l}Δw + T_{P_{22}}u,  (4.3)
∥Δ∥_{i2} ≤ γ,

if and only if there exists v ∈ ℝ^{ln_{v}} satisfying

T_{P_{21}}v = y − T_{P_{22}}u,
T_{v}^{T} T_{v} ≤ γ^{2}T_{P_{11}v+P_{12}u}^{T}T_{P_{11}v+P_{12}u}.

Proof: Given the P, γ, u and y in the proposition statement there exists a Δ ∈ T I_{∞}^{n_{v}×n_{w}} such that equations (4.3) have a solution for w if and only if there exists v ∈ ℝ^{ln_{v}}, and a corresponding w, where

w = T_{P_{11}}v + T_{P_{12}}u,

such that there exists a stable, causal, linear, time-invariant operator Δ satisfying

∥Δ∥_{i2} ≤ γ,

π_{l}Δ [ w_{0} ; w_{1} ; ⋮ ; w_{l−1} ] = [ v_{0} ; v_{1} ; ⋮ ; v_{l−1} ] .

From theorem 3.2, such a Δ exists if and only if

T_{v}^{T} T_{v} ≤ γ^{2}T_{w}^{T} T_{w}.

Remark 4.1 If I − T_{P_{11}}π_{l}Δ is invertible for all Δ ∈ T I_{∞}^{n_{v}×n_{w}}, which would be ensured if, for example, γ < (σ̄(T_{P_{11}}))^{−1}, then proposition 4.1 becomes:

Given a model-set LFT (P, T I_{∞}^{n_{v}×n_{w}}, γ) and a set of input-output data

u ∈ π_{l}Sn_{u}_{+} and y ∈ π_{l}Sn_{y}_{+}, there exists a P_{1} ∈ LFT (P, T I_{∞}^{n_{v}×n_{w}}, γ)

satisfying

y = π_{l}P_{1}u,

if and only if there exists v ∈ ℝ^{ln_{v}} satisfying

T_{P_{21}}v = y − T_{P_{22}}u,
T_{v}^{T} T_{v} ≤ γ^{2}T_{P_{11}v+P_{12}u}^{T}T_{P_{11}v+P_{12}u}.

In general the Δ that is implied by the existence of v ∈ ℝ^{ln_{v}} satisfying

T_{P_{21}}v = y − T_{P_{22}}u,
T_{v}^{T} T_{v} ≤ γ^{2}T_{P_{11}v+P_{12}u}^{T}T_{P_{11}v+P_{12}u},

is not guaranteed to ensure that I − T_{P_{11}}π_{l}Δ is invertible so F_{u}(T_{P} , π_{l}Δ) may not be well-defined. However, if the model is to be used in closed loop with a robustly stabilizing controller then the closed loop will be well-posed if we make assumptions about the existence of a robustly stabilizing controller. This is illustrated by the following proposition.

Proposition 4.2 Given a model P ∈ P^{(n_{w}+n_{y})×(n_{v}+n_{u})} and a bound on the infinity norm of the uncertainty, γ, then if there exists a strictly causal controller K such that I − KP_{22} is invertible and ∥F_{u}(P, K)∥_{i2} < γ^{−1} then the feedback loop in figure 4.2 is well-posed.

Δ P K d_{2}  e_{2} w  e_{1}  d_{1} v

Figure 4.2: Closed loop of uncertain model and controller

Proof: The closed loop equations of figure 4.2 are

[ \begin{matrix}I & −Δ & 0 & 0 \\ −P_{11} & I & 0 & −P_{12} \\ −P_{21} & 0 & I & −P_{22} \\ 0 & 0 & −K & I \end{matrix} ] [ \begin{matrix}v \\ w \\ e_{2} \\ e_{1} \end{matrix} ] = [ \begin{matrix}0 & 0 \\ 0 & 0 \\ 0 & I \\ I & 0 \end{matrix} ] [ \begin{matrix}d_{1} \\ d_{2} \end{matrix} ] .

Let A be the matrix on the left hand side of this equation and partition A into four submatrices so that

A_{22} = [ \begin{matrix}I & −P_{22} \\ −K & I \end{matrix} ] .

By assumption, A_{22} is invertible and using a Schur complement formula

det A = det(I − P_{22}K) det(I − ΔF_{l}(P, K)).

Hence if ∥F_{u}(P, K)∥_{i2} < γ^{−1} then the feedback loop in figure 4.2 is well-posed.

This result tells us that although we can use proposition 4.1 to test for the existence of Δ, ∥Δ∥_{i2} ≤ γ, for any γ the answer only makes sense if γ is small enough that a robustly stabilizing controller exists. This is true for all the results of this section. We continue by stating the result for validating models with timevarying uncertainty.

Proposition 4.3 Given a model P ∈ P^{(n_{w}+n_{y})×(n_{v}+n_{u})}, a bound on the induced 2norm of the uncertainty, γ, and a set of input-output data u ∈ π_{l}Sn_{u}_{+} and y ∈ π_{l}Sn_{y}_{+}, there exists Δ ∈ T V_{∞}^{n_{v}×n_{w}} such that the following equations have a solution for w

(I − T_{P_{11}}π_{l}Δ)w = T_{P_{12}}u,
y = T_{P_{21}}π_{l}Δw + T_{P_{22}}u,  (4.4)
∥Δ∥_{i2} ≤ γ,

if and only if there exists v ∈ ℝ^{ln_{v}} satisfying

T_{P_{21}}v = y − T_{P_{22}}u,  (4.5)
∥π_{k}v∥_{2} ≤ γ∥π_{k}(T_{P_{11}}v + T_{P_{12}}u)∥_{2},  (4.6)

for all k = 1, 2, …, l.

Proof: The proof is similar to the proof of proposition 4.1, using theorem 3.3 to obtain necessary and sufficient conditions for the existence of a suitable time-varying Δ.

4.1.2 Including noise

We state results for LFT model-set validation with a 2-norm bound on the noise, but the results are identical if the 2-norm is replaced by a different norm throughout. Results using an infinity-norm bound on the noise will be used in chapter 5.

Proposition 4.4 Given a model P ∈ P^{(n_{w}+n_{y})×(n_{v}+n_{n}+n_{u})}, a bound on the induced 2-norm of the uncertainty, γ, a bound on the 2-norm of the noise, δ, and a set of input-output data u ∈ π_{l}Sn_{u}_{+} and y ∈ π_{l}Sn_{y}_{+}, there exists Δ ∈ T I_{∞}^{n_{v}×n_{w}} such that the following equations have a solution for w and n

(I − T_{P_{11}}π_{l}Δ)w = T_{P_{12}}n + T_{P_{13}}u,
y = T_{P_{21}}π_{l}Δw + T_{P_{22}}n + T_{P_{23}}u,  (4.7)
∥Δ∥_{i2} ≤ γ,
∥n∥_{2} ≤ δ  (4.8)

if and only if there exists [ \begin{matrix}v \\ n \end{matrix} ] ∈ ℝ^{l(n_{v}+n_{n})} satisfying

[ T_{P_{21}} T_{P_{22}} ] [ \begin{matrix}v \\ n \end{matrix} ] = y − T_{P_{23}}u,
T_{v}^{T} T_{v} ≤ γ^{2}T_{P_{11}v+P_{12}n+P_{13}u}^{T}T_{P_{11}v+P_{12}n+P_{13}u},
∥n∥_{2} ≤ δ.

Proof: Given the P, γ, δ, u and y in the proposition statement there exists a Δ ∈ T I_{∞}^{n_{v}×n_{w}} such that equations (4.7) have a solution for w and n if and only if there exists v ∈ ℝ^{ln_{v}}, and a corresponding w, where

w = T_{P_{11}}v + T_{P_{12}}n + T_{P_{13}}u,

such that there exists a stable, causal, linear, time-invariant operator Δ satisfying

∥Δ∥_{i2} ≤ γ,

π_{l}Δ [ w_{0} ; w_{1} ; ⋮ ; w_{l−1} ] = [ v_{0} ; v_{1} ; ⋮ ; v_{l−1} ] .

From theorem 3.2, such a Δ exists if and only if

T_{v}^{T} T_{v} ≤ γ^{2}T_{w}^{T} T_{w}.

Proposition 4.5 Given a model P ∈ P^{(n_{w}+n_{y})×(n_{v}+n_{n}+n_{u})}, a bound on the induced 2-norm of the uncertainty, γ, a bound on the 2-norm of the noise, δ, and a set of input-output data u ∈ π_{l}Sn_{u}_{+} and y ∈ π_{l}Sn_{y}_{+}, there exists Δ ∈= T V_{∞}^{n_{v}×n_{w}} such that the following equations have a solution for w and n

(I − T_{P_{11}}π_{l}Δ)w = T_{P_{12}}n + T_{P_{13}}u,
y = T_{P_{21}}π_{l}Δw + T_{P_{22}}n + T_{P_{23}}u,  (4.9)
∥Δ∥_{i2} ≤ γ,
∥n∥_{2} ≤ δ  (4.10)

if and only if there exists [ \begin{matrix}v \\ n \end{matrix} ] ∈ ℝ^{l(n_{v}+n_{n})} satisfying

[ T_{P_{21}} T_{P_{22}} ] [ \begin{matrix}v \\ n \end{matrix} ] = y − T_{P_{23}}u,
∥π_{k}v∥_{2} ≤ γ∥π_{k}([ T_{P_{11}} T_{P_{12}} ] [ \begin{matrix}v \\ n \end{matrix} ] + T_{P_{13}}u)∥_{2},
∥n∥_{2} ≤ δ,

for all k = 1, 2, …, l.

Proof: The proof is again similar to the previous one, with theorem 3.3 used to obtain necessary and sufficient conditions for the existence of a suitable time-varying Δ.

4.1.3 Structured uncertainty

If the uncertainty in the model-sets has additional structure, such as that described in sub-section 2.3.1, the model-sets can be validated by partitioning v and w conformally with the block structure. The projection operator R_{i}, defined in subsection 2.3.1, can be used for this.

Abusing notation for the statement for the following result only, let v_{i} := R_{i}v and w_{i} := R_{i}w. Note that v_{i} and w_{i} are now vectors. For simplicity only the result for time-invariant uncertainty without noise is stated here. It is assumed that there are no repeated scalar complex blocks, so the block structure is defined by K(m_{r}, 0, m_{C}). This simplifies the problem, but still allows for real parameter variations and structured dynamic uncertainty. The result is:

Proposition 4.6 Given a model P ∈ P^{(n_{w}+n_{y})×(n_{v}+n_{u})}, a bound on the induced 2norm of the uncertainty, γ, and a set of input-output data u ∈ π_{l}Sn_{u}_{+} and y ∈ π_{l}Sn_{y}_{+}, there exists Δ ∈ ST I_{K}(∞) such that the following equations have a solution for w

(I − T_{P_{11}}π_{l}Δ)w = T_{P_{12}}u,
y = T_{P_{21}}π_{l}Δw + T_{P_{22}}u,
∥Δ∥_{i2} ≤ γ,

if and only if there exists v ∈ ℝ^{ln_{v}} and δ_{1}, δ_{2}, …, δ_{m_{r}} ∈ ℝ satisfying

T_{P_{21}}v = y − T_{P_{22}}u,
v_{1} = δ_{1}w_{1}, |δ_{1}| ≤ γ,
v_{2} = δ_{2}w_{2}, |δ_{2}| ≤ γ,
⋮
v_{m_{r}} = δ_{m_{r}}w_{m_{r}}, |δ_{m_{r}}| ≤ γ,
T_{v_{mr+1}}^{T}T_{v_{mr+1}} ≤ γ^{2}T_{w_{mr+1}}^{T}T_{w_{mr+1}},
T_{v_{mr+2}}^{T}T_{v_{mr+2}} ≤ γ^{2}T_{w_{mr+2}}^{T}T_{w_{mr+2}},
⋮
T_{v_{mr+m_{C}}}^{T} T_{v_{mr+m_{C}}} ≤ γ^{2}T_{w_{mr+m_{C}}}^{T} T_{w_{mr+m_{C}}} ,

where w = T_{P_{11}}v + T_{P_{12}}u.

Proof: The result can be proved by considering each section of v and w, partitioned conformally with the block structure, separately.

Note that this result can easily be generalized to include noise, using the method of proof for proposition 4.4.

These results have generalized the validation results in [PKT^{+}94, ZK92] to include LFT model-sets. Proposition 4.6, after being extended to include noise, solves a more general problem than is considered in [Smi90], since real uncertainty is included. The result is also less conservative than the results in [Smi90], since the uncertainty is guaranteed to be causal. We will now consider the computational complexity of the various generalizations.

4.2. COMPUTATIONAL COMPLEXITY

4.2 Computational complexity

The main result of this section is that validating a general LFT of P on Δ is at least as hard as solving a class of problems for which it is commonly believed no efficient algorithm exists. We prove this result using techniques from computational complexity theory. In the language of computational complexity we prove that solving an MVDP for a general LFT of P on Δ is NP-hard in the size of uncertainty and number of noise signals. An NP-hard problem is widely recognized to be computationally demanding.

Certain problems in mathematics, such as linear programming problems, can be solved quickly by computer even if the problem has a large number of variables. Problems like these have efficient algorithms that will solve any instance of the problem. For other problems however, it can be proved that no general algorithm exists. These problems are be said to be undecidable. An important example of this was shown by Turing [GJ79]; “that it is impossible to specify any algorithm which, given an arbitrary computer program and an arbitrary input to that program, can decide whether or not the program will eventually halt when applied to that input.”

Unfortunately proving that a given problem has an efficient algorithm, or that one does not exist, is often difficult. Usually, to prove that a problem has an efficient algorithm, it is necessary to construct the algorithm itself. There are many problems for which no efficient algorithm has been found, but proving that no algorithm exists has also proved unsuccessful. A large proportion of computational complexity theory is centred on proving that problems are equivalent in their complexity. This means that if an algorithm can be found for one, one can be found for the rest, and if no algorithm can be found for one, then no algorithm can be found for the remainder.

4.2.1 Central concepts

The ideas described above have been formalized mathematically, using Turing machines as the idealization of a computer (see [GJ79]), but this level of abstraction is not necessary to prove model validation is NP-hard. All the problems considered in complexity theory are decision problems, which means that have a binary answer “yes” or “no”. A famous problem studied in complexity theory is the Travelling Salesman problem. It is stated as, “Given n cities, and the distances between them, what is the shortest path such that every city is visited once only?” This is expressed as a decision problem by asking the question, “Is there a path passing through each city once whose total length is less than some given number k?” Answering either of these problems will enable the other to be answered. The relationship between these two versions of the Travelling Salesman problem is the same as between the MVDP and MVOP, described in section 2.4.

One of the key notions from complexity theory is that of a polynomial time algorithm and exponential time algorithm. A polynomial time algorithm is an algorithm whose execution time is bounded by a polynomial function of the number of variables, for example the time could be bounded by 10^{6}n^{5} + 10^{−3}n^{10}, where n is the number of variables. An exponential time algorithm is one whose execution time cannot be bounded by a polynomial function, for example bounded by 10^{−5}e^{n} + 10^{8}n^{42}. Polynomial time algorithms are considered to be efficient whereas exponential time algorithms are not. This notion is also used to define what is tractable and what is not. If there exists a polynomial time algorithm for a problem then it is said to be tractable^{1}. If a polynomial time algorithm does not exist then the problem is said to be intractable. No polynomial time algorithm has been found for the Travelling Salesman problem, but neither has it been proved that a polynomial time algorithm does not exist. However, it is widely believed that no polynomial time algorithm exists.

To describe the sets of problems that are considered intractable, let P denote the class of problems solvable in polynomial time, and NP the class of problems solvable by a non-deterministic Turing machine in polynomial time. Non-deterministic Turing machines are mathematical idealizations of computers that cannot be built. Roughly speaking the class NP contains all the problems whose solution can be checked in polynomial time. For example, the Travelling Salesman problem is an element of NP, since given a number k, and a sequence of n cities, it can clearly be “checked” in polynomial time whether or not the total distance is less then k. NP also contains many other problems for which it is widely believed no polynomial time algorithm exists.

The two concepts used to determine how different problems are related in their computational complexity are polynomial reducibility and polynomial equivalence. A polynomial reduction is a transformation from one problem to another, that takes polynomial time. Thus if problem A is polynomially reduced to problem B, and a polynomial time algorithm exists for problem B, then a polynomial time algorithm must also exist for problem A. If problem A can be polynomially reduced to problem B, and vice versa, then problems A and B are polynomially equivalent.

Within NP are the class of NP-complete problems. These are the class of problems that are polynomially equivalent to each other, and such that every problem in NP can be polynomially reduced to an NP-complete problem. Hence if a polynomial time algorithm is discovered for any NP-complete problem, then all the problems in NP can be solved in polynomial time. Also if any problem in NP is proved to be intractable, then so is every NP-complete problem. The Travelling Salesman problem is NP-complete.

Closely related to the class of NP-complete problems is the class of NP-hard problems. These are the problems that, whether or not they belong to NP, may be polynomially reduced from an NP-complete problem. Thus NP-hard problems are at least as hard as NP-complete problems. If a problem is proved to be NP-hard and a polynomial time algorithm is found for it, then there exists a polynomial time algorithm for every problem in NP. However, if a polynomial time algorithm is found for every element of NP, there is not necessarily a polynomial time algorithm for the NP-hard problem.

4.2.2 Model validation is NP-hard

We will now prove that the general model validation problem with unstructured uncertainty and noise bounded in its infinity-norm is NP-hard in the size of the uncertainty and number of noise signals. We prove the result in a similar way to which Braatz et al proved that calculating μ is NP-hard in the number of repeated real uncertainty blocks in [BYDM94]. They proved the result by showing that a general set of non-convex quadratic programming problems, which are known to be NP-hard from [MK87], are a special case of a μ problem with repeated real uncertainty blocks. We show that the special set of non-convex quadratic programming problems proved to be NP-hard in [MK87] are a special case of a model validation problem with unstructured uncertainty and sensor noise with bounded infinity-norm.

1This definition of intractable is taken from [GJ79]. It is only a rough approximation of the dictionary meaning.

The key result from [MK87] is

Theorem 4.1 ([MK87]) The following problem is NP-hard:

Given positive integers d_{0}; d_{1}, …, d_{n}, is

min_{0≤x_{i}≤1} (∑_{i=1}^{n} d_{i}x_{i} − d_{0})^{2} + ∑_{i=1}^{n} x_{i}(1 − x_{i}) ≤ 0?

With the following definitions,

d := [ d_{1} ; ⋮ ; d_{n} ] , c := [ 1 ⋮ 1 ] ∈ ℝ^{n},

z_{i} := x_{i} − \frac{1}{2}, z := [ z_{1} ; ⋮ ; z_{n} ] ,

we get the corollary

Corollary 4.1 The following problem is NP-hard:

Given positive integers d_{0}; d_{1}, …, d_{n}, is

min_{∥z∥_{∞}≤ \frac{1}{2}} z^{T} (dd^{T} − I)z + (c^{T} dd^{T} − 2d_{0}d^{T})z + (\frac{1}{2}d^{T} c − d_{0})^{2} + \frac{1}{4}c^{T} c ≤ 0?

The set of quadratic programs in this corollary are a special case of a model validation problem, which proves model validation is NP-hard. This is demonstrated by the main result of this section:

Theorem 4.2 Let Δ = ℝ^{m} then the following problem is NP-hard:

Given positive integers d_{0}; d_{1}, …, d_{m} let

P := [ \begin{matrix}P_{11} & P_{12} & P_{13} \\ P_{21} & P_{22} & P_{23} \end{matrix} ] = [ \begin{matrix}0 & I & 0 \\ −1 & d^{T} & \frac{1}{2}c^{T} d − d_{0} \end{matrix} ] ,

u := 1, y := 0,

does there exist n ∈ ℚ^{m}, Δ ∈ Δ, ∥Δ∥_{i2} ≤ 1, ∥n∥_{∞} ≤ \frac{1}{2} satisfying

y = F_{u}(P, Δ) [ \begin{matrix}n \\ u \end{matrix} ]?

4.3. CONVEXITY

Proof: Given the P, u and y of the theorem statement there exists n ∈ ℚ^{m}, Δ ∈ Δ, ∥Δ∥_{i2} ≤ 1, ∥n∥_{∞} ≤ \frac{1}{2} satisfying

y = F_{u}(P, Δ) [ \begin{matrix}n \\ u \end{matrix} ]

if and only if ∥v∥_{2}^{2} ≤ ∥w∥_{2}^{2} and ∥n∥_{∞} ≤ \frac{1}{2}, where w = n and v = d^{T} n + \frac{1}{2}c^{T} d − d_{0}. This inequality is true if and only if

min_{∥n∥_{∞}≤ \frac{1}{2}} n^{T} (dd^{T} − I)n + (c^{T} dd^{T} − 2d_{0}d^{T})n + (\frac{1}{2}d^{T} c − d_{0})^{2} ≤ 0

Solving this minimisation is clearly equivalent to the quadratic program in the corollary, which implies the result.

4.2.3 Comments

Theorem 4.2 is a powerful result and means any generalization of the model validation problem in theorem 4.2, which is a large class of model validation problems, is also NP-hard. Hence it is unlikely that a polynomial time algorithm can be found to solve a general LFT model validation problem. This means that if we want to solve model validation problems efficiently, we should look at approximate methods, or special cases of the general model validation problem. Some approximate methods for this problem, based on the function Ψ_{s}(M, χ) defined in section 3.2, are described in [Smi90].

We have already seen by the results in [PKT^{+}94] that there exist model validation problems that are equivalent to convex programming problems so are not NP-hard. We will concentrate on these model validation problems for the remainder of the thesis.

4.3 Convexity

In this section we derive conditions under which the results in section 4.1 are equivalent to convex feasibility problems. In the discussion following these results, we indicate why many model validation problems may automatically be convex. The motivation for using convexity as a method of assessing the computational properties of a problem was discussed in chapter 1.

Before deriving the results it is useful to define a convex set, function and problem. A set C is convex if for every x_{1}, x_{2} ∈ C, and every t ∈ [0, 1], the point tx_{1} + (1 − t)x_{2} ∈ C. A function f, defined on a convex set C, is convex if for every x_{1}, x_{2} ∈ C, and every t ∈ [0, 1],

f(tx_{1} + (1 − t)x_{2}) ≤ tf(x_{1}) + (1 − t)f(x_{2}).

A feasibility problem is convex if the set of feasible points is convex. Note that if f is a convex function, then the set of x satisfying f(x) ≤ c, for some constant c, is a convex set. Hence, a sufficient condition for a feasibility problem to be convex, is that the functions defining the constraints are convex. We will use this sufficient condition for the results in this section.

4.3.1 Unstructured uncertainty

Deriving conditions for the convexity of model validation problems is easier when the uncertainty is time-varying. The following result gives sufficient conditions for proposition 4.5, which includes proposition 4.3 as a special case, to be equivalent to a convex feasibility problem.

Theorem 4.3 Given a model-set LFT (P, T V_{∞}^{n_{v}×n_{w}}, γ, δ) and a set of input-output data u ∈ π_{l}Sn_{u}_{+} and y ∈ π_{l}Sn_{y}_{+}, define the matrix V_{k} as follows: If equation (4.11) has no solution let V_{k} be the k(n_{v} + n_{n}) × 1 matrix of zeros. Otherwise let V be the matrix whose columns span the nullspace of [ T_{P_{21}} T_{P_{22}} ].

Partition V as [ \begin{matrix}V_{1} \\ V_{2} \end{matrix} ] where V_{1} has ln_{v} rows and V_{2} has ln_{n} rows. Let V_{1k} be the

first kn_{v} rows of V_{1}, and V_{2k} be the first kn_{n} rows of V_{2}, then V_{k} := [ \begin{matrix}V_{1k} \\ V_{2k} \end{matrix} ] .

The the set of [ \begin{matrix}v \\ n \end{matrix} ] ∈ ℝ^{l(n_{v}+n_{n})} satisfying

[ T_{P_{21}} T_{P_{22}} ] [ \begin{matrix}v \\ n \end{matrix} ] = y − T_{P_{23}}u,  (4.11)
∥π_{k}v∥_{2} ≤ γ∥π_{k}([ T_{P_{11}} T_{P_{12}} ] [ \begin{matrix}v \\ n \end{matrix} ] + T_{P_{13}}u)∥_{2},  (4.12)
∥n∥_{2} ≤ δ,  (4.13)

for all k = 1, 2, …, l, is convex if

V_{k}^{T} [ \begin{matrix}γ^{−2}I_{kn_{v}} − π_{k}T_{P_{11}}^{T}π_{k}T_{P_{11}} & −π_{k}T_{P_{11}}^{T}π_{k}T_{P_{12}} \\ −π_{k}T_{P_{12}}^{T}π_{k}T_{P_{11}} & −π_{k}T_{P_{12}}^{T}π_{k}T_{P_{12}} \end{matrix} ] V_{k} ≥ 0,

for all k = 1, 2, …, l.

Proof: The result is proved by deriving necessary and sufficient conditions for the inequalities (4.12) to be convex on the space of all solutions to equation (4.11).

If equation (4.11) has no solution then no [ \begin{matrix}v \\ n \end{matrix} ] satisfies the feasibility problem.

In this case V_{k} is a matrix of zeros and the proposition is vacuously true as the set of feasible points is empty. Consequently assume that equation (4.11) has a solution

and let [ \begin{matrix}v_{0} \\ n_{0} \end{matrix} ] be the minimal norm solution that is orthogonal to the nullspace of

[ T_{P_{21}} T_{P_{22}} ]. Then all the solutions of equation (4.11) are of the form

[ \begin{matrix}v_{0} \\ n_{0} \end{matrix} ] ⊕ V_{l}x,

where x is a vector with dimensions the same as the nullspace of [ T_{P_{21}} T_{P_{22}} ].

Substituting the solution of equation (4.11) into inequalities (4.12), the quadratic term in the k-th inequality is

x^{T} [ V_{1k}^{T} V_{2k}^{T}] \{[ \begin{matrix}γ^{−2}I & 0 \\ 0 & 0 \end{matrix} ] − [ π_{k}T_{P_{11}}^{T} π_{k}T_{P_{12}}^{T} ] [ \begin{matrix}π_{k}T_{P_{11}} \\ π_{k}T_{P_{12}} \end{matrix} ]\} [ \begin{matrix}V_{1k} \\ V_{2k} \end{matrix} ] x.

The set of x satisfying each inequality is convex if the matrix sandwiched between the x^{T} and x is positive semi-definite (apply the results in [Lue84, p.176-178]). As the inequality (4.13) is always convex the result is proved.

The conditions for convexity given by this theorem appear to be very restrictive, but meaningful model validation problems automatically satisfy the sufficient conditions, as is demonstrated by the following example. Before stating the example it is necessary to recall two facts about block triangular matrices. Firstly, the inverse of a block triangular matrix is itself block triangular, ie if A ∈ ℝ^{n×n} and C ∈ ℝ^{m×m} are both invertible,

[ \begin{matrix}A & 0 \\ B & C \end{matrix} ]^{−1} = [ \begin{matrix}A^{−1} & 0 \\ −C^{−1}BA^{−1} & C^{−1} \end{matrix} ] .

Secondly, if X ∈ T I_{l}^{m×m}, then for k ≤ l, π_{k}X^{−1} = (π_{k}X)^{−1}.

Example 4.1 Consider N CF(Ñ, M̃, Δ, γ, δ), as defined in sub-section 2.3.2. This can be written as an LFT on [ ΔÑ ΔM̃ ] by taking

P = [ \begin{matrix}P_{11} & P_{12} & P_{13} \\ P_{21} & P_{22} & P_{23} \end{matrix} ] = [ \begin{matrix}0 & 0 & I \\ TM̃^{−1} & TM̃^{−1} & T_{M̃^{−1}Ñ} \\ TM̃^{−1} & TM̃^{−1} & T_{M̃^{−1}Ñ} \end{matrix} ] .

Using this P, equation (4.11) becomes

[ TM̃^{−1} TM̃^{−1} ] [ \begin{matrix}v \\ n \end{matrix} ] = y − T_{M̃^{−1}Ñ}u.

The matrix V, defined in theorem 4.3, spans the nullspace of [ TM̃^{−1} TM̃^{−1} ], so a possible V is

V = [ \begin{matrix}I \\ −I \end{matrix} ] .

Hence a simple calculation, using the two facts about the inverses of block triangular matrices, means the necessary and sufficient conditions of theorem 4.3 become γ^{2}I ≥ 0. This is clearly satisfied for any γ, so validating N CF(Ñ, M̃, Δ, γ, δ) can be achieved by solving a convex feasibility problem.

A simple sufficient condition for proposition 4.3, the case with time-varying uncertainty and no noise, to be equivalent to a convex feasibility problem is given by the following corollary of theorem 4.3.

Corollary 4.2 Given a model-set LFT (P, T V_{∞}^{n_{v}×n_{w}}, γ) and a set of input-output data u ∈ π_{l}Sn_{u}_{+} and y ∈ π_{l}Sn_{y}_{+}, the set of v ∈ ℝ^{ln_{v}} satisfying

T_{P_{21}}v = y − T_{P_{22}}u,
∥π_{k}v∥_{2} ≤ γ∥π_{k}(T_{P_{11}}v + T_{P_{12}}u)∥_{2},

for all k = 1, 2, …, l, is convex if σ̄(π_{l}T_{P_{11}}) < γ^{−1}.

The sufficient condition given by this corollary is also a sufficient condition for the LFT F_{u}(P, Δ) to be well-posed, and will be discussed later in this section.

Determining conditions for LFT model-sets with time-invariant uncertainty, stated in propositions 4.1 and 4.4, to be equivalent to a convex feasibility problem is more difficult. Recall the necessary and sufficient conditions in proposition 4.1

for validation; the existence of [ \begin{matrix}v \\ n \end{matrix} ] ∈ ℝ^{l(n_{v}+n_{n})} satisfying

[ T_{P_{21}} T_{P_{22}} ] [ \begin{matrix}v \\ n \end{matrix} ] = y − T_{P_{23}}u,  (4.14)
T_{v}^{T} T_{v} ≤ γ^{2}T_{P_{11}v+P_{12}n+P_{13}u}^{T}T_{P_{11}v+P_{12}n+P_{13}u},  (4.15)
∥n∥_{2} ≤ δ.  (4.16)

It is clear that this set of [ \begin{matrix}v \\ n \end{matrix} ] will be convex if and only if inequality (4.15) forms

a convex constraint on the set of solutions to equation (4.14). However, it is not clear how to express this condition as a simple condition in terms of V, as it was in theorem 4.3. Consequently only sufficient conditions are obtained for convexity of the result in proposition 4.1:

Proposition 4.7 Given a model-set LFT (P, T V_{∞}^{n_{v}×n_{w}}, γ) and a set of input-output data u ∈ π_{l}Sn_{u}_{+} and y ∈ π_{l}Sn_{y}_{+}, the set of v ∈ ℝ^{ln_{v}} satisfying

T_{P_{21}}v = y − T_{P_{22}}u,  (4.17)
T_{v}^{T} T_{v} ≤ γ^{2}T_{P_{11}v+P_{12}w}^{T}T_{P_{11}v+P_{12}w},  (4.18)

is convex if σ̄(π_{l}T_{P_{11}}) < γ^{−1}.

Proof: Assume σ̄(P_{11}) < γ^{−1}, then (γ^{−2}I − P_{11}^{T}P_{11}) is invertible and using Schur complements ([HJ85, p.472]) inequality (4.18) can be written as

[ \begin{matrix}T_{v}^{T} T_{P_{11}}^{T}T_{P_{12}}T_{u} − T_{u}^{T} T_{P_{12}}^{T}T_{P_{11}}T_{v} − T_{u}^{T} T_{P_{12}}^{T}T_{P_{12}}T_{u} & T_{v}^{T} \\ T_{v} & (γ^{−2}I − T_{P_{11}}^{T}T_{P_{11}})^{−1} \end{matrix} ] ≥ 0.

This is an LMI and forms a convex constraint on the set of possible v’s, [BEFB94].

The sufficient condition for convexity, σ̄(π_{l}T_{P_{11}}) < γ^{−1}, may often be assumed when the model P includes a robust controller, since F_{u}(P, Δ) is stable for all Δ, σ̄(Δ) ≤ γ, if and only if σ̄(P_{11}) < γ^{−1}. We would usually require that F_{u}(P, Δ) is stable for all Δ, σ̄(Δ) ≤ γ so one may argue that we would not want to validate a model when σ̄(π_{l}T_{P_{11}}) ≥ γ^{−1}.

When noise is included, assuming σ̄(π_{l}T_{P_{11}}) < γ^{−1} is no longer sufficient to ensure convexity of model validation. There is no easy solution to this problem and it also occurs when considering validation in the ν-gap metric in chapter 5. Convexity is ensured in this problem by assuming an LFT model where the equations are

w = P_{11}v + P_{11}n + P_{13}u,
y = P_{21}v + P_{21}n + P_{23}u,

and P_{21} is square and invertible. However, these assumptions may not be realistic for many problems.

4.3.2 Structured uncertainty

Allowing uncertainty with additional structure can make model validation difficult, as was demonstrated by proving that structured model validation is NP-hard. Not surprisingly the extra structure can easily destroy the convexity properties possessed by the unstructured case:

Proposition 4.8 There exists a model-set LFT (P, ST I_{K}(∞), γ), and a set of input-output data u ∈ π_{l}Sn_{u}_{+} and y ∈ π_{l}Sn_{y}_{+}, such that the set of v ∈ ℝ^{ln_{v}} satisfying the necessary and sufficient conditions of proposition 4.6 is not convex.

4.4. COMPUTATIONAL ISSUES

Proof: Let the model P be

P = [ \begin{matrix}P_{11} & P_{12} \\ P_{21} & P_{22} \end{matrix} ] = [ \begin{matrix}0 & 0 & 0 & 0 \\ 0.5 & 0 & 0 & 1 \\ 0 & 0.5 & 0 & 0 \\ 1 & 0 & 0 & 1 \end{matrix} ] ∈ T I_{1}^{(3+1)×(3+1)},

the uncertainty of the form diag(δ_{1}, δ_{2}, δ_{3}), for δ_{1}, δ_{2} and δ_{3} ∈ ℝ, and the input u = 1 and output y = 0. Let γ = 1 then the equality constraint T_{P_{21}}v = y − T_{P_{22}}u is

[ 1 0 0 ] [ \begin{matrix}v_{1} \\ v_{2} \\ v_{3} \end{matrix} ] = −1.

Hence v_{1} = 0 and v_{2} and v_{3} are not constrained. Calculating w,

w = T_{P_{11}}v + T_{P_{12}}u = [ \begin{matrix}0 \\ 0.5v_{1} + 1 \\ 0.5v_{2} \end{matrix} ] .

The model-set is valid, for this uncertainty structure, if |v_{i}| ≤ |w_{i}| for i = 1, 2, 3. For i = 1, 2 these inequalities form convex constraints but for i = 3 the constraint is

|v_{3}| ≤ |0.5v_{2}|.

A simple sketch shows that the set of v satisfying this constraint is not convex, which proves the result.

Note that in this proof σ̄(P_{11}) = 0.5, which is less than γ^{−1}, so the sufficient conditions for unstructured uncertainty in proposition 4.7 and corollary 4.2 are satisfied, yet the problem with structured uncertainty is still not convex.

4.4 Computational Issues

In this section we describe how the feasibility problems implied by different model validation problems can be implemented on a computer. The discussion is restricted to LFT model-sets. The special case of NCF model-sets is discussed in chapter 5 and further problem specific refinements are discussed in chapters 6 and 7.

Only convex feasibility problems are considered, for the reasons discussed in chapter 1. Hence the model-sets are assumed to have unstructured uncertainty and to satisfy the conditions for convexity described in the previous section. Despite this assumption, experience has shown that validating a model-set by solving a convex feasibility problem does not necessarily imply a solution can be computed with present computer resources. Fortunately, however, in many useful cases it does.

We defined the two fundamental validation problems, MVDP and MVOP, in section 2.4. All of programs used to calculate the numerical results in the thesis have solved the MVOP, the problem of finding the smallest γ. The smallest γ that solves an MVOP will be denoted by γ_{min}(T I), when the uncertainty is time-invariant, and γ_{min}(T V) when the uncertainty is time-varying. For LFT model-sets, γ_{min}(·) is calculated with a bisection search on γ, solving an MVDP at each iteration. For NCF model-sets it is possible to calculate γ_{min}(·) without iterating on γ, see sub-section 5.1.2.

4.4.1 Time-varying uncertainty

When the uncertainty is time-varying, LFT model-set validation (including noise)

can be accomplished by testing the existence of [ \begin{matrix}v \\ n \end{matrix} ] ∈ ℝ^{l(n_{v}+n_{n})} satisfying

[ T_{P_{21}} T_{P_{22}} ] [ \begin{matrix}v \\ n \end{matrix} ] = y − T_{P_{23}}u,  (4.19)
∥π_{k}v∥_{2} ≤ γ∥π_{k}([ T_{P_{11}} T_{P_{12}} ] [ \begin{matrix}v \\ n \end{matrix} ] + T_{P_{13}}u)∥_{2},  (4.20)
∥n∥_{2} ≤ δ,  (4.21)

for all k = 1, 2, …, l. Using the definition of the 2-norm, (4.20) can be expanded, on the space of all solutions to equation (4.19), into a quadratic inequality of the form

[ v^{T} n^{T} ] Q_{k} [ \begin{matrix}v \\ n \end{matrix} ] + 2p_{k}^{T} [ \begin{matrix}v \\ n \end{matrix} ] + r_{k} ≤ 0.  (4.22)

A sufficient condition for the feasibility problem to be convex is that all the Q_{k} are

positive semi-definite. If this is the case, the set of [ \begin{matrix}v \\ n \end{matrix} ] satisfying expressions of the

form of (4.22) are hyperellipsoids, ie ellipsoids in ℝ^{l(n_{v}+n_{n})}. The constraint on n, (4.21), is also restricts n to lie in a hyperellipsoid, but it is more naturally thought of as a ball^{2} in ℝ^{ln_{n}}. Hence the MVDP is equivalent to testing whether or not the intersection of a series of hyperellipsoids is empty.

This type of feasibility problem can be solved using a cutting plane algorithm [Lue84, p.416]. This is an iterative method that approximates the convex constraints, which are hyperellipsoids in the above problem, by hyperplanes and solves a linear program at each iteration. Under mild assumptions on the constraints, this algorithm is globally convergent [Lue84, p.419] and will converge whether the problem is feasible or not. This is ideal for the model validation problem.

A cutting plane algorithm has been implemented to solve MVDP’s for LFT model-sets without noise. This has worked well with simulated data. It requires few iterations to converge (typically less than 10) and will solve problems with up to 500 variables in under an hour when programmed in MATLAB on a SPARC-10. Using this algorithm γ_{min}(T V) can be calculated using a bisection search.

4.4.2 Time-invariant uncertainty

When the uncertainty is time-invariant, LFT model-set validation (including noise)

can be accomplished by testing the existence of [ \begin{matrix}v \\ n \end{matrix} ] ∈ ℝ^{l(n_{v}+n_{n})} satisfying

[ T_{P_{21}} T_{P_{22}} ] [ \begin{matrix}v \\ n \end{matrix} ] = y − T_{P_{23}}u,  (4.23)
T_{v}^{T} T_{v} ≤ γ^{2}T_{P_{11}v+P_{12}n+P_{13}u}^{T}T_{P_{11}v+P_{12}n+P_{13}u},  (4.24)
∥n∥_{2} ≤ δ.

2Note that if the 2-norm was replaced by the infinity-norm the same would be true. However, the set of n satisfying this constraint would be more intuitively thought of as a hypercube.

Let w(v, n) := T_{P_{11}v+P_{12}n+P_{13}u}, then it has been noted in [PKT^{+}94], that if the first element of w(v, n) ≠ 0, (4.24) is equivalent to

σ̄(T_{v}(T_{w(v,n)}^{T}T_{w(v,n)})^{− \frac{1}{2}} ) ≤ γ.  (4.25)

If this forms a convex constraint on the set of solutions to equation (4.24), then it is equivalent to the LMI

[ \begin{matrix}γ^{2}T_{w(v,n)}^{T}T_{w(v,n)} & T_{v}^{T} \\ T_{v} & I \end{matrix} ] ≥ 0.

Hence the MVDP is equivalent to finding a point that satisfies an LMI feasibility problem.

Both a cutting plane algorithm, and the MATLAB toolbox LMI-lab [GN93], have been used to solve MVOP’s for NCF model-sets with noise. The cutting plane algorithm works well for simulated data, with few iterations required to converge within 5 percent of γ_{min}(T I). It is however much slower than the cutting plane algorithm used to calculate γ_{min}(T V), as an eigenvalue decomposition is necessary at each iteration.

Using LMI-lab to solve MVOP’s for NCF model-sets with noise has also been implemented. This software uses interior point methods which are more efficient than cutting plane algorithms. However, the software itself requires the basis of the LMI to be stored in memory. Hence the memory required grows as l^{3}, where l is the length of the data. This means that length of data is restricted, on a 500Mb machine, to approximately 250 variables, which if often insufficient.

Gap Metric and ν-Gap Metric

In this chapter we show that left and right NCF model-sets are not equivalent for validation and show the equivalence between validating NCF model-sets and certain balls in the gap metric. We also derive necessary and sufficient conditions for validating model-sets defined as balls in the ν-gap metric. Validation in the ν-gap metric is important as the ν-gap metric provides an elegant characterization of the largest set of models that can a priori be guaranteed to be stabilized by a certain set of H_{∞} controllers.

The ν-gap metric is closely related to the gap metric, which was introduced into the control literature in [ZES80], and results in [GS90] provide an interesting interpretation of NCF model-set validation. In fact validation of NCF model-sets can be interpreted as model validation in the gap metric, where validation in the gap metric means answering:

Given a nominal model, a set of input-output data and a positive real number β, does there a exist a model in a gap ball, centred on the nominal model with radius β, that interpolates the data?

Having obtained a connection between NCF model-sets and the gap metric, it is natural to consider model validation in the ν-gap metric. Results for validating model-sets that are defined as balls in the ν-gap metric provide the best possible validation results for a certain class of H_{∞} controllers.

In addition to providing an interpretation of model validation in the gap metric, NCF model-sets are important in their own right. They are well suited to robust controller design, [MG90] and are the model-sets used to design controllers for the flexible beam and Harrier, studied in chapters 6 and 7.

5.1 Normalized coprime factor models

In this section we consider the validation results for different NCF model-sets, with the conclusion being that N CF(Ñ, M̃, Δ, γ, δ), defined in sub-section 2.3.2, is the easiest to validate. Necessary and sufficient conditions for validating this model-set are proved in [PKT^{+}92, PKT^{+}94]. However, many variations of this model-set are possible, such as taking a right coprime factorization, or assuming the noise enters at different points.

5.1. NORMALIZED COPRIME FACTOR MODELS

Left and right NCF model-sets, shown in figure 5.1, are equivalent for robust controller design [MG90, p.58], so it is surprising that they are not equivalent for model validation. They are equivalent in the sense that there exists a controller to stabilize every element of N CF(Ñ, M̃, Δ, γ), if and only if there exists a controller to stabilize every element of N CF(N, M, Δ, γ). They are not equivalent in the sense that validating left NCF model-sets (without noise) can always be accomplished by testing the positive definiteness of a matrix formed from the data, whereas validating right NCF model-sets (without noise) may require a search over a non-convex set. Also, there exists a model P, such that if (Ñ, M̃) is a left NCF of P and (N, M) is a right NCF of P, there exists P_{1} ∈ N CF(N, M, Δ, γ) such that y = π_{l}P_{1}u, but there does not exist P_{1} ∈ N CF(Ñ, M̃, Δ, γ) such that y = π_{l}P_{1}u.

u  y u  u w_{1}  v  w_{2}  v_{1}  w  v_{2}^{+} ^{+}^{+}  ^{+}^{+}^{+}^{−} (a)  (b)^{−}  N ΔÑ  ΔM̃ M̃^{−1}  Δ_{M}  Δ_{N} M^{−1} Ñ

Figure 5.1: Left (a) and right (b) normalized coprime factor models

5.1.1 Left and right factorizations

Necessary and sufficient conditions for validating left NCF model-sets are proved in [PKT^{+}92], but right NCF model-sets are not mentioned. The techniques described in [PKT^{+}92] can be used to derive necessary and sufficient conditions for validating right NCF model-sets, but these conditions may give non-convex feasibility problems.

Before stating the result for N CF(Ñ, M̃, Δ, γ) we define a permutation matrix U, to allow a concise statement of the result. Define matrices F_{i} ∈ ℝ^{n_{u}×ln_{u}} and G_{i} ∈ ℝ^{n_{y}×ln_{y}} by

F_{i} := [ 0_{n_{u}×(i−1)n_{u}} I_{n_{u}} 0_{n_{u}×(l−i)n_{u}} ] ,
G_{i} := [ 0_{n_{y}×(i−1)n_{y}} I_{n_{y}} 0_{n_{y}×(l−i)n_{y}} ] .

Then define the permutation matrix U, where

U := [ \begin{matrix}F_{1} & 0 \\ 0 & G_{1} \\ F_{2} & 0 \\ 0 & G_{2} \\ ⋮ & ⋮ \\ F_{l} & 0 \\ 0 & G_{l} \end{matrix} ] ∈ ℝ^{l(n_{u}+n_{y})×l(n_{u}+n_{y})}.

This matrix is used to permute the elements of a vector formed from the inputoutput data; given input-output data u ∈ π_{l}Sn_{u}_{+} and y ∈ π_{l}Sn_{y}_{+}, then

[ u_{0} ; y_{0} ; u_{1} ; y_{1} ; ⋮ ; u_{l−1} ; y_{l−1} ] = U [ u_{0} ; u_{1} ; ⋮ ; u_{l−1} ; y_{0} ; y_{1} ; ⋮ ; y_{l−1} ] .

Using this matrix, the validation result for left NCF model-sets is: Proposition 5.1 ([PKT^{+}92]) Given a model-set N CF(Ñ, M̃, Δ, γ), and a set of input-output data u ∈ π_{l}Sn_{u}_{+} and y ∈ π_{l}Sn_{y}_{+}, let w := U [ \begin{matrix}vec(u) \\ vec(y) \end{matrix} ] and v := TM̃vec(y) − TÑvec(u). Then

i) If Δ = T I_{∞}^{n_{y}×(n_{y}+n_{u})}, there exists P_{1} ∈ N CF(Ñ, M̃, Δ, γ) satisfying

y = π_{l}T_{P_{1}}u,

if and only if

T_{v}^{T} T_{v} ≤ γ^{2}T_{w}^{T} T_{w}.

ii) If Δ = T V_{∞}^{n_{y}×(n_{y}+n_{u})}, there exists P_{1} ∈ N CF(Ñ, M̃, Δ, γ) satisfying

y = π_{l}T_{P_{1}}u,

if and only if

∥π_{k}v∥_{2} ≤ γ∥π_{2k}w∥_{2},

for all k=1,2,…,l.

The situation for right NCF models is entirely different, which is illustrated by the following example^{1}: Example 5.1 Given P = ( \begin{matrix}2 \\ 1 \end{matrix} ) \frac{3z−1}{z−3}, l = 1, u = 1, y = [ \begin{matrix}3 \\ 6 \end{matrix} ], let (N, M) be a normalized right coprime factorization of P and (Ñ, M̃) be a normalized left coprime factorization of P. Then

1. ϵ_{max} = 0.9129.

2. There does not exist P_{1} ∈ N CF(Ñ, M̃, T I_{∞}^{2×3}, 0.56) satisfying y = π_{1}P_{1}u.

3. There does exist P_{1} ∈ N CF(N, M, T I_{∞}^{3×1}, 0.56) satisfying y = π_{1}P_{1}u.

1The model P in this example is obtained by a bilinear transform of an example system in [GS90].

4. If γ > 0.923, then the set of [ \begin{matrix}v_{1} \\ v_{2} \end{matrix} ] satisfying

y = [ −NM^{−1} I] [ \begin{matrix}v_{1} \\ v_{2} \end{matrix} ] + NM^{−1}u,
∥[ \begin{matrix}v_{1} \\ v_{2} \end{matrix} ]∥_{2}^{2} ≤ γ^{2}∥M^{−1}u − M^{−1}v_{1}∥_{2}^{2},

is not convex.

We can calculate that ϵ_{max} = 0.9129 using MATLAB and the MU-TOOLS toolbox [BDG^{+}]. To show (2) we will calculate the smallest γ such that there exists P_{1} ∈ N CF(Ñ, M̃, T I_{∞}^{3×1}, γ) satisfying y = π_{1}P_{1}u. Using theorem 3.2, the necessary and sufficient conditions for the existence of a suitable P_{1} are that

∥[ \begin{matrix}0.9839 \\ 3.9071 \end{matrix} ]∥_{2} ≤ γ ∥[ \begin{matrix}1 \\ 3 \\ 6 \end{matrix} ]∥_{2} .

Hence the smallest value of γ is 0.5941.

To show (3) and (4) we calculate the smallest γ such that there exists P_{1} ∈ N CF(N, M, T I_{∞}^{3×1}, γ) satisfying y = π_{1}P_{1}u. Necessary and sufficient conditions

for there to exist a suitable P_{1} are that there exists [ \begin{matrix}v_{1} \\ v_{2} \end{matrix} ] satisfying

y = v_{2} + NM^{−1}u − NM^{−1}v_{1},  (5.1)
∥[ \begin{matrix}v_{1} \\ v_{2} \end{matrix} ]∥_{2}^{2} ≤ γ^{2}∥M^{−1}u − M^{−1}v_{1}∥_{2}^{2}.  (5.2)

Solving equation (5.1) for v_{2} and substituting into inequality (5.2), the necessary and sufficient conditions become that there exists v_{1} satisfying

(46 − 54γ^{2})v_{1}^{2} + (108γ^{2} − 18)v_{1} + 18 − 54γ^{2} ≤ 0.  (5.3)

From this equation we can calculate the smallest γ, such that a suitable v_{1} exists, to be 0.5484. It is also clear that the set of v_{1} satisfying equation (5.3) is not convex for γ > 0.923.

The example illustrates two important points. Firstly, left and right NCF modelsets are not equivalent for validation, as they contain different elements. We will see in the next section that this is because of the difference between gap balls and T-gap balls. Secondly, validating right NCF model-sets may involve a search over a nonconvex set. Searching over non-convex sets causes computational problems since there may be non-trivial local minima, so we cannot guarantee global optimality. Consequently, we will only consider left NCF model-sets for the remainder of the thesis, which is acceptable since they are equivalent to right NCF model-sets for controller design.

5.1.2 Including noise

In any practical problem data is corrupted by noise, which means the set of possible v’s, in figure 5.2, may not be convex. It is stated in [PKT^{+}94] that if noise is considered to enter at position 3 of figure 5.2, then validating the model-set N CF(Ñ, M̃, Δ, γ, δ), with Δ = T I_{∞}^{p×m} or T V_{∞}^{p×m}, is equivalent to solving a convex feasibility problem. However, noise may also be considered at positions 1,2,4 and 5 of figure 5.2, with positions 1 and 5 the natural places to consider sensor noise at the inputs and outputs. The convexity of corresponding validation problems can be analysed using the results in section 4.3.

Ñ  M̃^{−1}^{+}  y u 1  2  3  4  5^{+}  ^{+}  ^{+}  ^{+}^{+}^{+} ^{+}  ^{+}  ^{+}  ^{+}^{+}  ^{−}^{+} v w_{1}  w_{2} ΔÑ  ΔM̃

Figure 5.2: Left NCF model-set with possible places for noise.

If the uncertainty is time-varying, ie the model-set is N CF(Ñ, M̃, T V_{,}^{p×m}γ, δ), we can use theorem 4.3 to determine sufficient conditions for convexity. The results are summarized in the following table:

Position of noise  Conditions for convexity
1  (π_{k}N)^{T} π_{k}N − γ^{2}I_{k} ≥ 0 ∀ k = 1, 2, …, l
2  Always convex
3  Always convex
4  Always convex
5  (π_{k}M)^{T} π_{k}M − γ^{2}I_{k} ≥ 0 ∀ k = 1, 2, …, l

The table shows that for the model validation problem to be convex, noise can only be allowed to enter at positions 2,3 or 4. However, to account for sensor noise, it is more natural to consider noise at positions 1 and 5. To try and balance these conflicting requirements we will make the following assumption:

Assumption 5.1 When considering a left NCF model-set for validation, the noise will be considered to enter at position 3 in figure 5.2.

This is an approximation but is not unreasonable since noise at position 3 can be considered as a sum of the noise at the inputs and outputs, after being filtered by Ñ and M̃ respectively.

Having discussed both left and right NCF model-sets, and the possible places noise can enter, the best NCF model-set for validation is a left NCF model-set with noise entering at position 3 of figure 5.2. This is the model-set used to model the Harrier in chapter 7 and, after a slight modification, the flexible beam in chapter 6. The set of models was defined in section 2.3.2 to be N CF(Ñ, M̃, Δ, γ, δ). It can be validated using the following result:

Theorem 5.1 ([PKT^{+}92]) Given a model-set N CF(Ñ, M̃, Δ, γ, δ) and a set

of input-output data u ∈ π_{l}Sn_{u}_{+} and y ∈ π_{l}Sn_{y}_{+}, let w = U [ \begin{matrix}vec(u) \\ vec(y) \end{matrix} ] and b :=

TM̃vec(y) − TÑvec(u). Then

i) If Δ = T I_{∞}^{n_{y}×(n_{y}+n_{u})}, there exists P_{1} ∈ N CF(Ñ, M̃, Δ, γ, δ) satisfying

y = π_{l}T_{P_{1}} [ \begin{matrix}n \\ u \end{matrix} ] ,

if and only if there exists n ∈ ℝ^{ln_{y}}, with ∥n∥_{2} ≤ δ, satisfying

T_{b−n}^{T}T_{b−n} ≤ γ^{2}T_{w}^{T} T_{w}.

ii) If Δ = T V_{∞}^{n_{y}×(n_{y}+n_{u})}, there exists P_{1} ∈ N CF(Ñ, M̃, Δ, γ, δ) satisfying

y = π_{l}T_{P_{1}} [ \begin{matrix}n \\ u \end{matrix} ] ,

if and only if there exists n ∈ ℝ^{ln_{y}}, with ∥n∥_{2} ≤ δ, satisfying

∥π_{k}(b − n)∥_{2} ≤ γ∥π_{2k}w∥_{2},

for all k=1,2,…,l.

5.1.3 Computational issues

Since NCF model-sets have more structure than LFT model-sets, we can exploit the additional structure to yield more efficient algorithms for validation. For left NCF-model sets both γ_{min}(T I) and γ_{min}(T V) can be calculated without iteratively solving an MVDP, and if the infinity-norm of the noise is bounded, γ_{min}(T V) has an explicit solution. We have used these observations to analyse larger sets of data for the flexible beam and Harrier than would be possible for more general model-sets.

The following proposition shows that γ_{min}(T I) and γ_{min}(T V) can be calculated without an iterative search over γ.

Proposition 5.2 Given a model-set N CF(Ñ, M̃, Δ, γ, δ) and a set of input-output data u ∈ π_{l}Sn_{u}_{+}, y ∈ π_{l}Sn_{y}_{+} then γ_{min}(T I) and γ_{min}(T V) can be calculated by finding the minimum of linear function subject to LMI constraints.

Proof: Firstly let w = U [ \begin{matrix}vec(u) \\ vec(y) \end{matrix} ] and b := TM̃vec(y) − TÑvec(u). To prove

γ_{min}(T I) can be calculated as described, by theorem 5.1 the model-set is validated if and only if there exists n ∈ ℝ^{ln_{y}}, with ∥n∥_{2} ≤ δ, satisfying

T_{b−n}^{T}T_{b−n} ≤ γ^{2}T_{w}^{T} T_{w}.  (5.4)

Using Schur complements (5.4) can be written as an LMI in [ \begin{matrix}γ^{2} \\ n \end{matrix} ];

[ \begin{matrix}γ^{2}T_{w}^{T} T_{w} & T_{b}^{T} − T_{n}^{T} \\ T_{b} − T_{n} & I \end{matrix} ] ≥ 0.

5.2. THE GAP METRIC

Hence the smallest γ^{2} can be calculated as linear function of [ \begin{matrix}γ^{2} \\ n \end{matrix} ], subject to LMI

constraints in [ \begin{matrix}γ^{2} \\ n \end{matrix} ].

Proving γ_{min}(T V) can be calculated without an iterative search is similar. By theorem 5.1 the model-set is validated if and only if there exists n ∈ ℝ^{ln_{y}}, with ∥n∥_{2} ≤ δ, satisfying

∥π_{k}(b − n)∥_{2} ≤ γ∥π_{2k}w∥_{2}, ∀k = 1, 2, …, l.  (5.5)

Using Schur complements (5.5) can be written as an LMI in [ \begin{matrix}γ^{2} \\ n \end{matrix} ];

[ \begin{matrix}γ^{2}(π_{2k}w)^{T} (π_{2k}w) & (π_{k}b)^{T} − (π_{k}n)^{T} \\ π_{k}b − π_{k}n & I \end{matrix} ] ≥ 0 ∀k = 1, 2, …, l.

Hence the result is proved.

Minimizing a linear function subject to LMI constraints is called Positive Definite Programming [BVG94]. It is a generalization of linear programming and interior point algorithms have been developed that solve large problems efficiently [BVG94]. There are also routines in LMI-lab that solve this type of problem.

Another useful observation is that γ_{min}(T V) can be calculated explicitly if a bound is known on the infinity-norm of n, which can be seen from (5.5). The right hand side of this inequality is independent of n, so γ_{min}(T V) can be calculated by minimizing ∥π_{k}(b − n)∥_{2}, subject to ∥n∥_{∞} < δ. Hence the i-th component of n can be chosen as

n_{i} = \{ \begin{matrix}b_{i}, & if |b_{i}| < δ, \\ δ, & if b_{i} ≥ δ, \\ −δ, & if b_{i} ≤ −δ. \end{matrix}

We will use this fact in the calculation of γ_{min}(T V) for the flight-test data in chapter 7.

5.2 The gap metric

In this section we will show that validating NCF model-sets is equivalent to validating balls defined in the gap metric. The gap metric, as described in [GS90], is defined for continuous time systems, whereas the results of this thesis are for discrete time systems. Hence the relevant results from [GS90] need to be restated in discrete time. This can be achieved with the bilinear transform that replaces s by \frac{1+λ}{1−λ}, which maps the closed right half plane into the closed unit disc. For more details see [Par88].

Before defining the gap between two systems it is necessary to define the graph of an operator. Let P ∈ P^{p×m} then M_{P}, the associated multiplication operator on H_{2}, is defined as

M_{P} : H_{2} → H_{2}
u ↦ P u.

The graph of the operator M_{P}, denoted by G_{P}, is defined as

G_{P} := [ \begin{matrix}M \\ N \end{matrix} ] H_{2},

where (M, N) is a right normalized coprime factorization of P. This is guaranteed to exist, as in the continuous time case, by results in [Vid87], and is a closed subspace.

The gap between two systems is defined as the distance between their respective graph spaces. If the gap between two subspaces is denoted by δ(K_{1}, K_{2}), then

δ(K_{1}, K_{2}) := ∥Π_{K_{1}} − Π_{K_{2}}∥,

where Π_{K} is the orthogonal projection operator onto K. Standard operator theoretic techniques can be used to show that this is the maximum of two directed gaps, δ⃗(K_{1}, K_{2}) and δ⃗(K_{2}, K_{1}), where

δ⃗(K_{1}, K_{2}) := ∥(I − Π_{K_{2}})Π_{K_{1}}∥.

This can be calculated for two given systems using a result in [Geo88]:

Proposition 5.3 ([Geo88]) : For i = 1, 2, let the system P_{i} have a transfer function with NCF P_{i}(z) = N_{i}(z)M_{i}(z)^{−1}. Then

δ⃗(P_{1}, P_{2}) = inf_{Q∈H_{∞}} ∥( \begin{matrix}M_{1} \\ M_{1} \end{matrix} ) − ( \begin{matrix}M_{2} \\ N_{2} \end{matrix} ) Q∥_{∞} .

5.2.1 Implications

The connection between model validation and the gap metric occurs when balls in the gap metric are considered. Define the directed gap ball, B⃗(P, β), and the gap ball, B(P, β), as

B⃗(P, β) := {P_{1} ∈ P_{p×m} : δ⃗(P, P_{1}) < β} ,
B(P, β) := {P_{1} ∈ P_{p×m} : δ(P, P_{1}) < β} .

Then the important result in [GS90], connecting gap balls and NCF model-sets, is:

Lemma 5.1 ([GS90]) Let P have a NCF P = NM^{−1}. Then for all 0 < β ≤ 1,

B⃗(P, β) = {P_{1} ∈ P_{p×m} : P_{1} = (N + Δ_{N})(M + Δ_{M})^{−1} where

Δ_{N}, Δ_{M} ∈ RH_{∞} and ∥( \begin{matrix}Δ_{M} \\ Δ_{N} \end{matrix} )∥_{∞} < β\} .

Hence asking if a right NCF model-set is consistent with a set of data, is the same as asking if a directed gap ball is consistent with the same set of data.

Similar results exist for left NCF model-sets by defining the transpose gap (Tgap);

δ⃗_{T} (P_{1}, P_{2}) := δ⃗(P_{1}^{T} , P_{2}^{T}),
δ_{T} (P_{1}, P_{2}) := δ(P_{1}^{T} , P_{2}^{T}).

The corresponding result is:

Lemma 5.2 ([GS90]) Let P have a NCF P = M̃^{−1}Ñ. Then for all 0 < β ≤ 1,

B⃗_{T} (P, β) = {P_{1} ∈ P_{p×m} : P_{1} = (M̃ + ΔM̃)^{−1}(Ñ + ΔÑ) where
ΔÑ, ΔM̃ ∈ RH_{∞} and ∥( ΔM̃ ΔÑ )∥_{∞} < β} .

Hence asking if a left NCF model-set is consistent with a set of data, is the same as asking if a directed T-gap ball is consistent with the same set of data.

It was demonstrated in section 5.1 that left and right NCF model-sets are not equivalent for model validation, so it is natural to ask about the relationship between gap balls and T-gap balls. The relationship is described by a result in [GS90]. Before stating the result it is necessary to define two quantities, λ(P) and b_{opt}(P). Firstly, define λ(P) and λ(P^{T}) by

λ(P) := inf_{|z|<1} σ_{min} ( \begin{matrix}M(z) \\ N(z) \end{matrix} ) ,
λ(P^{T}) := inf_{|z|<1} σ_{min} ( M̃(z) Ñ(z) ) .

Secondly, let b_{opt}(P) denote the supremum over all γ, such that there exists a controller C which stabilizes every element of the model-set N CF(Ñ, M̃, T I_{,}^{p×m}γ). Then a lemma relating these two quantities is

Lemma 5.3 ([GS90]) b_{opt}(P) ≤ λ(P).

The result showing the relationship between gap balls and directed gap balls is

Lemma 5.4 ([GS90]) If β < λ(P), then B⃗(P, β) = B(P, β).

The same result is also true in the T-gap, ie B⃗_{T} (P, β) = B_{T} (P, β) provided β < λ(P^{T}) [GS90]. Hence if the uncertainty in the NCF model-set is small enough to be stabilized, the directed gap ball and gap ball are equal.

The relationship between gap balls and T-gap balls is not as clear. The only result about this relationship in [GS90] is that for SISO systems B(P, β) = B_{T} (P, β), which is because for SISO systems the transfer functions commute. However, in the MIMO case gap and T-gap balls may be different, which is demonstrated by a counterexample in [GS90].

We can use these results to prove the main results of this section; that validation of left NCF models provides a lower bound for the directed T-gap between the model and system. A similar result also holds for right NCF models and directed gap balls.

Theorem 5.2 Given a “real system” P_{0}, a nominal model P_{1} and a set of inputoutput data from the “real system” u ∈ π_{l}Sn_{u}_{+} and y ∈ π_{l}Sn_{y}_{+} (y = π_{l}P_{0}u), let (M̃, Ñ) be a left NCF of P_{1} and γ_{min} be the smallest γ such that there exists P_{2} ∈ N CF(Ñ, M̃, T I_{∞}^{n_{y}×(n_{y}+n_{u})}, γ) satisfying y = π_{l}T_{P_{2}}u. Then γ_{min} ≤ δ⃗_{T} (P_{1}, P_{0}).

Proof: The result is a consequence of lemma 5.2.

Corollary 5.1 If γ_{min} < λ(P_{1}) then γ_{min} ≤ δ_{T} (P_{1}, P_{0}).

Corollary 5.2 If n_{u} = n_{y} = 1 and γ_{min} < λ(P_{1}) then γ_{min} ≤ δ(P_{1}, P_{0}).

The difference between γ_{min} and the gap between the model and system depends on the model and the length, and quality, of data used for validation.

Remark 5.1 For SISO systems, if the input to the uncertainty w_{0} is non-zero, γ_{min} is guaranteed to converge to the directed T-gap between the systems as the length of data tends to infinity. This is clear from the following equation, which relates the inputs and outputs of the uncertainty, w and v respectively, via the impulse response of the uncertainty (h_{0}, h_{1}, …, h_{l−1});

[ h_{l−1} h_{l−2} ⋯ h_{0} ] [ \begin{matrix}w_{0} & 0 & … & 0 \\ w_{1} & w_{0} & ... & ⋮ \\ ⋮ & ⋮ & ... & 0 \\ w_{l−1} & w_{l−2} & … & w_{0} \end{matrix} ] = [ v_{l−1} v_{l−2} ⋯ v_{0} ] .

Clearly, the impulse response of the uncertainty is uniquely determined by the v and w. Hence, as the length of this data tends to infinity, the infinite impulse response will be uniquely determined for all time, so γ_{min} will converge to the directed Tgap between the model and system. Note that this argument falls down for systems with multiple inputs as the uncertainty is no longer uniquely determined, although convergence may still occur.

The following simple example shows how γ_{min} varies with the length of data for a simple system with simulated data.

Example 5.2 Let P_{0} be the “real system” and P_{1} the model, where the transfer functions are P_{0} = \frac{7z−9}{3z−1} and P_{1} = \frac{19z−21}{3z−1}. Input-output data can be obtained using a random input signal u, with y being the corresponding output from P_{0}. Using this data and P_{1}, γ_{min}(T I) can be calculated and the results, for different lengths of data, are shown in figure 5.3.

γ_{min}(T I)  length of data  δ(P_{0}, P_{1})

Figure 5.3: Variation of γ_{min}(T I) with the length of data

5.2.2 Limitations

Results on the gap metric provide an interesting interpretation of model validation but, since gap balls and T-gap balls are not equal, problems may occur.

Let P ∈ P^{p×m}, C ∈ P^{m×p}, then say that [P, C] is stable if [ \begin{matrix}I & −C \\ −P & I \end{matrix} ]^{−1} ∈

RH_{∞}. This is the same as the standard feedback configuration of P and C being internally stable. Then a result of [GS90] is:

5.3. THE ν-GAP METRIC

Theorem 5.3 ([GS90]) The following statements are equivalent:

a) [P_{1}, C] is stable for all P_{1} such that δ(P, P_{1}) < β.

b) [P_{1}, C] is stable for all P_{1} such that δ_{T} (P, P_{1}) < β.

This theorem says that a controller, C, that stabilizes a gap ball around P will also stabilize a T-gap ball around P, and vice versa. However gap balls and T-gap balls are not necessarily the same. Hence a situation may arise, like in example 5.1, where validation of the left NCF model yields γ_{min}(T I) bigger than some number β. If validation of the right NCF model gives γ_{min}(T I) less that β then theorem 5.3 says that all left coprime factor uncertainty of size β are still stabilizable even though γ_{min}(T I) was greater than β. This is not acceptable if model validation is for robust controller design.

5.3 The ν-gap metric

In the previous sub-section we indicated a problem that may arise when validating model-sets in the gap metric. The problem arose because a gap ball does not contain every system that can be guaranteed to be stabilized by a certain set of H_{∞} controllers. The largest set of systems that can be guaranteed to be stabilized can be described neatly in the ν-gap metric [Vin93].

The ν-gap (in continuous time) is defined as follows:

Definition 5.1 ([Vin92]) Let P_{i} ∈ P^{p×m} have normalized right and left coprime factorizations denoted by (N_{i}, M_{i}) and (Ñ_{i}, M̃_{i}) respectively. Define a function δ_{ν}(·, ·) : P^{p×m} × P^{p×m} → ℝ as follows.

δ_{ν}(P_{1}, P_{2}) := \{ \begin{matrix}∥ − M̃_{2}N_{1} + Ñ_{2}M_{1}∥_{∞} & if det (N_{2}^{∗}N_{1} + M_{2}^{∗}M_{1})(jω) ≠ 0 ∀ω \\ & and wno det (N_{2}^{∗}N_{1} + M_{2}^{∗}M_{1}) = 0, \\ 1 & otherwise, \end{matrix}

where wno(g) denotes the winding number about the origin of g(s), as s follows the standard Nyquist D-contour.

To state the key result on the ν-gap metric, that shows that it describes the largest set of systems to be stabilized by a set of H_{∞} controllers, let b_{P,C} be defined as

b_{P,C} := \{ \begin{matrix}∥[ \begin{matrix}P \\ I \end{matrix} ] (I − CP)^{−1} [ −C I] ∥_{∞}^{−1} , & if [P, C] is stable \\ 0, & otherwise \end{matrix}

Then the result is:

Theorem 5.4 ([Vin92])

i) Given a nominal system P_{1} ∈ P^{p×q}, a compensator C ∈ P^{q×p} and a number β, then: [P_{2}, C] is stable for all systems P_{2} satisfying δ_{ν}(P_{1}, P_{2}) ≤ β if, and only if, b_{P_{1},C} > β. ii) Given a nominal system P_{1}, a perturbed system P_{2} ∈ P^{p×q} and a number β <

b_{opt}(P_{1}) then:

[P_{2}, C] is stable for all compensators C satisfying b_{P_{1},C} > β if, and only if, δ_{ν}(P_{1}, P_{2}) ≤ β.

Note that condition i) of this theorem also holds for the gap metric.

The importance of this result is succinctly described in [Vin92]:

So, as for the gap metric, any plant at a distance less then β from the nominal will be stabilized by any compensator stabilizing the nominal with a stability margin of β. Furthermore, any plant at a distance greater then β from the nominal will be destabilized by some compensator that stabilizes the nominal with a stability margin of at least β. The only metric that has both these properties is δ_{ν}(·, ·). So, a ball in this metric is the largest set of linear, time-invariant plants that can be guaranteed to be stabilized a priori by a compensator solving the following H_{∞} problem: find a C satisfying b_{P,C} > β.

Hence the controller C, satisfying b_{P,C} > β, is closely related to a ball in the νgap metric. In fact we will use a parameterization of all C satisfying b_{P,C} > β to parameterize a ball in the ν-gap.

A general H_{∞} problem has been solved, under mild assumptions, in [DGKF89, GD88]. Given a plant G, partitioned conformally with the block structure as follows;

G = [ \begin{matrix}G_{11} & G_{12} \\ G_{21} & G_{22} \end{matrix} ] ,

define a set of H_{∞} controllers as:

C(G, β) := {C : [G_{22}, C] is stable , ∥F_{l}(G, C)∥_{∞} < 1/β} .

Then from [GD88], if there exists a C such that F_{l}(G, C) ∈ H_{∞}, and G_{12}(jω) has full column rank, and G_{21}(jω) full row rank, for all ω, then

C(G, β) = {C : C = F_{l}(J, Q), Q ∈ RH_{∞}, ∥Q∥_{∞} < β^{−1}} .

The H_{∞} problem of finding a C satisfying b_{P,C} > β can easily be solved using this result for a suitable choice of G. An obvious choice for G is

[ \begin{matrix}0 & 0 & I \\ P & I & P \\ P & I & P \end{matrix} ] ,

but other G are possible, eg one corresponding to a right NCF model [GS90].

5.3.1 Implications

Using the observation that “… a ball in this (ν-gap) metric is the largest set of linear, time-invariant plants that can be guaranteed to be stabilized a priori by a compensator solving the following H_{∞} problem: find a C satisfying b_{P,C} > β,” and the parameterization of all sub-optimal H_{∞} controllers from [GD88], a parameterization of a ν-gap ball can be obtained. This provides a method of doing model validation in the ν-gap.

Firstly define a ν-gap ball,

B_{ν}(P, β) := {P_{1} ∈ P^{p×m} : δ(P, P_{1}) ≤ β} .

Then a parameterization of a ν-gap ball, recently obtained in [Vina], is:

Theorem 5.5 ([Vina]) Suppose that b_{opt}(P) > β and let J be the coefficient matrix in the parameterization of all sub-optimal H_{∞}controllers in [GD88] corresponding to the generalized plant, G, where G is partitioned as

G = [ \begin{matrix}[ 0 −M^{−1} ] & M^{−1} \\ [ I −P ] & P \end{matrix} ] ,

and (M, N) is a right normalized coprime factorization of P. Then

B_{ν}(P, β) = {F_{l}(J^{−1}, Δ) : Δ ∈ RH_{∞}, ∥Δ∥_{∞} < β} .

Proof: Suppose b_{opt}(P) > β, then the standard assumptions for G in [GD88, DGKF89] are satisfied and the set of all robustly stabilizing controllers can be parameterized as

{F_{l}(J, Q) : Q ∈ RH_{∞}, ∥Q∥_{∞} < β^{−1}} .

Furthermore, if J is partitioned as

J = [ \begin{matrix}J_{11} & J_{12} \\ J_{21} & J_{22} \end{matrix} ] ,

then from [GGLD90] it can be shown that J_{12} and J_{21} are invertible in RH_{∞}, and that the free parameter Q may be written as

Q = F_{u}(J^{−1}, C).

Hence, by application of the Small Gain Theorem, the set of all systems stabilized by all the controllers C satisfying b_{P,C} < β, can be parameterized as

{F_{l}(J^{−1}, Δ) : Δ ∈ RH_{∞}, ∥Δ∥_{∞} ≤ β} .

Theorem 5.4 states that this set of systems is the same as a ν-gap ball, of radius β, centred on P, so the result is proved.

This theorem can be used to do model validation in the ν-gap, which is the main result of this section:

Theorem 5.6 Given a model-set B_{ν}(P, γ) and a set of input-output data u ∈ π_{l}Sn_{u}_{+}, y ∈ π_{l}Sn_{y}_{+}, suppose b_{opt}(P) > γ and let J be the system defined in theorem 5.5 and K its inverse partitioned as

K := J^{−1} = [ \begin{matrix}K_{11} & K_{12} \\ K_{21} & K_{22} \end{matrix} ] .

Then there exists P_{1} ∈ B_{ν}(P, γ) satisfying

y = π_{l}T_{P_{1}}u,

if and only if

T_{v}^{T} T_{v} ≤ γ^{2}T_{w}^{T} T_{w},

where v = T_{K_{12}^{−1}}(y − T_{K_{11}}u) and w = T_{K_{21}}u + T_{K_{22}}v.

Proof: The result follows by a straightforward application of theorems 3.2 and 5.5.

Note that this result is only valid if the size of the uncertainty is less than b_{opt}(P). This is a limitation of the result, as there may exist a controller that stabilizes a model which is outside the ν-gap ball. However, if a model lies outside the ν-gap ball there exists an H_{∞} sub-optimal controller that will destabilize it. So, from an H_{∞} viewpoint, the ν-gap ball is the largest set of systems that would ever be required to be validated.

Theorem 5.5 is extremely powerful and provides an elegant parameterization of the largest set of linear, time-invariant systems that can be stabilized by a certain set H_{∞} controllers. However, a stronger result was also proved in [Vina], that the same set of controllers will also stabilize a set of possibly nonlinear/time-varying operators whose norm is bounded by β. Let B denote the space of possibly nonlinear/timevarying operators on l_{2+} with bounded incremental norm

∥B∥_{Δ} := sup_{u_{1}≠u_{2}} \frac{∥Bu_{1} − Bu_{2}∥_{2}}{∥u_{1} − u_{2}∥_{2}} .

Then from [Vina], with the same assumptions as theorem 5.5, any element of C(G, β) will also stabilize the set,

Bν^{NLT V} (P, β) := {F_{l}(J^{−1}, Δ) : Δ ∈ B, ∥Δ∥_{Δ} ≤ β} .

This model-set can be validated using interpolation results for time-varying operators:

Theorem 5.7 Given a model-set Bν^{NLT V} (P, γ) and a set of input-output data u ∈ π_{l}Sn_{u}_{+}, y ∈ π_{l}Sn_{y}_{+}, suppose b_{opt}(P) > γ and let J be the system defined in theorem 5.5 and K its inverse partitioned as

K := J^{−1} = [ \begin{matrix}K_{11} & K_{12} \\ K_{21} & K_{22} \end{matrix} ] .

Then there exists P_{1} ∈ Bν^{NLT V} (P, γ) satisfying

y = π_{l}T_{P_{1}}u,

if and only if

∥π_{k}v∥_{2} ≤ γ∥π_{k}w∥_{2} ∀k = 0, 1, …, l − 1,

where v = T_{K_{12}^{−1}}(y − T_{K_{11}}u) and w = T_{K_{21}}u + T_{K_{22}}v.

Proof: The result follows by a straightforward application of theorems 3.3 and the result from 5.5.

As for validation in the gap, the smallest γ that satisfies the necessary and sufficient conditions of theorem 5.6, is a lower bound on the gap between the model and real system. Denote this smallest γ by γ_{min}^{ν}(T I), then the result of repeating example 5.2 is shown in figure 5.4. This clearly shows γ_{min}^{ν}(T I) converging to the ν-gap between the nominal model and the model that produced the data. This is guaranteed for scalar systems for the same reasons as remark 5.1.

δ(P_{0}, P_{1}) δ_{ν}(P_{0}, P_{1}) γ_{min}(T I) γ_{min}^{ν}(T I)  length of data smallest γ
Figure 5.4: Variation of γ_{min}(T I) and γ_{min}^{ν}(T I) with the length of data

5.3.2 Including noise

Theorem 5.6 provides a method of doing model validation in the ν-gap metric, but it does not allow the data to be corrupted by noise. Allowing sensor noise on the inputs and outputs can make model validation equivalent to a non-convex feasibility problem. However, making similar assumptions to those in sub-section 5.1.2, the model validation problem is convex.

A sufficient condition for making the model validation problem convex is making the input to the uncertainty, w, independent of the noise. This can be accomplished, as in sub-section 5.1.2, by only allowing noise to enter the model at the same point as the uncertainty input, v, as in figure 5.5. This appears unnatural but is no more

y  u n v w  Δ K_{11} K_{21} K_{22} K_{12}

Figure 5.5: LFT parameterization of ν-gap ball with noise

unnatural than assumption 5.1. Assumption 5.1 considered the noise to be a sum of noise at the inputs and outputs, filtered by Ñ and M̃ respectively. The noise in figure 5.5 may also be considered as a sum of noise at the inputs and outputs, but filtered by K_{12}^{−1}K_{11} and K_{12}^{−1} respectively. This is equivalent to the case for left NCF model-sets, as if K is chosen so that figure 5.5 is a left NCF model-set, then K_{12}^{−1}K_{11} = Ñ and K_{12}^{−1} = M̃.

To state these ideas formally denote the set of systems in figure 5.5 by B̂_{ν}(P, β), where

B̂_{ν}(P, β) := \{F_{l} ([ \begin{matrix}K_{11} & K_{12} & K_{12} \\ K_{21} & K_{22} & K_{22} \end{matrix} ] , Δ) : ∥Δ∥_{∞} < β\} .

Then the result for validating model-sets of this form is:

Theorem 5.8 Given a model-set B̂_{ν}(P, γ) and a set of input-output data u ∈ π_{l}Sn_{u}_{+}, y ∈ π_{l}Sn_{y}_{+}, suppose b_{opt}(P) > γ and let J be the system defined in theorem 5.5 and K its inverse partitioned as

[ \begin{matrix}K_{11} & K_{12} \\ K_{21} & K_{22} \end{matrix} ] .

Then there exists P_{1} ∈ B̂_{ν}(P, γ) and a vector n ∈ ℝ^{ln_{v}} satisfying

y = π_{l}T_{P_{1}} [ \begin{matrix}u \\ n \end{matrix} ] ,
∥n∥_{2} ≤ δ,

if and only if there exists n ∈ ℝ^{ln_{v}} satisfying

[ \begin{matrix}γ^{2}T_{w}^{T} T_{w} & T_{b}^{T} − T_{n}^{T} \\ T_{b} − T_{n} & I \end{matrix} ] ≥ 0,
∥n∥_{2} ≤ δ,

where b := T_{K_{12}^{−1}}(y − T_{K_{11}}u) and w := T_{K_{21}}u + T_{K_{22}}b.

Proof: The proof follows easily from the structure of the elements in the set B̂_{ν}(P, γ) using Schur complements to expand the necessary and sufficient conditions.

This result clearly shows that the model validation for models of this kind results in a convex feasibility problem, as the necessary and sufficient conditions are LMI’s.

Flexible Beam

In this chapter, and the following one, we apply the theory developed in the previous chapters to model-sets and data from a flexible beam and Harrier VSTOL aircraft. A flexible beam and Harrier are very different, but both provide challenging control problems. We study the flexible beam first, since it has a less complicated model (having only one input and one output) than the Harrier.

The motivation for studying a flexible beam stems from interest in verifying theoretical models of violin bows [Woo93a, Woo93b]. To verify these models, experiments must be performed where the force between the string and bow is precisely controlled. Uncertainty about the true parameters of the system, and the precise dynamics, makes this is a challenging control problem, for which a robust controller is ideal. The design of robust controllers for this problem have been studied in [WWa, WWb], and this chapter is based on the apparatus, and models, described in [WWa, WWb].

Rigid Hub  Air table Electromagnet Flexible beam Motor housing  Permanent magnet Force transducer

Figure 6.1: Experimental setup of a flexible beam

6.1. MODEL

A diagram of the apparatus is shown in figure 6. It consists of a 0.67m long flexible beam, mounted on an aluminium hub, which is in turn attached to a brushless DC motor. The beam comes into contact with the force transducer through a low friction ball bearing, and the force transducer is free to move along the beam, since it is mounted on an airtable. The end of the beam is fixed using an electromagnet. The system is treated as having a single input; the current to the motor, and a single output; the force exerted on the force transducer at the point of contact between the beam and transducer. The controllers are implemented, and outputs measured, using a PC with dSPACE hardware and software [dSP]. This hardware samples the outputs and holds the inputs before being recorded.

6.1 Model

Hamilton’s principle was used in [WWb] to obtain a mathematical model of the beam. Several simplifying assumptions were made, such as the hub at the end of the beam was perfectly rigid, all displacements were small and that friction was negligible. Key parameters in this model, such as the hub inertia and modal damping factors, were then identified using data from the beam. Finally the infinite dimensional model was approximated by a finite dimensional model and written in state space form, with ten states, one input and one output, for control design.

The dynamics of the model vary significantly with the distance between the motor axis and the contact point between the beam and force transducer. Consequently the model was parameterized by the contact point. A Bode magnitude plot of the continuous time transfer function of the model, G(s), with the contact point at 0.35m, is shown in figure 6.2.

|G(jω)|  Frequency (rad/sec)

Figure 6.2: Bode magnitude plot of G(s)

6.1.1 Model-set

In [WWa], Wood describes how controllers were designed with the H_{∞}-loopshaping techniques described in [MG90], and extended in [Vin92]. These techniques require the model to include weighting functions, which increased the number of states in the model to fourteen. A Bode magnitude plot of the continuous time transfer function, with the contact point at 0.35m, including the weighting functions, G(s)W (s), is shown in figure 6.3.

To account for uncertainty about the true dynamics, the model-set was parameterized as perturbations in the normalized coprime factors of the nominal model. Hence the model-set chosen was N CF(Ñ, M̃, Δ, γ), with Δ = T I_{∞}^{p×m}

|G(jω)W (jω)|  Frequency (rad/sec)

Figure 6.3: Bode magnitude plot of G(s)W (s)

or Δ = T V_{∞}^{p×m}. Recall that an H_{∞}-controller, which is guaranteed to stabilize N CF(Ñ, M̃, T I_{∞}^{p×m}, γ), will also stabilize N CF(Ñ, M̃, T V_{∞}^{p×m}, γ). In fact it will also stabilize N CF(Ñ, M̃, Δ, γ), if Δ is a set of nonlinear operators, bounded in norm by γ, as a consequence of the Small Gain theorem [Zam81].

In N CF(Ñ, M̃, Δ, γ) the uncertainty is unstructured, which may be conservative if it is used to cover uncertainty in the parameters of the system. However knowledge of how parameter variations affect the model can be included by weighting the ball of uncertainty, as described in [Vin92]. This results in a left NCF model-set with an additional uncertainty weight, W_{u}, as shown in figure 6.4. Since W_{u} is invertible in

Ñ  ΔM̃ M̃^{−1}^{+}  ^{−} u  y w_{1}  ΔÑ  v n W_{u}  w_{2}^{+}^{+} ^{+}

Figure 6.4: Left NCF model-set with weighted ball of uncertainty and noise

RH_{∞}, the model validation problem is convex. Using the notation of section 5.1,

the equation for the uncertainty input, w, is unchanged (w = U [ \begin{matrix}vec(u) \\ vec(y) \end{matrix} ]) but the

uncertainty output, v, now becomes

v = T_{W_{u}^{−1}}(TM̃y − TÑu − n) .

As w is still independent of n, the model validation problem remains convex for any γ.

6.2. EXPERIMENT DESIGN

One feature not explicitly considered in Wood’s model was noise. If model validation is to be accomplished by solving a convex feasibility problem, the noise can only be allowed to enter at certain positions, as was noted in section 5.1. Consequently assumption 5.1 is made, and noise is only considered to enter at the same point as the uncertainty signal v. This is shown, together with the uncertainty weight in figure 6.4. This is an approximation but, as was noted in section 5.1, is not unreasonable.

In the previous chapter it was shown that testing the validity of NCF modelsets is equivalent to model validation in the gap metric, and that γ_{min}(T I) is a lower bound on the gap between the model and system. If γ_{min}(T I) is close to this gap, then a plot of how γ_{min}(T I) varies as the model changes with contact point, should be the similar to a plot of how the gap between different models varies with the contact point. A plot showing the variation of the gap, and ν-gap, between a weighted model and the weighted model at 0.35m contact point is shown in figure 6.5. This plot indicates that the beam may be difficult to robustly control,

ϵ_{max}  contact point (m) gap between models

Figure 6.5: Variation of gap and ν-gap (identical) between linear models

as a controller that stabilizes the model at 0.35m may destabilize a model with a contact point less than 2cm away. However, this is not the case as controllers have been designed that are robust over a much larger set of contact points [WWa]. The gap, and ν-gap metrics, require every element of their respective balls to be stabilized by every H_{∞} controller that satisfies a certain H_{∞} norm bound. However, for the beam, the H_{∞} central controller achieves poor robustness to variations in the contact point of the beam [WWa]. This is due to the transfer function zero appearing just before the first mode. Including W_{u} helps to overcome this problem, for further discussion see [Vinb].

6.2 Experiment design

The model used in [WWa] for controller design was in continuous time, but a discrete time model-set is required for validation. This was obtained by assuming the inputs and outputs were sampled and held (with a zero order hold) at 100Hz. This is justified as the hardware used to implement the controllers does exactly this.

To obtain data for model validation it was decided to input a swept sinusoid whose frequencies covered the first two modes of the beam. Due to the problems of analysing long lengths of data the sweeps were accomplished in 1 second. Hence the sinusoid swept from 0.1 to 30Hz in 1 second, with sampling at 100Hz. The actual input was

sin(2π(0.1t + 14.95t^{2})), 0 ≤ t ≤ 1.

6.2.1 Data length

To test how γ_{min}(T I) converges to the gap between the model and system for the chosen input, data was simulated using a linear model. The output of the discrete linear model at 0.35m (obtained by a sample and hold at 100Hz) was simulated with the input sin(2π(0.1t + 14.95t^{2})), 0 ≤ t ≤ 1. Portions of this data were then used to try and validate models at different contact points and the results are shown in figure 6.6. γ_{min}(T I) was calculated exactly as the data was not corrupted by noise and the results of section 5.1 could be applied. The plot shows good convergence

10 30 60 100 gap ϵ_{max}  contact point (m) γ_{min}(T I)

Figure 6.6: Variation of γ_{min}(T I) with data length and contact point

between γ_{min}(T I), and the gap between the model and “real system”, for the chosen input. This does not imply anything about the convergence between γ_{min}(T I) and the gap between the model and physical system, but gives some confidence that the input is suitable.

6.2.2 Computational issues

We calculated the results for γ_{min}(T I) and γ_{min}(T V) using the results in subsection 5.1.3. However, the noise was assumed to be bounded in the 2-norm so γ_{min}(T V) could not be calculated analytically. Instead a cutting plane algorithm was used, as described in section 4.4. This worked well on the 100 points of data analysed, converging in under 10 minutes on a SPARC-10 when coded in MATLAB.

We calculated γ_{min}(T I) without iteratively solving an MVDP for different γ, using the result in sub-section 5.1.3, by minimizing a linear function subject to LMI constraints. The MATLAB toolbox LMI-lab was used for this and the algorithm was far slower than the cutting plane algorithm used to calculate γ_{min}(T V). It typically took 3 hours on a SPARC-10, again using MATLAB, and consumed over 50Mb of RAM. Hence it was the speed, and memory requirements, of this algorithm that limited the number of data points to be 100.

6.3. EXPERIMENTAL DATA

6.3 Experimental data

The data, together with the simulated output from the discrete linear model, is shown in figures 6.7 and 6.8. Examining the plot of the actual output and the

voltage (V)  time (secs)

Figure 6.7: Input to power amplifier

voltage (V)  time (secs)

Figure 6.8: Output from force transducer (solid) and simulated output from discrete-time model (dash-dot)

simulated output reveals some marked differences. The model predicts that the first mode is excited after approximately 0.2s and the second mode after 0.7s. Around the first mode the simulated output is close to the actual output but this is not true for the second mode. This could be because the predicted position of the second mode is at a lower frequency than it really is.

This observation is borne out by doing another experiment, with the input a swept sinusoid sweeping over the same frequency range, but taking five seconds to do so. The plot of this data is shown in figure 6.9. As this input sweeps over frequency more slowly, there are more cycles at each frequency to excite the modes of the system. From this plot the predicted second mode appears approximately 0.25s before the actual second mode. This only corresponds to a difference in the frequencies of the modes of 1.5Hz. So although the actual and simulated outputs appear very different, the difference may be accounted for by a difference between the predicted and actual position of the second mode of 1.5Hz.

voltage (V)  time (secs)

Figure 6.9: Output from force transducer (solid), and simulated output from model (dotted), for 5 second sine sweep

It is misleading to compare the open loop responses of two systems when comparing their difference in the gap or ν-gap. It is shown in [Vin92] that two systems can have a small gap, or ν-gap, between them and yet have very different open loop responses. Also two systems that appear to have similar open loop responses may have a large gap between them. So although the actual and simulated outputs appear different, this does not imply that the validation results will be poor.

6.3.1 Measuring noise

To calculate γ_{min}(·) it is necessary to bound the norm of the noise. The norm bound was obtained by measuring the output of the force transducer, and power amplifier^{1}, with no input. The norm bound was then obtained by measuring the norm of these outputs after filtering through M̃ and Ñ respectively^{2}. Plots of the filtered signals are shown in figures 6.10 and 6.11.

We took the norm bound on the noise to be the sum of the 2-norms of the signals shown in figures 6.10 and 6.11, which for 500 data points was 0.2454. Hence we took the norm bound, for data of length l, to be

0.2454\sqrt{\frac{l}{500}}.

6.3.2 Invalidation

Using the noise bounds described, we calculated the smallest norm bound on the time-invariant uncertainty in the gap and ν-gap metrics, denoted by γ_{min}(T I) and γ_{min}^{ν}(T I) respectively, for different lengths of data and the model at 0.375m contact point. This was done because of the behaviour shown in figure 5.4, where γ_{min}(T I) tended towards the gap between the model and system as the length of

1The noise on the power amplifier was initially considered negligible. However the resulting invalidation of the model led to a careful measurement of this source of noise!

2It should be emphasized that the bound has been taken from a single set of data and that the noise is assumed to be independent of the input signal. Given more sets of data statistical techniques could be used to give greater confidence in the noise bound and to estimate the noise bound for the specific choice of input.

voltage (V)  time (secs)

Figure 6.10: Output noise after being filtered by M̃

voltage (V)  time (secs)

Figure 6.11: Input noise after being filtered by Ñ

data increased. Therefore to increase our confidence that γ_{min}(T I) is a good approximation to the gap between the model and system, and that γ_{min}^{ν}(T I) is a good approximation to the ν-gap between the model and system, we would like to see a convergence of γ_{min}(T I) as the length of data increases. This convergence appears to be happening in figure 6.12, but cannot be guaranteed. Note that ϵ_{max} is normalized to 1, as it is in all of the plots in this chapter. We repeated the calculation for time-varying uncertainty, and the results are shown in figure 6.13.

Figure 6.12 is interesting, showing that the model at 0.375m contact point is invalidated for time-invariant uncertainty in the gap metric, but validated for timeinvariant uncertainty in the ν-gap metric. Therefore, for this data, the difference between gap balls and ν-gap balls is crucial for the validity of the model-sets.

Figure 6.14 shows how γ_{min}(·) varies with the contact point of the beam in the gap metric^{3}. Note that γ_{min}(T I) is always greater than 1 whereas γ_{min}(T V) always less than 1. This poses a dilemma. If it is thought that the uncertainty is timevarying then the model cannot be invalidated for this set of data. However, if it is thought that the uncertainty is time-invariant then this set of data invalidates the model. There is no reason to consider uncertainty to be time-varying, but there are

3We did not obtain results for how γ_{min}(·) varies with the contact point in the ν-gap metric since validation in the ν-gap requires an iterative search over γ, which is prohibitively time consuming for a number of contact points.

γ_{min}^{ν}(T I) γ_{min}(T I) length of data smallest γ
Figure 6.12: Variation of γ_{min}(T I) and γ_{min}^{ν}(T I) with length of data

γ_{min}(T V) γ_{min}^{ν}(T V) length of data smallest γ
Figure 6.13: Variation of γ_{min}(T V) and γ_{min}^{ν}(T V) with length of data

certainly nonlinear effects that will not accurately be modelled by the linear model, such as back emf’s in the motor. It may therefore be appropriate to consider these effects as being modelled by linear time-varying uncertainty.

6.3.3 Modifying the model

Having invalidated the model-set for time-invariant uncertainty in the gap metric, the model was reviewed. The input and output from the power amplifier were measured and plots of these signals are shown in figure 6.15. This clearly shows that the power amplifier introduces a 20ms delay that is not included in the model. Using this observation the continuous time model of the beam was modified by premultiplying the model by a second order Pade approximation of the 20ms delay. The discrete model was then obtained using a sample and hold at 100Hz of the continuous time model.

6.3.4 Validation

With the modified model γ_{min}(T I) and γ_{min}(T V) were recalculated and the results shown in figures 6.16 and 6.17 . The plots show that both γ_{min}(T I) and γ_{min}(T V) are now both less than ϵ_{max} = 1 for a set of models. Hence the set of models with

ϵ_{max} γ_{min}(T V) γ_{min}(T I)  contact point (m) γ_{min}(·)

Figure 6.14: Variation of γ_{min}(·) with contact point

voltage (V)  time (secs)

Figure 6.15: Input (solid) to, and output (dash-dot) from, power amplifier

contact points between about 0.34m and 0.38m, and with time-invariant uncertainty, cannot be invalidated for the given set of data. It is interesting to note that the results indicated by the plots in figure 6.17 match closely the observed behaviour. It was observed in [WWb] that controllers designed to stabilize the beam at a contact point of 0.35m became unstable more quickly when the contact point was decreased, than when it was increased. This is to be expected from the plot in figure 6.17.

It is interesting to note the model that gives the smallest γ_{min}(T I) for the given set of data has a contact point of approximately 0.36m. This is 1cm longer than the actual contact point of the beam. It has already been noted that the second mode predicted by the model occurs at a frequency approximately 1.5Hz lower than is observed in the input-output data. As the contact point of the model is increased, the frequency of the first mode increases and the frequency of the second mode decreases. Hence it is not surprising that a model with a longer contact point than the actual contact point accounts for the data better, as the second mode will be modelled better.

6.4. CONCLUSIONS

γ_{min}(T V) γ_{min}(T I) length of data γ_{min}(·)

Figure 6.16: Variation of γ_{min}(·) with length of data, for modified model

ϵ_{max} γ_{min}(T V) γ_{min}(T I)  contact point (m) γ_{min}(·)

Figure 6.17: Variation of γ_{min}(·) with contact point, for modified model

6.4 Conclusions

We have demonstrated that the modified model of the beam was validated for the given set of data. Even though the model was relatively simple, having only one input and one output, and the data relatively short, the computations required for validation were very demanding. It is therefore reasonable to ask if the model can be validated by using standard identification techniques to obtain a model and seeing if this model lies inside a gap ball . To test this we used the experimental data to identify different models and calculated their gap, and ν-gap, from the the nominal model.

We used the identification toolbox in MATLAB to identify various parametric models [Lju86]. The results are summarized in the following table:

Model structure  Coeffs. ∥n∥_{2} Fit Gap ν-Gap
ARX  13,7,4 0.252 0.026 1  1
ARX  10,8,4 0.246 0.025 1  1
ARMAX  5,5,5,5 0.321 0.031 1  1
Box Jenkins  5,4,5,4,5 0.382 0.038 1  1
Output Error  9,5,5 0.286 0.028 1  1
State-space can. form  7  0.308 0.030 1  1
Optimized nominal model  16  0.227 0.020 1  1

The different model structure are those used in the identification toolbox, with the exception of the optimized nominal model. The optimized nominal model was the modified discrete time state space model, validated in the previous section, with the C and D matrices optimized to minimize the prediction error. The coefficients are the standard coefficients used in [Lju86] to define the model structure, which is the order of the model for state-space models. The noise, n, was taken to be the difference between the simulated output and the measured output and the fit is the mean square error, ie if y is the measured output (of length l) and y_{s} the simulated output then the fit is ∥y − y_{s}∥_{2}/\sqrt{l}.

The results in the table show that we were unable to obtain a valid model, or even get close. This is not surprising as the identification procedures attempted to minimize the prediction error, whereas we were looking for a model that minimizes the gap between the model and the nominal model. In fact the noise levels obtained by the identified models were considerably less than the levels we assumed, where ∥n∥_{2} on the output^{4} was approximately 0.8. Therefore we conclude that standard identification tools do not provide a solution to the computational difficulties of our model validation problem.

Validating and invalidating different models of the beam has illustrated some important points.

• It must be stressed that the γ_{min}(·) calculated are only lower bounds on the size of uncertainty in the system. Our confidence is possibly increased by taking as long a set of input-output data as possible, and checking that γ_{min}(·) appears to converge to a limit as the length of data is increased. However, it is still only a lower bound.

• The validation results obtained for the first model of the beam raise important philosophical questions. The model was modified as the smallest timeinvariant uncertainty, consistent with the data, was larger than could be stabilized. However for this model, the smallest time-varying uncertainty was considerably less than could be stabilized. This leaves the control designer with a dilemma. If he believes that the uncertainty is time-varying, or nonlinear, then the model will not be invalidated by the given set of data. However, if he believes the uncertainty is time-invariant then the model is invalidated by the data. It may be best to compromise and assume the true uncertainty lies between these two limits. It is probably unlikely that the uncertainty behaves as a linear, time-invariant system and it is equally unlikely that it behaves as a possibly arbitrarily fast time-varying or nonlinear system.

4This norm was obtained by filtering the input noise by the plant model P and adding it to the measured output noise.

• It must be emphasized that our validation results depend on the assumptions on the noise. We make the standard assumption in H_{∞} theory, that the noise is bounded in norm. This can be considered as the worst case, where the noise is allowed to act in the worst possible way on the system [Hja93]. Hence, in model validation we can choose the best possible noise, within the norm bound, that minimizes the bound on the uncertainty. Typically the resulting noise sequence does not “look like noise” or appear white. Assuming that the noise is white is a stronger assumption and requires stochastic assumptions on the noise.

Harrier VSTOL Aircraft

In this chapter we use flight test data to validate models of a Harrier VSTOL aircraft. Together with the results of the previous chapter, and results in [WWa] and [Hyd93], this demonstrates the applicability of robust control model-sets to physical systems.

The Harrier is the first production aircraft to have the ability to take-off and land vertically. It accomplishes this by directing the engine thrust through rotatable nozzles. While making the Harrier capable of incredible manoeuvres, the nozzles make it difficult to fly. For example, when decelerating from wingborne flight to the hover, the pilot has to control the pitch of the aircraft with his right hand, and both the engine thrust and nozzle angle with his left hand. Every time the nozzles are rotated downwards, there is a large kick in the pitch which must be compensated by the pilot.

The Harrier can be made easier to fly using feedback control, and in 1990 the Ministry of Defence sponsored the Defence Research Agency’s VAAC Harrier aircraft research programme, to assess different methods of improving the transition to hover of a Harrier like aircraft [SFAH94]. Several organizations designed controllers for this project but the only approach to use multivariable techniques is described in [Hyd91, HG93, Hyd93]. In this approach H_{∞} theory was used for the controller design and in this chapter we consider the validity of the model-sets used in [Hyd93].

7.1 Model

The controller was designed for the Harrier using linearizations of a nonlinear model. Only the longitudinal dynamics were controlled, as these were the critical dynamics in the transition to the hover. The linearizations were at airspeeds of 105, 175, 245 and 345 ft/sec, and the model had three inputs and three outputs. The inputs were longitudinal thrust, vertical thrust, and tailplane angle, and the outputs longitudinal acceleration, vertical acceleration and pitch rate^{1}. The vertical thrust component is obtained using the nozzles to direct the engine thrust downwards, but in level flight with the nozzles fully back, this is zero. The gains of the controllers designed for these linearizations were scheduled to obtain a controller valid for the whole flight envelope.

The controllers were designed in [Hyd93], as for the flexible beam in chapter 6, using H_{∞} loopshaping techniques. Consequently the linearized models, which are all in discrete time, were weighted to incorporate performance objectives. With these weighting functions each linearization had twenty states, and a plot of the singular values of the linearization at 380 ft/sec and 100 Hz sampling frequency is shown in figure 7.1. The smallest singular value is zero. This plot shows that the critical frequencies for the model lie between 0.1 and 20 rad/sec.

1For a discussion of why these variables were chosen for feedback see [Hyd91].

log of singular values frequency (rad/sec)

Figure 7.1: Two largest singular values of weighted linearized model

The model-sets used for controller design were N CF(Ñ, M̃, Δ, γ), with Δ = T I_{∞}^{p×m} or Δ = T V_{∞}^{p×m}, the same as for the flexible beam. We assumed noise entered the model at the same place as the uncertainty output, as in assumption 5.1. Hence we used the same model-set as for the flexible beam, shown in figure 6.4, but without the uncertainty weight^{2}, W_{u}.

We calculated the gap and ν-gap between linearizations at different airspeeds and the linearization at 380ft/sec, and the results are shown in figure 7.2. This plot is very different to the corresponding plot for the flexible beam, in figure 6.5. Whereas the flexible beam only had a small set of contact points where the gap between model was less than ϵ_{max}, the Harrier linearizations have gaps less than ϵ_{max} for a large set of airspeeds. This demonstrates that a certain set of H_{∞} controllers are guaranteed to be robust to changes in the airspeed.

7.2 Experiment design

The data we used for validation was obtained by adding additional inputs into the closed loop of the aircraft and controller. To do this we used a facility in the controller software designed for testing disturbance rejection. A block diagram of the interconnection between the model (or aircraft), represented by the weighted model W_{2}GW_{1}, and the controller K is shown in figure 7.3. This block diagram shows how disturbance signals, d1 and d2, can be injected into the closed loop, and the effects monitored by recording the signals e1 and e2. All of these signals have three separate channels.

2No uncertainty weight was used as its’ theory was not known when the controllers were initially designed.

airspeeds (ft/sec) gap between models ϵ_{max}

Figure 7.2: Variation of gap (solid) and ν-gap (dashed) with airspeed

W_{2}GW_{1} K d1(3)  d2(3) e1(3)  e2(3)

Figure 7.3: Interconnection of controller and model

To provide data for validation we decided to inject a swept sinusoid on each channel of d1 and d2 separately, whose frequencies covered the bandwidths of the model. The sinusoid swept from 0.1 to 20 rad/sec and was

sin(0.1t + \frac{19.9}{108000}t^{4}), 0 ≤ t ≤ 30.  (7.1)

This sinusoid was used, rather than a more conventional sinusoid of the form sin(at + bt^{2}), to ensure that the low frequencies were more fully excited. On a semilogarithmic plot, this means the frequency varies more linearly with time. This is demonstrated by the plots in figure 7.4, which shows how the frequency of two sinusoids, that cover the same frequency range, one of the form sin(a_{1}t + b_{2}t^{2}) and the other sin(a_{2}t + b_{2}t^{4}), vary with time on a semi-logarithmic plot. Both sinusoids begin at a frequency of 0.1 rad/sec and sweep up to 20 rad/sec, but it is clear that the sinusoid used to obtain data for validation, of the form sin(a_{2}t + b_{2}t^{4}), spends more time at low frequency.

We intended to implement the sine sweeps for each channel at three different airspeeds; 85, 255 and 425 ft/sec. With six channels for each airspeed this would have produced nine minutes of data. However, problems with the flight test meant that data was only obtained for two disturbance channels with the Harrier flying at 380 ft/sec.

frequency (rad/sec)  time (secs)

Figure 7.4: Variation of frequency of sin(0.1t+ \frac{19.9}{108000}t^{4}) (solid) and sin(0.1t+ \frac{19.9}{60}t^{2}) (dashed) with time

7.2.1 Computational issues

Computing γ_{min}(T I) proved difficult for the Harrier data because of the length of the data. A complete sine sweep at 100Hz produces 3000 data points, which cannot be handled by present software. Resampling the data at 6.25Hz, which is the smallest possible before aliasing occurs, reduces the number of data points to 190, on three channels. Hence it is necessary to solve LMI problems for a minimum of 570 data points. This proved impossible using LMI-lab, due to the memory required to store the LMI basis.

A cutting plane algorithm was implemented instead, which was far more efficient with memory storage, but for flight test data took prohibitively many iterations to converge. The problem was with the upper bound, which is used to decide when γ_{min}(T I) has been determined to a sufficient accuracy . The upper bound calculation requires the matrix T_{w}^{T} T_{w} to be inverted. For flight test data, the diagonal elements of T_{w} are often small relative to the off-diagonal elements, so T_{w}^{T} T_{w} is ill-conditioned. This makes the upper bound poor and convergence of the cutting plane algorithm slow.

Our solution to this problem was to calculate an upper bound for γ_{min}(T I), which is denoted by γ_{ub}(T I). Recall from section 4.4 that γ_{min}(T I) is the smallest

γ such that there exists [ \begin{matrix}v \\ n \end{matrix} ] satisfying

σ̄(T_{v}(T_{w(v,n)}^{T}T_{w(v,n)})^{− \frac{1}{2}} ) ≤ γ.

Let [ \begin{matrix}v_{0} \\ n_{0} \end{matrix} ] be the value of [ \begin{matrix}v \\ n \end{matrix} ] that minimizes the Frobenius norm^{3} of

T_{v}(T_{w(v,n)}^{T}T_{w(v,n)})^{− \frac{1}{2}} ,

7.3. SIMULATED DATA

subject to an infinity norm bound on n. This is a positive definite quadratic program, which can be solved efficiently [Lue84], and

3If A is a matrix then the Frobenius norm is defined as ∥A∥_{F} := \sqrt{trace(A^{T} A)}.

γ_{ub}(T I) := σ̄(T_{v}(T_{w(v_{0},n_{0})}^{T}T_{w(v_{0},n_{0})})^{− \frac{1}{2}} )

is clearly an upper bound on γ_{min}(T I).

The upper bound turns out to be surprisingly good. For short lengths of flight data, γ_{min}(T I) and γ_{ub}(T I) can be calculated (with an infinity-norm bound on n of 0.02) using a nonlinear optimization routine. The results are summarized in the following table:

Length of data γ_{ub}(T I) γ_{min}(T I)
5  0.1458  0.1365
10  0.1558  0.1510
15  0.2327  0.2218
20  0.3096  0.2761

Note that we can calculate γ_{min}(T V) analytically, using the observation in section 5.1.3, since we use the infinity-norm to bound the noise.

7.3 Simulated data

Before the test flight, we simulated the outputs for the swept sinusoid input using a nonlinear FORTRAN model implemented on a PC. There were several reasons for this:

1. We adjusted the amplitude of the injected sinusoid to give good data for validation. The amplitude was made as large as possible, to ensure good signal to noise ratios, but not as large as to affect the airworthiness of the aircraft, or to hit saturation and rate limits. In practice the maximum amplitude is that which the test pilot will tolerate.

2. We used data from the nonlinear model to calculate likely results from validation.

3. The nonlinear simulated data provides an indication of the validity of the linearization. If the linear models were invalidated using simulated data, then controller design based on linear models is unlikely to be successful.

7.3.1 Measuring noise

We estimated the noise, for both the nonlinear model and aircraft, by recording the signals e1 and e2 for the input on d1(2). The channel d1(2) corresponds to the vertical thrust of the engine, which with the nozzles fully back is zero. Hence the outputs are a realization of the noise on each channel. For the data from the nonlinear simulation, the appropriate signals are shown, after resampling at 12.5Hz (from 100Hz), in figure 7.5. The dimensions of the output are omitted, as in all the plots of this chapter, as the data is classified.

The amplitude of noise is surprising for data from a simulation model. The main source of this noise is thought to be the analogue joystick, used to input the pilot’s commands. Although no pilot commands are used during any of the flight

d1(2)

e1(1)

e1(2)

e1(3)

e2(1)

e2(2)

e2(3)

Figure 7.5: Input on d1(2) to nonlinear model, and corresponding outputs on e1(·) and e2(·) test, except for the initial flight to the suitable airspeed, the joystick produces a considerable amount of electrical noise.

The output e1(2) shows that the vertical acceleration output from the controller is zero. The signals e1(1), e1(3), e2(2) and e2(3) appear to contain significant amounts of noise but this cannot be said for e2(1). This is the longitudinal acceleration of the aircraft and would probably be more accurately modelled as an initial condition. However, the amplitude, relative to the amplitude of d1(2), is small and it will be considered as noise.

To obtain bounds on the noise, we consider noise at the input and output to be filtered by Ñ and M̃ respectively. The maximum amplitude of these filtered signals is then taken as an infinity-norm bound on the noise. The bounds are shown in the following table:

Channel number Infinity-norm bound
1  0.021
2  0.009
3  0.005

The relatively large bound on the noise on the first channel is due to the relatively large signal on e2(1).

7.3.2 Validation

Using the norm bounds on the noise described, we used the data from the nonlinear simulation to test the validity of linearizations as different airspeeds. The data and validation results are only shown for inputs on d1, as these correspond to the results from the flight test data. The data from a simulation, at an airspeed of 430 ft/sec, is shown in figures 7.5, 7.7 and 7.8, and the validation results in figures 7.6 and 7.9.

γ_{min}(T V) γ_{ub}(T I) length of data norm bound

Figure 7.6: Variation of γ_{min}(T V) and γ_{ub}(T I) with data length for data from nonlinear model

The plots in figure 7.6 shows how γ_{ub}(T I) and γ_{min}(T V) vary with data length. Both γ_{ub}(T I) and γ_{min}(T V) appear to be converging to an upper bound giving some

d1(1)

e1(1)

e1(2)

e1(3)

e2(1)

e2(2)

e2(3)

Figure 7.7: Input on d1(1) to nonlinear model, and corresponding outputs on e1(·) and e2(·)

d1(3)

e1(1)

e1(2)

e1(3)

e2(1)

e2(2)

e2(3)

Figure 7.8: Input on d1(3) to nonlinear model, and corresponding outputs on e1(·) and e2(·)

ϵ_{max} γ_{min}(T V) γ_{ub}(T I) norm bound  Validation for input on d1(1) γ_{min}(T V)γ_{ub}(T I) norm bound  Validation for input on d1(2) ϵ_{max} γ_{min}(T V) γ_{ub}(T I) norm bound  Validation for input on d1(3)

Figure 7.9: Variation of γ_{min}(T V) and γ_{ub}(T I) with airspeed for data from nonlinear model confidence that γ_{ub}(T I) and γ_{min}(T V) are good lower bounds on the gap between plant and system.

In all of the plots in figure 7.9, both γ_{ub}(T I) and γ_{min}(T V) are considerably smaller than ϵ_{max} for most of the airspeeds considered. This a good validation result and indicates that there is no evidence that the linear model-sets to not adequately describe the behaviour of the nonlinear model.

The plot in figure 7.9 for an input on d1(2) is to be expected, with many linearizations matching the data. This is because the outputs were assumed to be noise.

The plot in figure 7.9 for an input on d1(3) is more interesting. The shape of the plot seems reasonable, with a linearization close to 430 ft/sec best accounting for the data, but one may expect γ_{min}(T V) to get closer to zero at 430 ft/sec. The difference between γ_{min}(T V) and zero may be due to the nonlinearities. If this is the case then one would expect the difference to get smaller if the size of the inputs is decreased. A plot of the variation of γ_{min}(T V) for an input of half the size is shown in figure 7.10, which shows a gap between γ_{min}(T V) and zero of about half the size of the original input. This indicates that the discrepancy may be due to the nonlinearities in the model.

ϵ_{max} γ_{min}(T V) γ_{ub}(T I)  airspeed (ft/sec) norm bound  Validation for input on d1(3)

Figure 7.10: Variation of γ_{min}(T V) and γ_{ub}(T I) for smaller input signal

The plot in figure 7.9 for the input on d1(1) is a little surprising. Although the true airspeed is 430 ft/sec the results in this plot indicate that the 480 ft/sec linearization accounts for the data better. The reason for this poor match may be the poor signal to noise ratio on this set of data. The approximate size of the noise can be estimated from figure 7.5 and this can be compared with the signals in figure 7.7. This shows that noise can account for between 10% and 50% of the signal, which is clearly poor. However, with the signals shown in figure 7.8 the noise only accounts for between 1% and 20% of the signal. Hence results based on the data shown in figure 7.8 should be expected to be more reliable than for the data in figure 7.7.

Together these results support the design of controllers using linearizations of a nonlinear model. The model-sets derived from the linearizations have been shown to contain elements that exactly describe the data obtained from the nonlinear model. This supports the approach to controller design in [Hyd91] and the subsequent work.

7.4. FLIGHT-TEST DATA

7.4 Flight-test data

Several Harrier test flights have been flown at the Defence Research Agency, Bedford, England between December 1993 and March 1994. The final one of these test flights included the validation tests described in the previous section. Unfortunately not all of the planned tests were completed and only two complete sets of data were obtained; for inputs on d1(2) and d1(3), and at an airspeed of 380 ± 10 ft/sec.

7.4.1 Measuring noise

We again assumed that the outputs on e1 and e2 for the input on d1(2) were a realization of the noise, and the signals are shown in figure 7.11. It is clear from these plots that there is considerably more noise in the flight-test data than in the data taken from the nonlinear model. This is not surprising as there were many disturbances acting on the plane that were not modelled in the nonlinear model, for example wind gusts. Again we filtered the signals e1 and e2 by Ñ and M̃ respectively, to obtain bounds on the noise. The norm bounds are shown in the following table:

Channel number Infinity-norm bound
1  0.043
2  0.064
3  0.035

These bounds are all at least twice as large as the norm bounds assumed for the nonlinear model.

7.4.2 Validation

The only successful flight test was for the input on channel d1(3). Only the first 27 seconds of this sweep was recorded, but this still covers frequencies up to approximately 15 rad/sec. The original data was sampled at 102.4Hz and this was resampled at 12.8Hz to provide a manageable length of data. This resampled data is shown in figure 7.12.

The results of validation are shown in figures 7.13 and 7.14. The results are good with a large range of airspeeds giving γ_{ub}(T I) and γ_{min}(T V) less than ϵ_{max}. The minimum of γ_{ub}(T I) does occur at an airspeed of 405 ft/sec when the true airspeed is closer to 380 ft/sec. However, γ_{ub}(T I) is only an upper bound and the number of linearizations used was not great enough to determine the airspeed to a greater accuracy than plus or minus 35 ft/sec.

d1(2)

e1(1)

e1(2)

e1(3)

e2(1)

e2(2)

e2(3)

Figure 7.11: Input on d1(2) to the aircraft, and corresponding outputs on e1(·) and e2(·)

d1(3)

e1(1)

e1(2)

e1(3)

e2(1)

e2(2)

e2(3)

Figure 7.12: Input on d1(3) to the aircraft, and corresponding outputs on e1(·) and e2(·)

γ_{min}(T V) γ_{ub}(T I) length of data norm bound

Figure 7.13: Variation of γ_{min}(T V) and γ_{ub}(T I) with data length

ϵ_{max}  γ_{ub}(T I) γ_{min}(T V)  airspeed (ft/sec) norm bound

Figure 7.14: Variation of γ_{min}(T V) and γ_{ub}(T I) with airspeed

Conclusions

We conclude the thesis by summarising the contributions and suggesting possibilities for further research. The contributions fall into three main areas; extending time domain validation results to include more general model-sets, providing a new interpretation of existing validation results in the gap metric and deriving results for validation in the ν-gap metric, and applying the theory to the model and data from a flexible beam and Harrier.

8.1 Contributions

• The time domain results described in [PKT^{+}94, ZK92] have been extended to include LFT model-sets with noise and structured uncertainty. LFT modelsets with structured uncertainty were considered in the frequency domain approach of Smith and Doyle, but the uncertainty was not necessarily causal. In our time domain results the uncertainty is guaranteed to be causal.

• We have proved that the general model validation problem is NP-hard in the number of uncertainty blocks. This indicates that it is unlikely that we will find a computationally tractable algorithm that solves the general model validation problem.

• We have derived conditions for various model validation problems to be convex. This is useful when deciding the model validation problems that are likely to be computationally tractable, especially as realistic problems are likely to have a lot of data.

• In chapter 5 we analysed left and right NCF models and showed that they are not equivalent for model validation. This was shown to be because of the difference between balls defined in the gap and T-gap metrics. The gap metric also provides an interesting interpretation of NCF model-set validation.

• The most important theoretical result of the thesis is the extension of the validation results to include model-sets defined as balls in the ν-gap metric. Results for validating balls in the ν-gap can be considered as the best possible validation results for robust control, since a ν-gap ball characterizes the largest set of systems that can a priori be guaranteed to be stabilized by a certain set of H_{∞} controllers.

8.2. FURTHER RESEARCH

• The validation of models of the flexible beam and Harrier, in chapters 6 and 7, demonstrates the applicability of the model-sets used for robust design to real systems. Validation of the flexible beam revealed deficiencies in the model which were overcome by a modifying the model. The behaviour predicted by validation was also born out in practice.

8.2 Further Research

There is much scope for further research, with the two obvious areas being in the application of current validation techniques, and the theory. There is much that can be done in the algorithms for model validation. The amount of data that was analysed in chapters 6 and 7 was severely restricted by computational limitations, but this should not be the case. The validation problems solved for the model-sets and data from the flexible beam and Harrier were all convex feasibility problems, for which interior point algorithms are claimed to work very well. It is claimed in [BVG94] that interior point algorithms have been developed for which

… we are able to solve convex optimization problems with over 1000 variables and 10000 constraints in around 10 minutes on a workstation.

We have achieved no where near this level of performance for the model-sets and data for the flexible beam and Harrier, but this work indicates that significant improvements are possible.

Recent work in [RP] has solved model validation problems for sample-data uncertainty models. The models included additive uncertainty and additive noise, so a natural extension of this work is to include more general uncertainty descriptions.

A useful theoretical result would be necessary and sufficient conditions for validating model-sets where the uncertainty has a bounded rate of time variation. Theoretical results for systems where the uncertainty has a bounded rate of variation have recently been obtained in [WP, PT], and it may be possible to apply these results to model validation problems. A way of measuring the rate of variation of an operator Δ is to measure the norm of the commutant of Δ and the right shift S, ie ∥SΔ − ΔS∥. Therefore the interpolation result that would allow results similar to this thesis is:

Given sequences u ∈ π_{l}Sm_{+}, y ∈ π_{l}Sn_{+}, and positive real numbers γ and δ, when does there exists a stable, causal, linear, time-varying operator Δ satisfying

∥Δ∥_{i2} ≤ γ,
∥SΔ − ΔS∥ ≤ δ,
y = π_{l}Δu?

Bibliography

[ÅE71] K J Åström and P Eykhoff. System identification – a survey. Automatica, 7:123–162, 1971.

[Bai94] N T J Bailey. Prediction and validation in the public-health modeling of HIV/AIDS. Statistics in Medicine, 13:1933–1943, 1994.

[BBFE93] S Boyd, V Balakrishnan, E Feron, and L ElGhaoui. Control system analysis and synthesis via linear matrix inequalities. In Proc. American Cont. Conf., pages 2147–2154, 1993.

[BDG^{+}] G Balas, J C Doyle, K Glover, A Packard, and R Smith. MU-TOOLS Manual. MUSYN INC. and The MathWorks Inc.

[BEFB94] S Boyd, L Elghaoui, E Feron, and V Balakrishnan. Linear Matrix Inequalities in System and Control Theory. SIAM, 1994.

[BT86a] M H Butterfield and P J Thomas. Methods of quantitative validation for dynamic simulation models - part 1: Theory. Transactions of the Institute of Measurement and Control, 8(4):182–200, October 1986.

[BT86b] M H Butterfield and P J Thomas. Methods of quantitative validation for dynamic simulation models - part 2: Examples. Transactions of the Institute of Measurement and Control, 8(4):201–219, October 1986.

[BVG94] S. Boyd, L. Vandenberghe, and M. Grant. Efficient convex optimization for engineering design. In Proceedings IFAC Symposium on Robust Control Design, pages 14–23, 1994.

[BYDM94] R D Braatz, P M Young, J C Doyle, and Manfred M. Computational complexity of μ calculation. IEEE Trans. Automat. Contr., 39(5):1000– 1002, 1994.

[Cor89] G O Corrêa. A system identification problem motivated by robust control. Int. J. Control, 50(2):575–602, 1989.

[DFT92] J Doyle, B A Francis, and A R Tannenbaum. Feedback Control Theory. Mcmillan Publishers, 1992.

[DGKF89] J C Doyle, K Glover, P P Khargonekar, and B A Francis. State-space solutions to standard H_{2} and H_{∞} control problems. IEEE Trans. Automat. Contr., 34(8):831–847, August 1989.

[Doy82] J C Doyle. Analysis of feedback systems with structured uncertainty. IEE Proc., 129(6), November 1982.

[DP87] M Dahleh and J B Pearson. l^{1}-optimal controllers for MIMO discretetime systems. IEEE Trans. Automat. Contr., 32:314–323, 1987.

[dSP] dSPACE. dSPACE manual. dSPACE GmbH, Germany.

[Dym89] H Dym. J-Contractive Matrix Functions, Reproducing Kernel Hilbert Spaces and Interpolation. American Math. Soc., 1989.

[FF90] C Foias and A E Frazho. The Commutant Lifting Approach to Interpolation Problems. Birkhäuser, 1990.

[GD88] K Glover and J C Doyle. State-space formulae for all stabilizing controllers that satisfy an H_{∞}-norm bound and relations to risk sensitivity. Systems and Control Letters, 11:167–172, 1988.

[Geo88] T T Georgiou. On the computation of the gap metric. Systems and Control Letters, 11:253–257, 1988.

[GGLD90] M Green, K Glover, D Limebeer, and J Doyle. A J-spectral factorization approach to H_{∞} design. SIAM J. Control and Optimization, 28(6):1350–1371, 1990.

[GGN92] G C Goodwin, M Gevers, and B Ninness. Quantifying the error in estimated transfer functions with application to model order selection. IEEE Trans. Automat. Contr., 37(7):913–928, 1992.

[GJ79] M R Garey and D S Johnson. Computers and Intractability - A Guide to the Theory of NP-Completeness. W H Freeman and Co., 1979.

[GK92] G Gu and P P Khargonekar. Linear and nonlinear algorithms for identification in H_{∞} with error bounds. IEEE Trans. Automat. Contr., 37(7):953–963, 1992.

[Glo86] K Glover. Robust stabilization of linear multivariable systems: relations to approximation. Int. J. Control, 43(3):741–766, 1986.

[GN93] P Gahinet and A Nemirovskii. LMI-lab: A package for Manipulating and Solving LMI’s, 1993.

[GS90] T T Georgiou and M S Smith. Optimal robustness in the gap metric. IEEE Trans. Automat. Contr., 35(6):673–686, 1990.

[HG93] R A Hyde and K Glover. The application of scheduled H_{∞} controllers to a VSTOL aircraft. IEEE Trans. Automat. Contr., 38(7):1021–1039, 1993.

[HJ85] R A Horn and C R Johnson. Matrix Analysis. Cambridge University Press, 1985.

[Hja93] H Hjalmarsson. Aspects on Incomplete Modeling in System Identification. Phd thesis, Linköping University, 1993.

[HJN91] A J Helmicki, C A Jacobsen, and C N Nett. Control oriented system identification: A worst-case/deterministic approach in H_{∞}. IEEE Trans. Automat. Contr., 36(10):1163–1176, 1991.

[Hyd91] R A Hyde. The Application of Robust Control to VSTOL Aircraft. Phd thesis, University of Cambridge, 1991.

[Hyd93] R A Hyde. CONTROL LAW 005: A robust multivariable control law for the DRA VAAC programme. Technical report, Cambridge University Engineering Department CUED/F-INFENG/TR159, UK, 1993.

[Kai80] T Kailath. Linear Systems. Prentice-Hall, 1980.

[KGED94] L W Keene, D B Gibbons, D J Evans, and T W Davies. Computer modeling validation for shell-and-tube heat-exchangers. Chem. Eng. Res. and Des., 72:611–615, 1994.

[LDF94] M M Livstone, M A Dahleh, and J A Farrell. A framework for robust control based model invalidation. In Proc. American Cont. Conf., pages 3017–3020, 1994.

[Lju86] L Ljung. System Identification Toolbox: User’s Guide. The Mathworks Inc., 1986.

[Lju87] L Ljung. System Identification: Theory for the User. Prentice Hall, 1987.

[Lue84] D G Luenberger. Linear and Nonlinear Programming. Addison-Wesley, 1984.

[Mäk91] P M Mäkilä. Laguerre methods and H_{∞} identification of continuous time systems. Int. J. Control, 53(3):689–707, 1991.

[MFM94] R Malewski, M A Franchek, and V T McWhirter. Experimental validation of a computer-model simulating an impulse voltage distribution in HV transformer windings. IEEE Trans. Power Deliv., 9:1789–1798, 1994.

[MG90] D C McFarlane and K Glover. Robust Controller Design Using Normalized Coprime Factor Plant Descriptions. Springer-Verlag, 1990.

[MK87] K G Murty and S N Kabadi. Some NP-complete problems in quadratic and nonlinear programming. Mathematical Programming, 39:117–127, 1987.

[MR85] M Milanese and Tempo R. Optimal algorithms theory for robust estimation and prediction. IEEE Trans. Automat. Contr., 30:730–738, 1985.

[NS91] M P Newlin and R Smith. Model validation and a generalization of μ. In Proc. IEEE Conf. Dec. and Cont., pages 1257–1258, 1991.

[Pac88] A Packard. What’s New with μ: Structured Uncertainty in Multivariable Control. Phd thesis, University of California, Berkeley, 1988.

[Par88] J R Partington. An introduction to Hankel operators. Cambridge University Press, 1988.

[Par91] J R Partington. Robust identification and interpolation in H_{∞}. Int. J. Control, 54(5):1281–1290, 1991.

[Par93] J R Partington. Algorithms for identification in H_{∞} with unequally spaced function measurements. Int. J. Control, 58(1):21–31, 1993.

[PKT^{+}92] K Poolla, P P Khargonekar, A Tikku, J Krause, and K M Nagpal. A time-domain approach to model validation. In Proc. Amer. Cont. Conf., pages 313–317, 1992.

[PKT^{+}94] K Poolla, P Khargonekar, A Tikku, J Krause, and K Nagpal. A timedomain approach to model validation. IEEE Trans. Automat. Contr., 39(5):951–959, May 1994.

[PT] K Poolla and A Tikku. Robust performance against time-varying structured perturbations. To appear in IEEE Trans. on Aut. Cont.

[RL92] N P Rubin and D J N Limebeer. System identification for H_{∞} control. In Proc. IEEE Conf. Dec. and Cont., pages 1694–1695, 1992.

[Roc93] R T Rockafellar. Lagrange multipliers and optimality. SIAM Review, 35:183–283, 1993.

[RP] S Rangan and K Poolla. Time-domain validation for sample-data uncertainty models. Submitted to IEEE Trans. on Aut. Cont.

[SD90] R S Smith and J C Doyle. Towards a methodology for robust parameter identification. In Proc. American Cont. Conf., pages 2394–2399, 1990.

[SD92] R S Smith and J C Doyle. Model validation: A connection between robust control and identification. IEEE Trans. Automat. Contr., 37(7):942–952, July 1992.

[SFAH94] G T Shanks, C Fielding, S J Andrews, and R A Hyde. Flight control and handling research with the VAAC harrier aircraft. Int. J. Control, 59(1):291–319, 1994.

[Smi90] R Smith. Model Validation for Uncertain Systems. Phd thesis, California Institute of Technology, 1990.

[SVdH92] R J P Schrama and P M J Van den Hof. An iterative scheme for identification and control design based on coprime factorizations. In Proc. American Cont. Conf., pages 2842–2846, 1992.

[vdB93] T van den Boom. MIMO system identification for H_{∞} robust control. Phd thesis, Technical University of Delft, 1993.

[Vid87] M Vidyasagar. Control System Synthesis: A Factorization Approach. MIT Press, 1987.

[Vina] G Vinnicombe. The best possible results for H_{∞} controllers, with application to identification. Submitted to IEEE Trans. on Aut. Cont.

[Vinb] G Vinnicombe. The robustness of feedback systems with bounded complexity controllers. Submitted to IEEE Trans. on Aut. Cont.

[Vin92] G Vinnicombe. Measuring the Robustness of Feedback Systems. Phd thesis, University of Cambridge, 1992.

[Vin93] G Vinnicombe. Frequency domain uncertainty and the graph topology. IEEE Trans. Automat. Contr., 38(9):1371–1383, 1993.

[WAK94] G W Wall, J S Amthor, and B A Kimball. Cotco2 - a cotton growth simulation-model for global change. Agricultural and Forest Meteorology, 70:289–342, 1994.

[WL92] B Wahlberg and L Ljung. Hard frequency-domain model error bounds from least-squares like identification techniques. IEEE Trans. Automat. Contr., 37(7):900–912, 1992.

[Woo93a] J Woodhouse. On the playability of violins. part i: Reflection functions. Acustica, 78:125–136, 1993.

[Woo93b] J Woodhouse. On the playability of violins. part ii: Minimum bow force and transients. Acustica, 78:137–153, 1993.

[WP] G Wolodkin and K Poolla. A power distribution theorem for robust control. Submitted to Automatica.

[WWa] G D Wood and J Woodhouse. Experimental validation of robust control schemes for a flexible manipulator. To be submitted to Int. Jnl. Cntrl.

[WWb] G D Wood and J Woodhouse. Robust force control through a flexible manipulator. To be submitted to Int. Jnl. Cntrl.

[You88] N Young. An Introduction to Hilbert space. Cambridge University Press, 1988.

[You93] P M Young. Robustness with Parametric and Dynamic Uncertainty. Phd thesis, California Institute of Technology, 1993.

[Zam81] G Zames. Feedback and optimal sensitivity: Model reference transformations, multiplicative semi-norms, and approximate inverses. IEEE Trans. Automat. Contr., 26(2):301–320, 1981.

[ZDG] K Zhou, J Doyle, and K Glover. Robust and optimal control. Draft.

[ZES80] G Zames and A K El-Sakkary. Unstable systems and feedback: The gap metric. In Proc. Allerton Conf., pages 380–385, 1980.

[ZK92] T Zhou and H Kimura. Minimal H_{∞}-norm of transfer functions consistent with prescribed finite input-output data. In Proc. SICE, pages 1079–1082, 1992.

[ZK93] T Zhou and H Kimura. Time domain identification for robust control. Systems and Control Letters, 20:167–178, 1993.

[ZK94] T Zhou and H Kimura. Simultaneous identification of nominal model, parametric uncertainty and unstructured uncertainty for robust control. Automatica, 30(3):391–402, 1994.