A More Practical Formulation for Robustness Analysis

Copyright 10 IFAC System Structure and Control, Nantes, France, 1998

A MORE PRACTICAL FORMULATION FOR ROBUSTNESS ANALYSIS Gerard Scorletti *,1

* CESAME,4 av. G. Lemaitre, B1348 Louvain-La-Neuve,

Belgium e-mail:

scorlett~ensta.fr

.

Abstract: The structured singular value J.L has proved its usefulness to assess the stability of a linear time invariant model submitted to structured perturbations (real parametric perturbations or model perturbations). To the purpose, the model must be cast into a special form: a Linear Time Invariant system connected with a block diagonal uncertainty operator whose gain is less than one. For practical issue, such a formulation is nevertheless awkward when considering some problems of interest such that one sided stability analysis or distance to passivity. We propose a new formulation, more general and flexible, which allows to address these problems in a direct way. This extension relies on the sector notion developed by Safonov in the beginning of the eighties. For the computation of the margin associated to our formulation, we propose an extension of the J.L-upper bound. Copyright© 1998 IFAC Keywords: Uncertain linear systems, Robust stability, Stability analysis, Multipliers, Convex optimisation

1. INTRODUCTION

lower estimates of km (Fan et al. 1991) can be efficiently computed, in the form of a generalized eigenvalue problem, which is a quasi convex OJr timization problem, efficiently computed (Boyd and El Ghaoui 1993). Performance, measured by a weighted Hoo norm (Zames 1981), is analyzed by introducing a fictitious uncertainty and applying robustness analysis framework (Packard and Doyle 1993). From a practical point of view, it may be interesting to make a difference between the ''level'' of performance and the ''level'' of uncertainty that the system can achieve. This leads to the introduction of a related tool to the J.L-analysis: the v-analysis (Fan and Tits 1992), or skewed J.Lanalysis (Ferreres and Fromion 1996). This problem can be computed by solving a sequence of J.L analysis problems. A more efficient approach is to directly formulate it as one generalized eigenvalue problem (Fan and Tits 1992, Ferreres and Fromion 1996).

A physical system is classically modeled as a nominal linear time invariant system submitted to sector uncertainties a (Doyle et al. 1982, Safonov and Athans 1981). Under certain conditions, such a system can be modeled by a causal linear time invariant system connected by feedback with a (block) diagonal causal linear time invariant oJr erator or a nonlinear operator. Using loop shifting (Zames 1966, Willems 1969), the uncertain operator is ''normalized'' in such a way that the transformed operator has its Hoo norm less or equal one. For this formulation, one may assess the stability of the interconnection or determine the maximum "size" (termed km or 1/ J.L) of the perturbation which destabilizes the system. The exact calculation of this size is known to be an NPhard problem (Braatz et al. 1994). Nevertheless, 1 Part of this research was done when the author was with LMA, ENSTA and Service Auto, ESE, France. At now, he is with CESAME.

107

10 this paper, we propose to consider that the uncertainty ~ is directly described by a sector, linearly dependent on a parameter k. The parameter k can be interpreted as its size. Our purpose is to extend both J.L and v analysis to address related problems such as one-sided multi variable stability margin (Tekawy et al. 1992) or distance to passivity for uncertain systems (Wen 1988). These problems could be reduced to solve a sequence of J.L analysis problems, that is, for the lower bound computation of k, a sequence of generalized eigenvalue problems. Nevertheless, we propose a direct formulation in the form of one generalized eigenvalue problem. To the best of our knowledge, this practical computation was never considered.

two matrices A and B, diag(A, B) denotes the matrix [ ~ ~ ] .

2. AN UNCERTAIN SYSTEM

10 this section, we propose an approach to the modeling of uncertain systems. Let an uncertain system be modeled by the connection of a causal linear time invariant operator M connected with a time invariant causal operator ~

q(s) = M(s)p(s) and p(s) =

~(s)q(s)

(1)

The operator M models the known part of the system. 10 the sequel, this system is referred as the .M-~ system. We assume that the closed loop I. is well-posed. The ~ models the uncertain ' part of the system and it is characterized by a dissipative-like property. More precisely, ~(s) is an n x n transfer matrix belonging to the set a defined as: a = {~(s) ~(s) ~ diag (~l (s), . .. , ~i(S), ... , ~r(s)) } where ~i(S) is characterized by bounded measurable functions Xi, Yi and Zi from jR to en.: x n.: (with Xi(jW) = Xi(jW)* and Zi(jW) = Zi(jW)* ::; 0) in the following way: with Pi(jW) =

