Technical report from Automatic Control at Linköpings universitet
Employing Kronecker Canonical Form for
H
∞
Synthesis Problems
Anders Helmersson
Division of Automatic Control
E-mail: andersh@isy.liu.se
23rd March 2011
Report no.: LiTH-ISY-R-3007
Submitted to IEEE Transactions of Automatic Control, 2010-10-20
Address:
Department of Electrical Engineering Linköpings universitet
SE-581 83 Linköping, Sweden
WWW: http://www.control.isy.liu.se
AUTOMATIC CONTROL REGLERTEKNIK LINKÖPINGS UNIVERSITET
Technical reports from the Automatic Control group in Linköping are available from http://www.control.isy.liu.se/publications.
Abstract
The Kronecker canonical form (KCF) can be employed when solving H∞
synthesis problem. The KCF structure reveals variables that can be elimi-nated in the semidefinite program that defines the controller. The structure can also be used to remove states in the controller without sacrificing per-formance.
In order to find the KCF structure, we can transform the relevant matrices to a Guptri (Generalized upper triangular) form using orthogonal trans-formations. Thus we can avoid finding the KCF structure explicitly, which is a badly conditioned problem.
Keywords: KCF, GUPTRI, LMI, semidefinite programming, controller synthesis, model reduction
Employing Kronecker Canonical Form for H
∞
Synthesis Problems
Anders Helmersson
2011-03-23
Abstract
The Kronecker canonical form (KCF) can be employed when solving H∞synthesis problem. The KCF structure reveals variables that can be
eliminated in the semidefinite program that defines the controller. The structure can also be used to remove states in the controller without sac-rificing performance.
In order to find the KCF structure, we can transform the relevant ma-trices to a Guptri (Generalized upper triangular) form using orthogonal transformations. Thus we can avoid finding the KCF structure explicitly, which is a badly conditioned problem.
1
Introduction
In the control community, there was a shift in the 1980s from the more tradi-tional Linear Quadratic Gaussian (LQG) control towards H∞ control. Zames’
article from 1981 [15] was an important influence for the development of H∞
synthesis algorithms, initially based on Riccati equations [8, 3]. Soon after, Lin-ear Matrix Inequalities (LMIs) were used to solve the same synthesis problem [4, 9].
In the original algorithms for solving the H∞synthesis problem using Riccati
equation, there exist a conditions on the rank on some of the state-space ma-trices, in particual the D12and D21matrices. The LMI approach, on the other
hand, these conditions are not present, and a controller can be found without modifying the problem as needed in the Riccati approach.
In the general case the controller will have the same number of states as the system it is designed to control. This includes any states added to the plant in order to specify the performance criteria of the closed loop system. As we will see in this paper, the rank deficiencies in D12and D21can actually be employed
to reduce the order of the controller even without sacrificing performance. This paper discusses the structure of the LMI formulation, which can be described using Kronecker Canonical Forms (KCF) [7]. The KCF structure is used to partition the variables in the LMI program into two parts: unbounded and bounded variables. In order to find the minimum H∞ gain, γ, of the
closed-loop system, the unbounded variables can be removed, while the bounded variables must be kept. This can be used to reduced the complexity of the LMI problem and it can also improve the numerical properties. In addition, the KCF
structure can be used to reduce the order of the controller, without reducing its performance.
Finding the KCF of a system is a badly conditioned problem. However, we do not necessarily need the explicit KCF structure. Instead we can use nu-merically well-behaved algorithms like Guptri (Generalized Upper Triangular) decomposition [1, 2].
When unbounded elements are present in the LMI problem, this also indi-cates that the optimal controller may have very high gains (increasing to infinity when γ approaches its minimum value). Also, such a controller may be fragile in the sense that it becomes sensitive to parameter variations [10]. Due to these aspects, we could even say that such a synthesis problem is not properly defined since some properties are missing in the specification.
The paper is organized as follows. We start by making a short recapitulation of the Real Bounded Lemma and the LMI approach to solve H∞ synthesis
problems. We also discuss how rank deficiencies in D21 can be used to reduce
the order of the controller. We then give an overview of the Kronecker Canonical Form and the Guptri form.
Before presenting the main result we use the KCF structure on an example from the COMPleib’s benchmark suite. We also discuss how the controller can
be recovered when infinite variables are present in the LMIs. Finally, we present the results from the COMPleib benchmarks followed by conclusions.
1.1
Notations
For symmetric or Hermitian matrices X ≺ () 0 denotes negative (semi-)definite matrices. In sums, ()T denotes the transpose of the previous term. X∗ de-notes the complex conjugate transpose of a complex matrix and X† denotes the pseudo-inverse of X.
2
LMI formulation
We start by taking a quick review of the LMI formulation used to measure the H∞performance of a system. We then proceed to the LMI synthesis technique
used for finding an optimal H∞ controller. Finally, we discuss how rank
defi-cencies in the D12 and D21 matrices can be used for reducing the order of the
controller.
2.1
The Bounded Real Lemma
Consider the Bounded Real Lemma (BRL) [13, 5], which can be used to deter-mine the H∞ gain, γ, of a continuous-time system, G(s) = C(sI − A)−1B + D,
P A + ATP P B CT BTP −γI DT C D −γI ≺ 0, (1) where P 0.
Using Schur complement the BRL (1) can also be written as 0 0 0 −γI + γ−1 CT DT C D + I 0 P A B + T ≺ 0,
where ()T denotes the transpose of the previous term. This LMI can have
unbounded solutions in P , if there are strictly stable, uncontrollable modes. To see this, let x belong to the null space of
A − λI B ∗
, and let P = P0+αxx∗,
where P0 0 is any solution to (1). If Re λ ≤ 0 then P is also a solution for
any α ≥ 0, since I 0 xx∗ A B + ∗ = 2(Re λ) xx∗ 0 0 0 0.
In addition we must assure that A is strictly Hurwitz (Re λ < 0); otherwise (1) is not feasible.
2.2
Synthesis
We will now consider the system, G, which takes the inputs w and u and pro-duces the outputs z and y:
G := A B1 B2 C1 D11 D12 C2 D21 D22 to be controlled by K := ˆ A Bˆ ˆ C Dˆ ,
which connects to y and u of G, see Figure 1. The closed loop system, G from w to z, is given by G := A B C D = A 0 B1 0 0 0 C1 0 D11 + 0 B2 I 0 0 D12 ˆ A Bˆ ˆ C Dˆ 0 I 0 C2 0 D21 . (2) K G -- -w z u y
Figure 1: The synthesis problem consists of the design task of finding a con-troller K to the system G such that the H∞ norm of the closed-loop system,
from w to z, is minimized.
The Bounded Real Lemma (1) applied to G(s) = C(sI −A)−1B +D, together with the Elimination Lemma [4, 9] can be used to derive the necessary and sufficient conditions for the existence of a controller with a closed-loop H∞
norm less than γ. The conditions are FX(γ, X) := NX 0 0 I T XA + ATX XB 1 C1T BT 1X −γI DT11 C1 D11 −γI NX 0 0 I ≺ 0, (3a) FY(γ, Y ) := NY 0 0 I T AY + Y AT Y CT 1 B1 C1Y −γI D11 B1T D11T −γI NY 0 0 I ≺ 0, (3b) where NX and NY designate any bases of the null spaces of
C2 D21 and
BT
2 DT12 , respectively. In addition X and Y are connected by
X I
I Y
0 (4)
and rank(XY − I) ≤ r, where r is the order of the controller.
When X and Y are found, we proceed by finding a P 0 such that P = X ∗ ∗ ∗ = Y ∗ ∗ ∗ −1 , (5)
which together with the Bounded Real Lemma of the closed loop system gives an LMI condition for the controller.
2.3
Rank deficiencies in D
21In order to illustrate how the structure of the synthesis problem can be used to reduce the order of the controller, we assume that D21 in the augmented plant
is rank deficient. Rank deficencies in the D12 matrix can be handled similarly.
Without loss of generality, we can assume that
C2 D21 has full rank,
which means that the number of outputs to the controller, m2= rank
C2 D21 .
If this is not the case we can reduce the number of outputs to the controller. When D21∈ Rm2×p1 has full row rank, then we can choose
NX= NX1 NX2 = I −D†21C2 ,
which means that NX1∈ Rn×n has full rank. When D21 looses row rank then
also NX1 will do so, which in turn means that some elements in X become
free and will not affect FX(γ, X) in (3a). Note that D21 looses row rank when
p1< m2 (more control outputs than disturbance inputs) or when D21is zero or
has linearly dependent rows.
Suppose that NX1looses rank. We assume that NX2=
0 ∗
, possibly after a similarity transformation. In this case X can be partitioned into
X = X11 X12 X21 X22 , where X11does not affect FX(γ, X).
Using Schur complement, (4) can be rewritten as X − Y−1 0 and the rank constraint becomes rank(X − Y−1) ≤ r. Note that (4) always implies that both X 0 and Y 0 so ¯Y = Y−1 is well defined. Then (4) is equivalent to
X11− ¯Y11 X12− ¯Y12
X21− ¯Y21 X22− ¯Y22
0.
When X22− ¯Y22is nonsingular, this condition is equivalent to
X11− ¯Y11− X12− ¯Y12 X22− ¯Y22 −1
X21− ¯Y21 0,
X22− ¯Y22 0.
Otherwise (when X22− ¯Y22 is singular) we need to replace the inverse by the
pseudo inverse. In this case the null space of X22− ¯Y22belongs to the null space
of X12− ¯Y12.
We now assume that X11 can be chosen freely without affecting FX(γ, X).
Let
X11= ¯Y11+ X12− ¯Y12 X22− ¯Y22 †
X21− ¯Y21 ,
which implies that rank X − Y−1 = rank X22− ¯Y22 ≤ rank NX1.
This means that when D21 looses row rank, we can reduce the number of
states in the controller accordingly (without affecting the performance, γ). For instance if D21 = 0 we can remove m2 states from the controller, from n to
n − m2states. In the state feedback case, when m2= n, the controller becomes
static (state feedback).
This illustrates that the structure of the H∞ problem sometimes can be
employed for both reducing the number of variables in the LMI problem and also for reducing the number of states in the optimal controller. We will generalize this in the sequel.
3
Kronecker Canonical Form (KCF)
We will now take a more general view of the structure of the H∞ synthesis
problem.
3.1
The Structure of the Synthesis LMIs
Using Schur complement we can rewrite (3a) asFX(γ, X) := NXT XA + ATX XB 1 BT 1X −γI + γ−1 CT 1 DT 11 C1 D11 NX≺ 0 (6) This matrix inequality is on the form
Q(γ) + U XVT + V XUT ≺ 0, (7)
with
UT =
A B1 NX and VT = I 0 NX
In the discrete-time case the matrix inequalities will be of the type
Q(γ) + U XUT− V XVT ≺ 0. (8)
We can treat the discrete-time case using similar techniques as the for continuous-time case, but we will not pursue that track here.
3.2
The KCF structure
We will now consider (7) and (4) where we assume that Y is given. Q(γ) + U XVT + V XUT ≺ 0,
X Y−1 0.
In order to solve this type of LMIs we can employ the Kronecker Canonical Form (KCF). The KCF is a generalization of the Jordan form and describes the structure of the generalized eigenvalues and eigenvectors of the pencil U − λV . By applying nonsingular P and Q, we can partition the pencil into a block-diagonal structure:
P−1(U − λV )Q = diag(Ui− λVi),
where each block Ui− λVimust be on one of the following four forms [7]: Jj(α),
Nj, Lj, and LTj. The Jordan j-by-j block is defined by
Jj(µ) := µ − λ 1 . .. . .. . .. 1 µ − λ and Nj:= 1 −λ . .. . .. . .. −λ 1 ∈ Rj×j
corresponds to Jj(∞). The remaining two blocks, Lj ∈ Rj×(j+1) and LTj ∈
R(j+1)×j, are rectangular blocks defined by
Lj:= −λ 1 . .. . .. −λ 1 .
The null space of Lj is given by
1 λ · · · λj T for any value of λ. We illustrate the KCF structure on the following example:
U = 1 2 3 4 5 7 , V = 7 8 9 10 11 13 ,
where we can use P = 6 1 6 0 and Q = −1 −5 −2 5 1 3 −3 3 −1 , to obtain P−1(U − λV )Q = −λ 1 0 0 0 1 − λ = L1⊕ J1(1).
Here we use the notation A ⊕ B = diag[A, B]. Thus, this eigenvalue problem has a finite eigenvalue in 1 and an undefined eigenvalue (which can take any value).
3.3
The Generalized Upper Triangular (Guptri) Form
In order to solve these types of generalized eigenvalue problems we employ algorithms like Guptri (Generalized upper triangular) form [1, 2]. Unitary matrices, P and Q, are used here, which makes the algorithm numerically well-conditioned. For the previous example we can use
P = −0.7071 −0.7071 −0.7071 0.7071 and Q = 0.1690 0.8407 0.5145 −0.8452 −0.1449 0.5145 0.5071 −0.5218 0.6860 to obtain PT(U − λV )Q = 0 1.4349 −9.2164 0 0 4.1231 − λ 1.4343 −0.0410 −23.7685 0 0 4.1231 . The KCF structure can be derived from the staircase block structure of the Guptri form. Due to numerical tolerances in the computations and the finite representation, the elements below the staircase may not be exactly equal to zero but they may differ with some given tolerance. This also means that the Guptri form is not unique in the general case: for a given tolerance we can sometimes find several forms. As the tolerance is increased more forms may become acceptable.
The basic algorithm for computing the Guptri form is based on the singular value decomposition (SVD). For the example, we first apply the SVD on U to find obtain a Q1 such that the first column of U Q1 becomes zero. Next, find a
P1 such that P1TV Q1 becomes zero in its left, bottom element. We have then
found the first step and can proceed to the next step of the staircase structure. When ready we have separated the Lj and Jj(0) blocks.
The LT
j and Njblocks can be extracted using the same algorithm but applied
on VT and UT. Finally, the remaining blocks, which are J
j(µ) blocks, are
obtained by using a standard eigenvalue decomposition. The order of extraction of blocks can be tailored by applying the algorithm on (U, V ), (UT, VT), (V, U ),
4
The ROC7 Example
In order to illustrate how the LMI structure can be used for solving H∞synthesis
problems, we have taken the COMPleib’s benchmarks as examples [11]. We start
by considering the ROC7 benchmark. Here we have removed the augmented state before applying the H∞ synthesis. This means that the reduced ROC7
plant has four states, two outputs to the controller, y, and one control input, u.
A B1 B2 C1 D11 D12 C2 D21 D22 = 0 1 0 0 0 0 −1 0 0 0 1 −0.2 0 0 0 1.02 0 0 0.2 0 0 0 −0.2 1 0.1 0 0 0 0 0 0 0 0.1 0 0 0 0 0 0 0 0 0.2 1 0 0 0 0 0 0 0 1 0 0 0
The KCF structure of FX(γ, X) is L1⊕ N2, which can be written as:
FX(γ, X) := Q(γ) + 0 1 0 0 0 0 0 1 0 0 0 0 X 1 0 0 0 0 0 1 0 0 0 0 1 T + T = Q(γ) + 2x12 x23+ x14 x24 x23+ x14 2x34 x44 x24 x44 0 ≺ 0 (9) where X = x11 x12 x13 x14 x12 x22 x23 x24 x13 x23 x33 x34 x14 x24 x34 x44 Y−1. (10)
Note that FX(γ, X) only has five linearly independent variables: x12, x23+
x14, x24, x34 and x44.
Here we can modify x12 to decrease without bounds and at the same time
increasing x11and x22in order to satisfy (10). This means that x12will dominate
the first row and column in (9). After removing the first row and column we get q22(γ) q23(γ) q23(γ) q33(γ) + 2x34 x44 x44 0 ≺ 0. (11)
In addition (10) can be reduced to
x33 x34
x34 x44
¯Y22, (12)
by removing the two first rows and columns. We can go one step further since both x44 and x34 can grow without bounds. If there exists a solution to (11)
then the lower right element, q33(γ), must be negative, say −ε < 0. For instance,
let x33 x34 x34 x44 = ε 5t3 2t2 2t2 t ε 6 t3 0 0 t ¯Y22
by choosing t large enough. This allows the 3-4 block to dominate both in X and FX.
Using this argument we can get rid of X in (9) altogether. What remains is just a lower bound on γ given by q33(γ) < 0, which here is 0. In this way we
have reduced the number of variables from 10 + 10 + 1 = 21 to 10 + 1 = 11. Next, looking at the structure of the FY(γ, Y ) ≺ 0, which is LT2 ⊕ LT2, or
FY(γ, Y ) := QY(γ) + 0 0 0 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 0 0 0 0 1 Y 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 0 0 0 0 1 0 0 0 0 T + T ≺ 0.
We first look at the upper left 3 × 3 block, which must be negative definite, q11 q12 q13 q12 q22 q23 q13 q23 q33 + 0 y11 y12 y11 2y12 y22 y12 y22 0 ≺ 0.
This block bounds all elements in the upper left block of Y . To see this we first consider y12in the upper right corner. This element is bounded by the first and
last elements on the diagonal as |q13+ y12|2< q11q33. Using the same argument
y12in turn bounds y11and y22. Similarly, the lower right block puts bounds on
y33, y34 and y44. Finally, the elements in the off-diagonal blocks are bounded
by the diagonal blocks. Thus all elements in Y are bounded.
Solving FY(γ, Y ) ≺ 0 with Y 0 yields γ = 1.12033. Using the full set of
LMIs gives a slightly higher value, γ = 1.12037. Elements in X become large (∼ 105) and the solver also gives warnings about problems when finding the
solution.
When recovering a controller that achieves the optimal gain we can observe that the controller has very high gains:
ˆ A Bˆ ˆ C Dˆ ≈ 2858.981 982.834 −41590.878 8451.580 −8656.141 −2970.635 125899.336 −25584.559 1649.430 566.424 −23991.035 4875.370 . Even worse, the controller is fragile, since small changes in its parameters may cause large increase in γ. When the parameters are rounded to three decimals as above we obtain γ = 1.1233; when rounded to one decimal the gain more than doubles, γ = 2.3921, for the closed loop system.
However, if we relax the performance by fixing γ to say 1.13, we obtain a controller with much less gain and higher robustness for parameter variations:
ˆ A Bˆ ˆ C Dˆ ≈ −1.3 −12.0 −4.7 −1.3 1.8 −15.7 −12.6 −0.1 −1.3 −18.2 −7.7 −1.8 .
The conclusion of this elaborated example is that depending on the KCF struc-ture it is sometimes possible to remove elements from X or Y when solving the LMIs during the synthesis step. For KCF blocks like Lj, Jj(0) and Nj the
kinds of KCF blocks are present, the order of the controller can also be reduced without sacrificing performance.
However, when the corresponding controller is recovered it may become frag-ile (sensitive to parameter variations) and may have unrealistically high gains when the γ approaches its optimal value. We can even argue that the synthesis problem is not well formulated since there are no bounds on the controller gains or on its sensitivity to parameter variations.
5
KCF Partitioning
We will now take a more general look at how the various types of KCF blocks affect the solution of the LMIs. The matrices U and V are on a block diagonal KCF structure and the Q and X = X∗matrices have a corresponding full block structure. Note the use of complex conjugate transpose ()∗, since the eigenvalues and eigenvectors may be complex.
FX(γ, X) = Q(γ) + U XV∗+ V XU∗= Q11 Q12 Q∗ 12 Q22 + U1X11V1∗ U1X12V2∗ U2X12∗ V1∗ U2X22V2∗ + ∗ = Q11+ U1X11V1∗+ V1X11U1∗ Q12+ U1X12V2∗+ V1X12U2∗ Q∗12+ U2X12∗V1∗+ V2X12∗U1∗ Q22+ U2X22V2∗+ V2X22U2∗
5.1
L
jblocks
We first consider the case when (U1, V1) =
0 I , I 0 , which corre-sponds to Lj. The 1-1-block is dominated by
U1X11V1∗+V1X11U1∗= x12+ x21 x13+ x22 . . . x1n+ x2(n−1) x22+ x31 x23+ x32 . . . x2n+ x3(n−1) .. . ... . .. ... xn1+ x(n−1)2 xn2+ x(n−1)3 . . . xn(n−1)+ x(n−1)n . (13) To see this we first let x21to go to −∞ and at the same time letting x11approach
∞ in order to satisfy X Y−1. Then the first row and column will dominate.
We continue with the second row and column in the same way until the final row and column. Thus only
Q22(γ) + U2X22V2∗+ V2X22U2∗≺ 0 (14)
remains.
5.2
J
jand N
jblocks
We next consider the case when (U1, V1) = (µI + T, I), which corresponds to
Jj(µ), where T = U1− µV1= 0 1 . .. . .. . .. 1 0 .
Then FX(γ, X) = Q11+ 2(Re µ)X11+ T X11+ X11T∗ Q12+ X12(U2+ ¯µV2)∗+ T X12V2∗ Q∗12+ (U2+ ¯µV2)X12∗ + V2X12∗T∗ Q22+ U2X22V2∗+ V2X22U2∗ . If Re µ < 0 then X11 = X11∗ dominates in the first j rows and columns. The
remaining condition in FX(γ, X) ≺ 0 becomes (14).
Next, assume that Re µ = 0 by letting µ = jω where ω ∈ R. Then
T X11+ X11T∗= x12+ x21 . . . x1n+ x2(n−1) x2n .. . . .. ... ... xn1+ x(n−1)2 . . . xn(n−1)+ x(n−1)n xnn xn2 . . . xnn 0 .
Here we can identify the framed block as (13), which according to previous discussion, can be made arbitrarily negative definite. Thus it will dominate and only the last row and column remain. Consequently, the first j − 1 rows and columns can be removed and the term T X12V2∗disappears altogether, since the
last row of T is zero.
We can now use the Elimination lemma [4, 9] when solving for the last row of X12. Let N = diag[N1, N2] be defined as a matrix spanning the null space of
(U − jωV )∗=
(U1− jωV1)∗ (U2− jωV2)∗ . The two equivalent conditions
of the elimination lemma become N∗Q(γ)N ≺ 0 and (14), where we have used the fact that
N2∗(U2X22V2∗+ V2X22U2∗) N2
= N2∗(jωV2X22V2∗− jωV2X22V2∗) N2= 0,
since (U2− jωV2)∗N2 = 0. Note that the condition N∗Q(γ)N ≺ 0 does not
involve X and can be replaced by a lower bound on γ. When Re µ > 0, the 1-1-block becomes
Q11+ 2(Re µ)X11+ T X11+ X11T∗≺ 0, or Q11+ (2 Re µ) | {z } >0 x11 . . . x1(n−1) x1n .. . . .. ... ... x(n−1)1 . . . x(n−1)(n−1) x(n−1)n xn1 . . . xn(n−1) xnn + x12+ x21 . . . x1n+ x2(n−1) x2n .. . . .. ... ... xn1+ x(n−1)2 . . . xn(n−1)+ x(n−1)n xnn xn2 . . . xnn 0 ≺ 0.
Looking at the rightmost, bottom element, we conclude that xnnis bounded by
Q11 and X Y−1 0.
Next, looking at the second last element along the diagonal, the elements x(n−1)(n−1)and xn(n−1)are bounded, since |xn(n−1)|2< xnnx(n−1)(n−1)(X 0)
and q(n−1)(n−1) + (2 Re µ)x(n−1)(n−1) < 2|xn(n−1|. We can now proceed by
taking x(n−2)(n−2) and so on until x11. We conlude that all elements in X11
5.3
L
Tj
blocks
We next look at the last type of KCF block: LT
j for which U1 = 0 I and V1= I 0 : Q11+ 0 x11 . . . x1(n−1) x1n x11 x12+ x21 . . . x1n+ x2(n−1) x2n .. . ... . .. ... ... x(n−1)1 x(n−1)2+ xn1 . . . x(n−1)n+ xn(n−1) xnn xn1 xn2 . . . xnn 0 ≺ 0 together with X11 0.
Let kI Q11 be an upper bound of Q11 and use Re x ≤ |x|. Considering
the elements on the diagonal and sub and super-diagonals we get x211< k (k + 2|x12|) |x12|2< x11x22 x222< (k + 2|x12|)(k + 2|x23|) . . . |x(n−1)n|2< x(n−1)(n−1)xnn x2nn< k k + 2|x(n−1)n| .
Replacing the xi(i+1)terms with their upper bounds gives
x211< k (k + 2√x11x22) x222< (k + 2√x11x22) (k + 2 √ x22x33) . . . x2(n−1)(n−1)< k + 2√x(n−2)(n−2)x(n−1)(n−1) k + 2√x(n−1)(n−1)xnn x2nn< k k + 2√x(n−1)(n−1)xnn .
Next, choose a κ such that κ√xiix(i+1)(i+1)≥ k holds for all i = 1, . . . , n − 1:
x411< k2(κ + 2)2x11x22
x422< (κ + 2)4x11x222x33
. . .
x4(n−1)(n−1)< (κ + 2)4x(n−2)(n−2)x2(n−1)(n−1)xnn
x4nn< k2(κ + 2)2x(n−1)(n−1)xnn.
Further we can show that
x311< k2(κ + 2)2x22
x522< k2(κ + 2)4×22−2x333 . . .
x2n−1(n−1)(n−1)< k2(κ + 2)4(n−1)2−2x2n−3nn x3nn< k2(κ + 2)2x(n−1)(n−1).
Table 1: Classification of KCF structure.
Lj unbounded
Jj(µ), Re µ < 0 unbounded
Nj ∼ Jj(∞) mixed, elimination lemma
Jj(µ), Re µ = 0 mixed, elimination lemma
Jj(µ), Re µ > 0 bounded
LT
j bounded
Using the two last inequalities, we get xnn < k(κ + 2)n−1 and x(n−1)(n−1) <
k(κ+2)3n−5, which also means that all x
ii< k(κ+2)(2i−1)(n+1−i)−iare bounded
(using κ =√2 − 1 can be used for all n ≥ 2). Further, all off-diagonal xij
ele-ments are bounded by the xii elements. Thus, all elements in X11are bounded.
5.4
Summary
We summarize this section by stating the following lemma.
Lemma 1. Consider the problem of minimizing γ with respect to Q(γ)+U XVT+
V XUT < 0 and X > 0. Partition (U, V ) into KCF blocks and partition Q and
X accordingly.
(i) Blocks corresponding to Lj and Jj(µ) where Re µ < 0, can be truncated.
(ii) Blocks of the type Nj and Jj(µ) where Re µ = 0, can be truncated and
replaced by a additional condition N∗Q(γ)N < 0, where N designate any basis of the null space of (U − µV )∗.
(iii) The remaining blocks in X corresponding to LT
j and Jj(µ) where Re µ > 0
are bounded.
For some examples in our test suite, both the FX and FY inequalities can
be removed and replaced by a lower bound on γ. These are AC1, AC12, DLR1, UWV, PAS and ROC5.
When solving the LMIs (FX and FY), we start by finding the KCF structure
using the Guptri algorithm. For finding the optimal γ we sort the blocks in the following order: Lj, Jj(0), Nj, Jj(µ), LTj. The Jj(µ) blocks are sorted so
that the eigenvectors corresponding to Re µ ≤ 0 come first. After the optimal γ is obtained, it can be relaxed to a slightly higher value and the controller can be computed. The blocks are then sorted in the order Lj, Jj(µ), Jj(0), Nj, LTj.
The Jj(µ) blocks are sorted so that the eigenvectors corresponding to Re µ ≤ −ε
come first. These, together with the Lj blocks, are truncated in the LMI while
the others are kept. The truncated blocks will correspond to uncontrollable or unobservable modes in the closed-loop system. The parameter ε should selected to get sufficient closed-loop decay rate.
6
Recovering X, Y and the controller, K
After finding X and Y we need to compute the controller, K. This can be done in two steps: first, find a P that satisfies (5) and, secondly, compute a controller such that (1) holds for the closed-loop system.
When X and Y are bounded, we can perform a similarity transformation of the system so that X = Y = Σ I becomes diagonal. Then we can use
P = Σ Γ Γ Σ = Σ −Γ −Γ Σ −1 , where Γ2= Σ2− I 0. When Σ = diag[Σ
r, In−r] contains ones in its diagonal,
the corresponding states can be removed in the controller by using
P = Σr 0 Γr 0 In−r 0 Γr 0 Σr .
If there are unbounded elements in X and Y then this scheme needs to be modified accordingly. The unbounded elements correspond to uncontrol-lable and unobservable modes in the closed-loop system. The closed-loop state space can be divided into four subspaces corresponding to uncontrollable and unobservable modes, uncontrollable modes, unobservable modes, and both con-trollable and observable modes:
P = ∞ 0 0 0 0 ∞ 0 0 0 0 ∞ 0 0 0 0 ∞ 0 0 0 0 I 0 0 0 0 I 0 0 0 0 Σ 0 0 0 0 Γ 0 0 0 0 I 0 0 0 0 ∞ 0 0 0 0 ∞ 0 0 0 0 ∞ 0 0 0 0 ∞ 0 0 0 0 I 0 0 0 0 I 0 0 0 0 Γ 0 0 0 0 Σ = ∞ 0 0 0 0 −∞ 0 0 0 0 I 0 0 0 0 −I 0 0 0 0 ∞ 0 0 0 0 −∞ 0 0 0 0 Σ 0 0 0 0 Γ 0 0 0 0 I 0 0 0 0 −∞ 0 0 0 0 ∞ 0 0 0 0 −I 0 0 0 0 I 0 0 0 0 −∞ 0 0 0 0 ∞ 0 0 0 0 −Γ 0 0 0 0 Σ −1 .
There are still some questions that remain: how do we recover the removed variables in an efficient way and can we even find a reduced order controller in the same process? There even may be no need to recover the complete X and Y , and we can possibly find the controller parameters anyway.
Even if the LMIs don’t bound some elements in X and Y this does not necessarily mean the they should approach infinity. Often, we can find bounded solutions and in that case X and Y can be fully reconstructed and K can be searched for. This is also needed when Jj(µ) where Re µ = 0 are involved since
P in (1) cannot be unbounded for these elements, since this would imply an uncontrollable mode on the imaginary axis.
On the other hand, we can always let the unbounded elements in X and Y approach infinity. Then they will correspond to unbounded KCF blocks (Lj
and Jj(µ) where Re µ < 0) in FX(γ, X) or FY(γ, Y ).
If the generalized eigenvalue, Re λ ≤ 0, corresponds to an unbounded element in FX(γ, X) < 0, then there exists a vector xλ such that
xTλ
A − λI B1 NX = 0. (15)
Suppose there exists an uncontrollable mode in the closed-loop system (2). Then, there exist x 6= 0 and λ such that
xT
A − λI B = 0.
Note that the eigenvalues corresponding to the uncontrollable modes must be defined since there cannot exist LT
Expanding this we get xT A − λI B = x1 x2 T A − λI 0 B1 0 −λI 0 + 0 B2 I 0 ˆ A Bˆ ˆ C Dˆ 0 I 0 C2 0 D21 = xT1(A − λI) −xT2λI xT1B1 + xT2 xT1B2 ˆ A Bˆ ˆ C Dˆ 0 I 0 C2 0 D21 = 0,
This has a solution, K(λ) = ˆC(λI − ˆA)−1B + ˆˆ D, if and only if xT1 A − λI B1 NX = 0 NT ZxT1 A − λI B1 = 0
where NZ designates any base of the null space of
BT 2x1 x2 .
Note that NZ only exists if both B2Tx1 and x2 are zero; otherwise the
second condition disappears. When NZ exists the second condition becomes
xT 1
A − λI B1 = 0, which makes the first condition redundant. This case
corresponds to an uncontrollable mode in λ for the original system described by the pair (A,
B1 B2 ). If we assume that the original system is controllable,
NZ disappears and only the first condition remains, which is equivalent to (15).
What happens if Re λ = 0? In this case we can choose the corresponding element in P arbitrarily large, which drives λ towards the imaginary axis. How-ever, we must choose a finite P in order to obtain a stable closed loop system. Thus, we need to treat the case when Re µ = 0 specially, compared to Jj(µ)
blocks with Re µ < 0.
In some cases the performance bounds obtained using the reduced LMIs are theoretical lower bounds. For instance for AC1 the lower bound is zero, but this would require a controller with infinite gain. In order to obtain a practical, finite-gain controller, we must settle for some γ > 0, which gives reasonable gains in the controller. This can also be interpreted as if the problem is not well formulated and more constraints should be added to limit the gain of the controller.
6.1
Rank Reduction
The KCF structure reveals if we can reduce the order of the controller by choos-ing the free elements in X or Y appropriately. In general, we can reduce the rank of XY − I by one for each block of the type Lj. In addition we can further
reduce the rank by one if there exists a block of the type Jj(0) or Nj. If there
exist several such blocks they can only be combined if the eigenvalues are equal. For instance, for the KCF structure defined by L1⊕ L2⊕ J1(0) ⊕ N2 the rank
can be reduced by 3 (from 8 to 5); the J1(0) block cannot be combined with
the N2 block, but it can be combined with both L1 and L2. Another example
is DLR1, in which FX has the structure J1(0) ⊕ N1⊕ 8J1(−). The controller
can be reduced by 1 from 10 to 9 according to the rule above.
For some block structures we can reduce the rank even more. In many cases the possibility to reduce the rank depends on the values of the fixed elements (in
Table 2: Jet engine (JE) and reactor (REA) models Case n X Y r γ JE1 30 5N1, 25LT1 5LT6 25 3.8500 (3.8702) JE2 21 3N1, 18LT1 15LT1, 3LT2 18 42.348 (42.508) JE3 24 24J1(−) 6LT4 24 2.8833 (2.8860) REA1 4 3N1, LT1 2LT2 1 0.86172 REA2 4 2N1, 2LT1 2L T 2 2 1.1341 REA3 12 3N1, 9LT1 10L T 1, L T 2 9 74.251
X − Y−1or Y − X−1). The conditions are in general non-convex and algebraic. For instance, if x1, x2 and x3 are free variables and a, b and c are fixed in
x1 a b − x2 a b + x2 c b − x2 c x3 0,
we can obtain rank 1, either when b2≥ ac 6= 0, when a = c = 0, or when b > 0
and ac = 0. This example corresponds to an L2block.
7
Examples
In order to illustrate the method using KCF, we have taken the COMPleib’s
benchmarks as examples [11]. In COMPleib the ROC examples (all but ROC3)
have been augmented with internal states (1, 2 or 3), which corresponds to the minimum order of a stabilizing controller. These extra states have been removed before the controller synthesis. If they are kept addition KCF blocks of the type L0 appear, one for each added state.
In the Tables 2-10, the order of the original system and the optimal controller are denoted by n and r, respectively. The H∞ norm of the closed loop system
is given by γ. The existence of an optimal controller has been verified for most systems, but not all (unverified systems are marked by * in the γ-column). For some systems the gain can be made arbitrarily small, denoted by → 0, but as γ approaches zero the controller gain is likely to increase to infinity. If a γ-value is given within paranthesis, this gives the best performance obtained for the closed loop system using the synthesized controller.
The LMIs have been solved in Matlab using the Yalmip interpreter [12] together with SeDuMi 1.3 [14] and LMILAB [6] solvers.
8
Conclusions
In many control synthesis problems, the structure of the plant augmented with weighting functions can be employed to reduce the complexity of the LMI prob-lem and to reduce the number of states in the controller. There are still some details to be elaborated, such that recovering the removed variables in the LMIs in an efficient way.
On the other hand, it can also be argued that if the LMIs are associated with unbounded elements in X or Y (that is if Lj and some Jj blocks are
Table 3: 2D heat flow problems Case n X Y r γ HF2D10 5 3N1, 2LT1 5LT1 2 79536 HF2D11 5 3N1, 2LT1 5L T 1 2 76443 HF2D12 5 4N1, LT1 5L T 1 1 1037666 HF2D13 5 4N1, LT1 5L T 1 1 101548 HF2D14 5 4N1, LT1 5LT1 1 526354 HF2D15 5 4N1, LT1 5LT1 1 173042 HF2D16 5 4N1, LT1 5LT1 1 444144 HF2D17 5 4N1, LT1 5LT1 1 300236 HF2D18 5 2N1, 3LT1 5LT1 3 93.153 HF2D_CD4 7 2N1, 5LT1 7LT1 5 24.458 HF2D_CD5 7 2N1, 5LT1 7LT1 5 25.771 HF2D_CD6 7 2N1, 5LT1 7LT1 5 15.922 HF2D_IS5 5 4N1, LT1 5LT1 1 1131511 HF2D_IS6 5 4N1, LT1 5LT1 1 132263 HF2D_IS7 5 2N1, 3LT1 5LT1 3 36.550 HF2D_IS8 5 2N1, 3LT1 5LT1 3 91.297
Table 4: Decentralized interconnected systems (DIS) and Euler-Bernoulli beams (EB) Case n X Y r γ DIS1 8 L0, N1, 2L2 LT1, 2LT2, LT3 4 4.1593 DIS2 3 2N1, LT1 LT3 1 0.94741 DIS3 6 4N1, 2LT1 LT1, LT5 2 1.0424 DIS4 6 6N1 2LT3 0 0.73147 DIS5 4 LT 4 LT4 3 666.27 EB1 10 LT10 LT10 9 3.1041 EB2 10 LT10 LT10 9 1.7676 EB3 10 LT10 LT10 9 1.7976 EB4 20 LT20 LT20 19 1.7973 EB5 40 LT40 LT40 39 1.7973
Table 5: Helicopter problems
Case n X Y r γ HE1 4 N1, LT3 J1(−), 3J1(+) 3 0.073615 HE2 4 2N1, 2LT1 2LT2 2 2.4181 HE3 8 2L0, 3L1 4LT1, 2LT2 3 0.79897 HE4 8 6N1, 2LT1 8LT1, 2LT2 2 22.838 HE5 8 LT 8 6J1(−), 2J1(+) 7 1.7650 HE6 20 18J1(−), 2J1(+) 4LT1, 8LT2 19 2.3876 HE7 20 J1(−), 2LT6, L7T 4LT1, 8LT2 19 2.6130
Table 6: Aircraft models Case n X Y r γ AC1 5 2N1, N2, J1(−) L4 2 → 0 AC2 5 2N1, N2, J1(−) J1(0), 2LT2 2 0.11149 AC3 5 4N1, LT1 LT1, 2LT2 1 2.9674 AC4 4 3J1(−), J1(+) LT4 3 0.55729 AC5 4 2N1, 2LT1 2LT2 2 657.06 AC6 7 4N1, 3LT1 3L T 1, 2L T 2 3 3.4273 AC7 9 N1, LT0, L T 8 7J1(−), 2J1(+) 8 0.037862 (0.037865) AC8 9 LT0, LT1, LT2, 2LT3 LT9 8 1.6165 AC9 10 2LT1, LT2, 2LT3 L3, L4, N1 7 0.99998 AC10 55 N1, N2, L51T 2LT18, LT19 53 3.2324 * AC11 5 4N1, LT1 LT1, 2LT2 1 2.8079 AC12 4 L1, 2N1 2L1 1 → 0 (0.13) AC13 28 4N1, 24LT1 22LT1, 3LT2 24 156.23 AC14 40 39J1(−), J1(+) 2LT2, LT3, 2LT6, 3LT7 39 100.00 AC15 4 3N1, LT1 4LT1 1 14.863 AC16 4 4N1 4LT1 0 14.856 AC17 4 2N1, 2LT1 2LT1, LT2 2 6.6124 AC18 10 2N1, LT8 2LT3, LT4 8 5.3839 (5.3991)
Table 7: Second order problems
Case n X Y r γ CM1 20 L19 LT1, LT19 19 0.81650 CM2 60 L59 LT1, LT59 59 0.81650 CM1_IS 20 L19 LT1, LT19 19 0.81650 CM2_IS 60 L59 LT1, LT59 59 0.81650 (0.81817) TMD 6 3L0, J2(0), J1(−) 2J2(0), LT2 2 2.1207 FS 5 3N1, 2LT1 3LT1, LT2 2 77919 DLR1 10 J1(0), N1, 8J1(−) 10J1(−) 9 0.061930 DLR2 40 L37, N2 40LT1 38 82.679 * DLR3 40 L38, N1 40LT1 38 82.679 * LAH 48 48J1(−) LT0, LT48 47 5.3727 × 10−5
Table 8: Problems for various physical systems Case n X Y r γ TG1 10 2N1, 8LT1 6LT1, 2LT2 8 3.4652 AGS 12 2N1, 10LT1 8L T 1, 2L T 2 10 8.1732 WEC1 10 4N1, 6LT1 4L T 1, 3L T 2 6 3.6363 WEC2 10 4N1, 6LT1 4L T 1, 3L T 2 6 3.5981 WEC3 10 4N1, 6LT1 4LT1, 3LT2 6 3.7685 BDT1 11 L0, L4, N5 LT1, 2LT5 8 0.26621 MFP 4 2N1, 2LT1 LT4 2 4.1865 UWV 8 N2, N6 L7 6 → 0 IH 21 10N1, 11LT1 21J1(−) 11 → 0 CSE1 20 9L1, N1, J1(−) J1(−), LT1, 9LT2 10 0.019871 CSE2 60 29L1, N1, J1(−) J1(−), LT1, 29LT2 30 0.019996 (0.02445) PAS 5 L0, N1, N2, J1(−) N4, J1(−) 2 → 0 (0.000765) TF1 7 2L1, L2 LT3, LT4 4 0.24636 TF2 7 L0, L1, J1(0), 3J1(−) LT3, LT4 4 5200 TF3 7 L0, L3, 2J1(−) LT3, LT4 5 0.24636 PSM 7 L0, 2N3 LT1, 2LT3 4 0.92024
Table 9: Academic test problems
Case n X Y r γ NN1 3 2N1, LT1 L T 1, L T 2 1 13.130 NN2 2 N1, LT1 LT2 1 1.7644 NN3 4 N2, 2J1(+) N2, J1(−), J1(+) 3 7.9019 (7.97) NN4 4 3N1, LT1 LT1, LT4 1 1.2862 NN5 7 2N1, 5LT1 5LT1, LT2 5 237.87 NN6 9 4N1, 5LT1 7LT1, LT2 5 126.58 NN7 9 4N1, LT5 LT1, LT8 5 33.039 (33.046) NN8 3 2N1, LT1 LT3 1 2.3573 NN9 5 N1, 2J1(−), 2J1(+) L0, N1, 2LT1, J1(+) 3 13.639 (13.641) NN10 8 L1, 2L2, 3LT0 N1, L2, L3, LT0 5 → 0 NN11 16 3N1, 13J1(−) N1, N3, 10J1(−), J1(+) 13 0.011994 (0.01760) NN12 6 2N1, 4LT1 2LT1, 3LT2 4 6.2830 NN13 6 LT 6 LT6 5 10.184 NN14 6 LT 6 LT6 5 9.4314 NN15 3 N1, L1 LT1, LT2 1 0.097741 NN16 8 4N1, 4LT1 2J1(+), 6J1(±j) 4 0.95556 NN17 3 N1, J1(−), J1(+) N1, J1(−), J1(+) 2 2.6364
Table 10: Reduced order controller models Case n X Y r γ ROC1 8+1 LT8 LT8 7 1.1267 ROC2 9+1 2LT4, N1 8J1(−), J1(+) 8 0.041030 (0.041032) ROC3 11 4N1, 7LT1 4L1T, 2LT2, LT3 7 42.647 ROC4 8+1 3J1(−), 2J1(+), N3, L0T 7J1(−), N1, LT0 7 7.6111 (12.56) ROC5 6+1 L0, 2N2, N1 2N2, 2J1(−) 2 → 0 (1 × 10−6) ROC6 3+2 LT 1, LT2 LT1, LT2 2 21.528 ROC7 4+1 L1, N2 2LT2 2 1.1203 ROC8 6+3 6J1(±j) 6LT1 5 3.4870 ROC9 4+2 4J1(±j) 4LT1 3 2.2361 ROC10 5+1 2L1, N1, LT0 N2, J1(−), J1(+), LT1 2 0.071440 (0.07188)
present), this would indicate that the synthesis problem is badly formulated since constraints limiting the controller gain could be missing. If so, the KCF structure will readily provide a tool for diagnosing this kind of problems.
References
[1] J. Demmel and B. Kågström. The generalized Schur decomposition of an arbitrary pencil A − λB: Robust software with error bounds and applica-tions. Part I: Theory and algorithms. ACM Transactions of Mathematical Software, 19(2):160–174, 1993.
[2] J. Demmel and B. Kågström. The generalized Schur decomposition of an arbitrary pencil A − λB: Robust software with error bounds and applica-tions. Part II: Software and applicaapplica-tions. ACM Transactions of Mathemat-ical Software, 19(2):175–201, 1993.
[3] J. C. Doyle, K. Glover, P. Khargonekar, and B. A. Francis. State-space solutions to the standard H2and H∞control problems. IEEE Transactions
on Automatic Control, 34(8):831–847, August 1989.
[4] P. Gahinet and P. Apkarian. A linear matrix inequality approach to H∞
control. International Journal of Robust and Nonlinear Control, 1:656–661, October 1992.
[5] P. Gahinet and P. Apkarian. An LMI-based parametrization of all H∞
controllers with applications. In IEEE Proceedings of the 32nd Conference on Decision and Control, volume 1, pages 656–661, San Antonio, Texas, December 1993.
[6] P. Gahinet and A. Nemirovsikii. LMI-Lab – A Package for Manipulating and Solving LMIs, version 2.0, 1993.
[7] F. R. Gantmacher. The Theory of Matrices, volume II. New York: Chelsea Publishing Company, 1959, 1971.
[8] K. Glover and J. C. Doyle. State-space formulae for all stabilizing con-trollers that satisfy an H∞ norm bound and relations to risk sensitivity.
Systems and Control Letters, 11:167–172, 1988.
[9] T. Iwasaki and R. E. Skelton. All controllers for the general H∞ control
problem: LMI existence conditions and state space formulas. Automatica, 30:1307–1317, August 1994.
[10] L. H. Keel and S. P. Bhattacharyya. Robust, fragile or optimal. IEEE Transactions on Automatic Control, 42(8):1098–1105, August 1997. [11] F. Leibfritz. Compleib: Constraint matrix-optimization problem library
– a collection of test examples for nonlinear semidefinite programs, con-trol system design and related problems. Technical report, Department of Mathematics, University of Trier, Trier, Germany, 2004.
[12] J. Löfberg. YALMIP : A toolbox for modeling and optimization in MAT-LAB. In Proceedings of the CACSD Conference, Taipei, Taiwan, 2004. [13] C. Scherer. The Riccati Inequality and State-Space H∞-Optimal Control.
PhD thesis, Universitat Wurtzburg, Wurtzburg, Germany, 1990. [14] Jos F. Sturm. SeDuMi, 2003.
[15] G. Zames. Feedback and optimal sensitivity: Model reference transforma-tions, multiplicative seminorms, and approximative inverses. IEEE Trans-actions on Automatic Control, 26(2):301–320, February 1981.
Avdelning, Institution Division, Department
Division of Automatic Control Department of Electrical Engineering
Datum Date 2011-03-23 Språk Language 2 Svenska/Swedish 2 Engelska/English 2 Rapporttyp Report category 2 Licentiatavhandling 2 Examensarbete 2 C-uppsats 2 D-uppsats 2 Övrig rapport 2
URL för elektronisk version http://www.control.isy.liu.se
ISBN — ISRN
—
Serietitel och serienummer Title of series, numbering
ISSN 1400-3902
LiTH-ISY-R-3007
Titel Title
Employing Kronecker Canonical Form for H∞Synthesis Problems
Författare Author
Anders Helmersson
Sammanfattning Abstract
The Kronecker canonical form (KCF) can be employed when solving H∞synthesis problem.
The KCF structure reveals variables that can be eliminated in the semidefinite program that defines the controller. The structure can also be used to remove states in the controller without sacrificing performance.
In order to find the KCF structure, we can transform the relevant matrices to a Guptri (Generalized upper triangular) form using orthogonal transformations. Thus we can avoid finding the KCF structure explicitly, which is a badly conditioned problem.