Such an approach avoids the normalization step: it is then more practical to use. The uncertainties are considered in their "original" form, that is, modeled as a (sector bounded) operator described by a quadratic constraint. 10 this way, it relies on the notion of sector introduced by Zames (Zames 1966) and developed by Safonov in the beginning of the eighties (Safonov 1980). The main idea is to deal with the ''natural'' form of the model and not with a mathematical equivalent (small gain or passivity) form without explicit physical meaning. Uncertainty description based on quadratic constraints was already used for improving robustness analysis with uncertain time delays (Scorletti 1997). The main result of this article is to establish a more adequate stability margin extending the J.L/km stability margin (Doyle et al. 1982, Safonov 1982). The paper is organized as follows . The problem under consideration is presented in the section 2. 10 contrast with usual robustness analysis, the uncertainties are modeled by quadratic constraints on their input-output signals instead of an Hoc norm. 10 this framework, the usual J.L upper bound is reinterpreted in the section 3 and in the section 4. These two sections focus on the following problem: "Is the uncertain system stable ?". 10 the sequel, this is referred as the stability analysis problem. The main result is stated in the section 5, where the following problem is considered: "what is the maximum "size" of one or several uncertainties for which the uncertain system is stable ?". 10 the sequel, this is referred as the stability margin analysis problem. The computational aspects are briefly discussed in the section 6. The approach is finally illustrated by a numerical example in the following section.

I

~i (jW)qi(jW)

qi(jW)]* [Xi(jW) Yi(jW)] [qi(jW)] >0(2) [ Pi(jW) Pi(jW) Yi(jw)* Zi(jW) This is referred as a quadratic constraint. 10 the sequel, ~ belongs to a set defined by the following quadratic inequality: with p(jw) = ~(jw)q(jw)

q(jW)]* [X(jW) Y(jW)] [q(jW)] >0(3) [ p(jw) Y(jw)* Z(jw) p(jw) where

X(jw) = diag(Xl (jw), .. . , Xi(jW), ... , Xr(jw)) , Y(jw) = diag (Y1 (jw), ... , Yi(jw) , . .. , Yr(jw)) , Z(jw) = diag (Zl(jW), .. . , Zi(jW), ... , Zr(jW)). 10 fact, ~i(S), Xi(jw), Yi(jw) and Zi(jW) may have a special structure, depending on the nature of the uncertainty under consideration. We are specially interested in the two following cases: ~i is a (unstructured) model uncertainty: so it is a multi-input multi-output linear time invariant system only known to satisfy (2), and Xi(jw), Yi(jw) and Zi(jW) are matrices of same dimensions, • ~i is a parameter uncertainty: ~i = bi1n.: where bi is a real time invariant scalar belonging to some interval and Xi(jw), Yi(jw) and Zi (jw) are also repeated scalars on the diagonal of the same dimension. •

Notations. 11 .11 denotes the £2-norm. Ir and Or denote the identity and the zero matrices of Rrxr . The subscript is omitted when it is evident from the context. it denotes R U { -00, +oo}. For the

Note that X = I, Y = 0 and Z = -I express that the operator ~ has an £2-gain less or equal 1 and

108

that X = 0, Y = I and Z operator Ll is passive.

= 0 express that the

• the sum of the numbers of non stable poles of M (s) and of Ll (s) is constant for every Ll (s) in the set A, • every Ll(s) in the set A have the same poles on the imaginary axis, • there exists a Llo in the set A such that the closed loop system is stable.

Now, we propose several models Xi(jW), Yi(jw) and Zi (jw) of the perturbation Lli with an associated "size" ~. This parameter is a scaling factor allowing to narrow or to extend the set of perturbations to which ~ belongs. Its introduction is convenient to state in the section 5 a precise definition of our stability margin. To the set A is then associated the vector: k~ = [kl ... kr ] . We consider the case when Lli = Oiln. with Oi is a real time invariant scalar, that is, a real parametric uncertainty. The case when Lli is a model perturbation can be considered as an extention. Real parametric perturbations are considered such that:

Then the system (1) is stable if there exist a strictly positive t: such that for all w in it: [

4. CONSTRAINT PARAMETERIZATION AND THE STABILITY ANALYSIS PROBLEM

The characterization of the linear operator Ll by the Quadratic Constraint (3) is not unique. By this characterization, the set of signals P, q such that p = Llq are "embedded" in a bigger set of signals Pe, qe related by (3) . As a matter offact, we obtained a very crude description of the operator Ll. Some important informations about Ll are lost: (i) informations about the structure of Ll: it is block diagonal with (possibly) repeated elements, (ii) information about the nature of Ll: some block operators Lli are the multiplication by a repeated real scalar. To avoid this loss of information, a (partial) parameterization of the Quadratic Constraints checked by Ll is proposed. A rule of thumb is to think that the more Quadratic Constraints you consider, the more precise description of the operator Ll you get. To this purpose, scaling functions S and G from jR to c nxn (Doyle 1982, Doyle 1985, Fan et al. 1991) are introduced. Both of them give informations about the (block diagonal with possibly repeated blocks) structure of Ll. Furthermore, the function G gives information about the real nature of some repeated blocks. More precisely, S is such that S(jw) = S(jw)* > 0 and Ll(jw) = S(jW)-1/2Ll(jw)S(jW)1/2 and the function G is such that with p(jw) = Ll(jw)q(jw),

In this section, we give a precise definition of the stability analysis problem. For this problem,

[1 ... 1].

Definition 3.1. The system M -Ll is stable if the

[_~

-;W]

[X(jW) Y(jW)] [M(jW)] . . . Y(jw)* Z(jw) I

The conditions proposed in this section are sufficient. Tighter conditions can be obtained by a (partial) parameterization of the Quadratic Constraints describing Ll as seen in the next section.

3. THE STABILITY ANALYSIS PROBLEM

operator

*

... ::; -tM(jw)* M(jw)

• if Oi E le - r,c + r] with c and r two constant real scalars then Xi = (r2 - c2)I... , Yi = cl... , Zi = -I.... In the sequel, it is referred as a parameter conic sector perturbation (c, r) . This perturbation set can be narrowed or extended by considering the perturbation (c, y'k;r) where ~ is a positive scalar. • if Oi E la, b] with a and b two constant real scalars then Xi = -2abI",., Yi = (a + b)I... , Zi = - 2In.. In the sequel, it is referred as a parameter interval sector perturbation (a, b). This perturbation set can be narrowed or extended by considering the perturbation (a+ (1- ~)(b-a), b) or (a, b+ (1- ~)(a - b» where k i is a positive scalar. • if Oi E la, 001 with a a constant real scalar then Xi = -2aI",., Yi = In., Zi = On•. In the sequel, it is referred as a parameter semiinterval sector perturbation (a). This perturbation set can be narrowed or extended by considering the perturbation (a + (1 - ki)b) where b and k i are positive scalars.

k~ =

M(jW)] I

has a stable inverse in £2.

p(jw)*G(jw)q(jw) + q(jw)*G(jw)*p(jw) =

o.

The previous properties of the functions S and G imply that (1) if Lli is a multi-input multi-output linear time invariant system then they are of the form Di(jW) = d;(jw)I... where d;(jw) is a scalar function, (2) if ~ is a Oil... where Oi is a real time invariant scalar, then they are of the form Di(jW) is an 11i x 11i matrix function.

The stability analysis problem consists in checking the stability of the system (1). Theorem 3.1. Suppose that the uncertain operator Ll belongs to the set A. It thus verifies the Quadratic Constraints (2). Furthermore, suppose that:

109

Theorem 4.1. Suppose that the uncertain operator ~ is characterized by (3) . Then, it is also characterized by the Quadratic Constraint:

In this case, the definition 5.1 reduces to the definition of 1/v (Fan and Tits 1992). At last, consider a system M -~ with ~ an interval sector perturbation (0, I) and with k4 defined by as above. In this case, the definition 5.1 reduces to the definition of the one-sided multi variable stability margin (Tekawy et al. 1992). We just illustrate this margin in the case of a linear time invariant system with one uncertain parameter. In the p,-analysis approach, after normalization, the biggest k can be estimated such that the system is stable for an uncertain parameter in an interval [-k, k]. In fact, it may be much more interesting to compute the biggest uncertain scalar k lying in an interval, say, [0, k] such that the system is still stable. Addressing this skewed p,-analysis problem with a standard approach is quite intricate.

q(jW)]' [ s{jw)x{jw) Y(jw)' s{jw)* + G{jw)' . .. [ p(jw) ... S(jW)Y(jW)+G(jW)] [q{jW)] >0 s{jw)z(jw) p{jw)-

for all functions S and G defined above. The main idea is to exhibit from the set of Quadratic Constraints parameterized by the previous theorem 4.1 a Quadratic Constraint for which the sufficient condition of the theorem 3.1 holds. The main issue is thus to find the functions S and G which characterize this Quadratic Constraint (if they exist) . The theorem 3.1 is combined with the theorem 4.1 to obtained the following theorem:

The exact resolUtion of all these problems is known to be NP-har-rl. So, we propose a numerical tractable lower bound of the kmargin. It is based on the Theorem 4.2: it consists on analyzing the stability of the uncertain system with the partially scaled perturbation. This perturbation is then characterized by the matrices X(jw), Y(jw) and Z(jw) which can be decomposed such that:

Theorem 4.2. Suppose that the uncertain operator ~ is such that the assumptions of the Theorem 3.1 hold. Then the system (1) is stable if there exist two functions S and G defined as above and an E > such that Vw E R

°

Xo(jw) YO(jW)] [ Yo(jw)' Zo{jw)

[ MVw) ]. M [ MVw) ] ::; -!M(jw), M(jw) M = [

S(jw)X(jw) S(jw)Y(jw) + G(jW)] Y(jw)' S(jw)' + G(jw)' S(jw)Z(jw)

To sum up, proving the stability of the uncertain system (1) reduces to find two functions S and G such that the condition (4.2) holds. As a matter of fact, this problem is an infinite dimensional linear optimization problem. Computational aspects are discussed in the section 6.

+ k .... argin

[Xk{jW) Yk(jW)] Yk(jW), Zk(jW)

Theorem 5.1. Suppose that the uncertain operator ~ is such that the assumptions of the Theorem 3.1 hold. A lower bound on the stability margin is then given by the optimization problem: minimize S{jw), G(jw) subject to

kmargin -1

kmargin -1 B -

A

2: 0

(4)

B>O

with B:

5. THE STABILITY MARGIN ANALYSIS PROBLEM

- [

M{jW)]' [ Xo(jw)S(jw) I Yo(jw)' S{jw) + G{jw)* . . . Yo{jw)S{jw) + G(jW)] [M{jW)] Zo{jw)S(jw) I

In this section, we consider the system M -~. Its associated size vector k4 has two kind of components: either k i = 1 or ~ = kmargin where kmargin is a positive scalar. We are now able to propose a definition for the stability margin:

and with A:

Definition 5.1. The stability margin ofthe system M -~ with the associated size vector k4 defined as above is the largest kmargin such that the system M -~ is stable with ~ extended by the vector k 4 .

Remark: The condition B > 0 means that for 0, the system M-t::. is stable.

[

Consider a system M -~ with ~ a conic sector perturbation (0, I, 1) and with k4 = [kmargin . . . kmargin]. In this case, the definition 5.1 reduces to the definition of km or lip,. Consider a system M -~ with ~ a conic sector perturbation (0, I, I) and with k4 defined as above.

M(jW)] * [Xk(jW)S(jW) Yk{jW)S{jW)] [M(jW)] I Yk{jW)* S(jw) Zk(jW)S(jW) I

kmargin

=

6. COMPUTATIONAL ASPECTS The problem (4) is a quasi-convex optimization problem. However, it is infinite dimensional: the variables are the unknown functions S and G. Two approaches are possible.

110

Ibwl ::; 1 and ';0 = ';n + ';gbe , with Ibel ::; 1 where Wn and ';n are the nominal values of Wo and ';0. By "separating" bw and be from the overall system, we obtain the following representation q(s) = Mqp(s)p(s) and p(s) = Ll(s)q(s) , with Mqp(s) as

• The first one consists in using a frequency gridding (Doyle 1982) . For a fixed frequency w, this problem reduces to a quasi-convex optimization problem over (complex) Linear Matrix Inequalities (LMIs). It is equivalent to a problem over twice as big real LMIs (Boyd and El Ghaoui 1993). Several polynomial time algorithms have been devised for solving these problems, see e.g. (Boyd et al. 1994) . • The second approach consists in restricting the functions S and G to be the sum of N fixed scalings SI and GI : S(jw) = N

0 Wg o 0 [

o

0

€g 0

1+

0

with d(s) =

10.-(.+10)4(.) '4{')(4{') 10)

=

1=1

+ 2';nwns + w; W

.6. =

L:7JIG I(jW),

where 01 and 7Jl are the decision variables. The next step consists in considering the state space representation of the transfer functions M , S and G. Then, the frequency domain inequalities of (4) are equivalent to (real) Linear Matrix Inequalities (for more details see (Ly et al. 1994, Balakrishnan 1995)). One ofthe drawback of this approach is that it leads to the resolution of large Linear Matrix Inequality problems.

and .6.(5) =

-

k

1 .

[~o ~o ~ 1which is simpler than the o

1

60 6.., 0 0 0 [ o 0 6~

0

€o

previous one.

7. NUMERICAL EXAMPLE

s2+~owos +w~

and

Despite the simplicity of this example, the obtained representation is quite intricated. By "separating" ';0 and wo , avoiding the normalization step, as it is allowed by our approach, we obtain the following form:

1=1

J

1

5

S2

N

L:0ISI(jW) and G(jw)

[2{n 51+ Wn

Now, we focus on the following problem: what is the biggest kmarge such that, for the parameters Wo and ';0 in the intervals [8 + 2(1 - kmarge) , 10J and [0.1 + 0.2(1 - kmarge) , 0.3], the closed loop system remains stable ? For kmarge = 1, we get the nominal intervals, whose vertices are represented by a star, figure 2. This particular problem

I--

Fig. 1. Closed loop system We want to investigate the stability of the closedloop system, represented figure 1, for variations of the damping ';0 and the frequency wo . The gain k is set to -1. Applying the well-known RouthHurwitz criterium, this system is proved to be stable if ';0 > 0 and Wo > 0 and ';0 > ~ . The o stability domain, in the parameter plane (wo , ';0), is represented figure 2. We now investigate this

.,.... . .

•." ,

. ..........

... ... ... .

u

.... --..

.....

...

...... ...... :. ........ .

Fig. 3. Plot of k versus frequency w

." .~1.... ;-. .....•.,:..................................... + ........,..........•.........H

_\

can not be directly addressed using classical f.L analysis. In this framework, for a given frequency, several Generalized Eigenvalue problems must be solved. Our approach addresses this problem by solving one Generalized Eigenvalue problem, for a given frequency. To the purpose, the optimization problem, defined theorem 5.1 is coded, using the Generalized Eigenvalue problem solver of the Matlab Tool-box LMI Control Toolbox (Gahinet et al. 1995). The problem is solved, using a frequency

:1

•.• .• \\~ ........•..•.....•..........•.........'-....... .t .. .....J.•.......+ ........-tj

. .rN: ... .:..... ...

"

Fig. 2. Stability domain (x-axis: Wo and y-axis ';0) domain, using the robustness analysis approach. To apply standard f.L analysis, parameters Wo and ';0 need to be scaled: Wo = Wn + Wgbw , with

111

sweeping. For a given value of w, the corresponding k{w) is obtained; kmarge is then estimated by the minimum value of k(w), for the tested set of frequencies. For w = 0, the obtained k(O) is 5 X 105 • For the other tested frequency, k(w) is represented figure 3. The obtained kmarge is 1.49. The corresponding box is plotted figure 2.

Doyle, J.C. (1982). Analysis of feedback systems with structured uncertainties. IEE Proc. 129D(6), 242-250. Fan, M. K. H., A. 1. Tits and J. C. Doyle (1991). Robustness in the presence of mixed parametric uncertainty and unmodeled dynamics. IEEE Trans. Aut. Control 36(1), 25-38. Fan, M. K. H. and A. L. Tits (1992). A measure of worst-case Hoo performance and of largest acceptable uncertainty. Syst. Control Letters 18, 409-42l. Ferreres, G. and V. Fromion (1996) . Robustness analysis using the 11 tool. Proc. of the IEEE CDC 1996 4, 4566-4570. Ferreres, G., G. Scorletti and V. F'romion (1996). Advanced computation of the robustness margin. In: Proc. IEEE Conf. on Decision and Control. pp. ·4~80-4584. Gahinet, P., A. NEmlirovsky, A. L. Laub and M. Chilali (1995): LMI Control Toolbox. The Mathworks Inc. Ly, J .H., M.G. Safonov and R.Y. Chiang (1994). Real / complex multivariable stability margin computation via generalized popov multiplier - LMI approach. In: Proc. American Control Conf. Packard, A. and J.C. Doyle (1993). The complex structured singular value. Automatica 29(1),71-109. Safonov, M. G. (1980). Stability and Robustness of Multivariable Feedback Systems. MIT Press. Cambridge. Safonov, M. G. (1982) . Stability margin of diagonaly perturbed multi variable feedback systems. IEE Proc., Part D 129(6), 251-256. Safonov, M. G. and M. Athans (1981) . A multiloop generalization of the circle criterion for stability margin analysis. IEEE Trans. Aut. Control AC-26(2), 415-422. Scorletti, G. (1997). Robustness analysis with time delays. In: Proc. IEEE Conf. on Decision and Control (IEEE, Ed.). San Diego, California. pp. 3824-3829. Tekawy, J . A., M. G. Safonov and R. Y. Chiang (1992). Convexity property of the one-sided multi variable stability margin. Syst. Control Letters 19(2), 131-137. Wen, J.T. (1988). Robustness amalysis based on passivity. In: Proc. IEEE Conf. on Decision and Control. Willems, J. C. (1969). The Analysis of Feedback Systems. Vol. 62 of Research Monographs. MIT Press. Cambridge, Massachusetts. Zames, G. (1966). On the input-output stability of time-varying nonlinear feedback systemsPart I, IT. IEEE Trans. Aut. Control. Zames, G. (1981). Feedback and optimal sensitivity: Model reference transformations, multiplicative seminorms, and approximate inverses. IEEE Trans. Aut. Control AC26(2),301-320.

8. CONCLUSION In this paper, we propose an extension of the multi-variable stability margin, the one-sided multi-variable stability margin and the skewed I-l margin (see the definition 5.1). The theorem 5.1 establishes a lower bound whose estimation reduces to solve LMI based problems. The main idea is to avoid the normalization step which is always prerequisite to the stability analysis. This allows a far less intricate resolution of the stability problems. An illustration of these advantages is provided in the paper (Ferreres et al. 1996) . Acknowledgments. This paper presents research results of the Belgian Programme on Interuniversity Poles of Attraction initiated by the Belgian State, Prime Minister's Office for Science, Technology and Culture. The scientific responsability rests with the author. The author wishes to thank Vincent Fromion, Stephane Font from, ESE, Paris, Gilles Ferreres from the CERT, Toulouse, Laurent El Ghaoui ENSTA and Gilles Due, ESE, Paris for valuable help and discussions.

9. REFERENCES

Balakrishnan, V. (1995). Linear Matrix Inequalities in robustness analysis with multipliers. Systems and Control Letters. Boyd, S. and L. El Ghaoui (1993). Method of centers for minimizing generalized eigenvalues. Linear Algebra and Applications, special issue on Numerical Linear Algebra Methods in Control, Signals and Systems. Boyd, S., L. El Ghaoui, E. Feron and V. Balakrishnan (1994). Linear Matrix Inequalities in Systems and Control Theory. Vo!. 15 of Studies in App!. Math .. SIAM. Philadelphia. Braatz, R. D., P. M. Young, J. C. Doyle and M. Morari (1994). Computational complexity of I-l calculation. IEEE Trans. Aut. Control AC-39(5), 1000-1002. Doyle, J. C. (1985). Structured uncertainty in control system design. Proc. IEEE Conf. on Decision and Control pp. 260-265. Doyle, J., J. E . Wall and G. Stein (1982). Performance and robustness analysis for structured uncertainties. In: Proc. IEEE Conf. on Decision and Control. pp. 629-636.

112