## Advanced LMI Techniques

This last section gives a few hints for making the most out of the LMI Lab. It is directed toward users who are comfortable with the basics, as described in Tools for Specifying and Solving LMIs.

### Structured Matrix Variables

Fairly complex matrix variable structures and interdependencies can be specified
with `lmivar`

. Recall that the symmetric
block-diagonal or rectangular structures are covered by Types 1 and 2 of `lmivar`

provided that the matrix
variables are independent. To describe more complex structures or correlations
between variables, you must use Type 3 and specify each entry of the matrix
variables directly in terms of the free scalar variables of the problem (the
so-called decision variables).

With Type 3, each entry is specified as either 0 or
±*x*_{n} where
*x*_{n} is the *n*-th
decision variable. The following examples illustrate how to specify nontrivial
matrix variable structures with `lmivar`

. The following examples
show variable structures with uncorrelated and interdependent matrix
variables.

#### Specify Matrix Variable Structures

Suppose that the variables of the problem include a 3-by-3 symmetric matrix *X* and a 3-by-3 symmetric Toeplitz matrix, *Y*, given by:

$$Y=\left(\begin{array}{ccc}{y}_{1}& {y}_{2}& {y}_{3}\\ {y}_{2}& {y}_{1}& {y}_{2}\\ {y}_{3}& {y}_{2}& {y}_{1}\end{array}\right).$$

The variable *Y* has three independent entries, and thus involves three decision variables. Since *Y* is independent of *X*, label these decision variables *n *+ 1, *n *+ 2, and *n* + 3, where *n* is the number of decision variables involved in *X*. To retrieve this number, define the Type 1 variable *X*.

setlmis([]) [X,n] = lmivar(1,[3 1]); n

n = 6

The second output argument `n`

gives the total number of decision variables used so far, which in this case is `n`

= 6. Given this number, you can define *Y. *

Y = lmivar(3,n+[1 2 3;2 1 2;3 2 1]);

An equivalent expression to define Y uses the MATLAB® command `toeplitz`

to generate the matrix.

Y = lmivar(3,toeplitz(n+[1 2 3]));

To confirm the variables, visualize the decision variable distributions in *X* and *Y* using `decinfo`

.

lmis = getlmis; decinfo(lmis,X)

`ans = `*3×3*
1 2 4
2 3 5
4 5 6

decinfo(lmis,Y)

`ans = `*3×3*
7 8 9
8 7 8
9 8 7

#### Specify Interdependent Matrix Variables

Consider three matrix variables, *X*, *Y*, and *Z*, with the following structure.

$$X=\left(\begin{array}{cc}x& 0\\ 0& y\end{array}\right),\phantom{\rule{1em}{0ex}}Y=\left(\begin{array}{cc}z& 0\\ 0& t\end{array}\right),\phantom{\rule{1em}{0ex}}Z=\left(\begin{array}{cc}0& -x\\ -t& 0\end{array}\right),$$

where *x*, *y*, *z*, and *t* are independent scalar variables. To specify such a triple, first define the two independent variables, *X* and *Y*, which are both Type 1.

setlmis([]); [X,n,sX] = lmivar(1,[1 0;1 0]); [Y,n,sY] = lmivar(1,[1 0;1 0]);

The third output of lmivar gives the entry-wise dependence of *X* and *Y* on the decision variables $\left({\mathit{x}}_{1},{\mathit{x}}_{2},{\mathit{x}}_{3},{\mathit{x}}_{4}\right):=\text{\hspace{0.17em}}\left(\mathit{x},\mathit{y},\mathit{z},\mathit{t}\right).$

sX

`sX = `*2×2*
1 0
0 2

sY

`sY = `*2×2*
3 0
0 4

Using lmivar, you can now specify the structure of the Type 3 variable *Z* in terms of the decision variables ${\mathit{x}}_{1}=\mathit{x}$ and ${\mathit{x}}_{4}=\mathit{t}$.

[Z,n,sZ] = lmivar(3,[0 -sX(1,1);-sY(2,2) 0]);

Because `sX(1,1) `

refers to ${\mathit{x}}_{1}$ and `sY(2,2)`

refers to ${\mathit{x}}_{4}$, this expression defines the variable `Z`

as:

$$Z=\left(\begin{array}{cc}0& -{x}_{1}\\ -{x}_{4}& 0\end{array}\right)=\left(\begin{array}{cc}0& -x\\ -t& 0\end{array}\right).$$

Confirm this result by checking the entry-wise dependence of `Z`

on its decision variables.

sZ

`sZ = `*2×2*
0 -1
-4 0

### Complex-Valued LMIs

The LMI solvers are written for real-valued matrices and cannot directly handle
LMI problems involving complex-valued matrices. However, complex-valued LMIs can be
turned into real-valued LMIs by observing that a complex Hermitian matrix
*L*(*x*) satisfies

*L*(*x*) < 0

if and only if

$$\left(\begin{array}{cc}\mathrm{Re}\left(L\left(x\right)\right)& \mathrm{Im}\left(L\left(x\right)\right)\\ -\mathrm{Im}\left(L\left(x\right)\right)& \mathrm{Re}\left(L\left(x\right)\right)\end{array}\right)<0.$$

This suggests the following systematic procedure for turning complex LMIs into real ones:

Decompose every complex matrix variable

*X*as*X*=*X*_{1}+*jX*_{2}where

*X*_{1}and*X*_{2}are realDecompose every complex matrix coefficient

*A*as*A*=*A*_{1}+*jA*_{2}where

*A*_{1}and*A*_{2}are realCarry out all complex matrix products. This yields affine expressions in

*X*_{1},*X*_{2}for the real and imaginary parts of each LMI, and an equivalent real-valued LMI is readily derived from the above observation.

For LMIs without outer factor, a streamlined version of this procedure consists of
replacing any occurrence of the matrix variable *X* =
*X*_{1} +
*jX*_{2 }by

$$\left(\begin{array}{cc}{X}_{1}& {X}_{2}\\ -{X}_{2}& {X}_{1}\end{array}\right)$$

and any fixed matrix *A = A*_{1}*
+ jA*_{2}, including real ones, by

$$\left(\begin{array}{cc}{A}_{1}& {A}_{2}\\ -{A}_{2}& {A}_{1}\end{array}\right).$$

For instance, the real counterpart of the LMI system

M < ^{H}XMX,
X = X > ^{H}I | (1) |

reads (given the decompositions
*M* = *M*_{1} +
*jM*_{2} and *X* =
*X*_{1} +
*jX*_{2} with
*M*j*, X**j*
real):

$$\begin{array}{c}{\left(\begin{array}{cc}{M}_{1}& {M}_{2}\\ -{M}_{2}& {M}_{1}\end{array}\right)}^{T}\left(\begin{array}{cc}{X}_{1}& {X}_{2}\\ -{X}_{2}& {X}_{1}\end{array}\right)\left(\begin{array}{cc}{M}_{1}& {M}_{2}\\ -{M}_{2}& {M}_{1}\end{array}\right)<\left(\begin{array}{cc}{X}_{1}& {X}_{2}\\ -{X}_{2}& {X}_{1}\end{array}\right)\\ \left(\begin{array}{cc}{X}_{1}& {X}_{2}\\ -{X}_{2}& {X}_{1}\end{array}\right)<I.\end{array}$$

Note that *X =
X** ^{H}* in turn requires
that $${X}_{1}={X}_{1}^{H}$$ and $${X}_{2}+{X}_{2}^{T}=0$$. Consequently,

*X*

_{1}and X

_{2}should be declared as symmetric and skew- symmetric matrix variables, respectively.

Assuming, for instance, that *M* ∊ **C**^{5}^{×}^{5},
the LMI system (Equation 1) would be specified as
follows:

M1=real(M), M2=imag(M) bigM=[M1 M2;-M2 M1] setlmis([]) % declare bigX=[X1 X2;-X2 X1] with X1=X1' and X2+X2'=0: [X1,n1,sX1] = lmivar(1,[5 1]) [X2,n2,sX2] = lmivar(3,skewdec(5,n1)) bigX = lmivar(3,[sX1 sX2;-sX2 sX1]) % describe the real counterpart of (1.12): lmiterm([1 1 1 0],1) lmiterm([-1 1 1 bigX],1,1) lmiterm([2 1 1 bigX],bigM',bigM) lmiterm([-2 1 1 bigX],1,1) lmis = getlmis

Note the three-step declaration of the structured matrix variable
`bigX`

,

$$bigX=\left(\begin{array}{cc}{X}_{1}& {X}_{2}\\ -{X}_{2}& {X}_{1}\end{array}\right).$$

Specify

*X*_{1}as a (real) symmetric matrix variable and save its structure description`sX1`

as well as the number`n1`

of decision variables used in*X*_{1}.Specify

*X*_{2}as a skew-symmetric matrix variable using Type 3 of`lmivar`

and the utility`skewdec`

. The command`skewdec(5,n1)`

creates a 5-by–5 skew-symmetric structure depending on the decision variables`n1 +`

1,`n1 +`

2,...Define the structure of

`bigX`

in terms of the structures`sX1`

and`sX2`

of*X*_{1}and*X*_{2}.

See Structured Matrix Variables for more details on such structure manipulations.

### Specifying c^{T}x Objectives for mincx

The LMI solver `mincx`

minimizes linear objectives
of the form *c*^{T}*x* where *x* is the
vector of decision variables. In most control problems, however, such objectives are
expressed in terms of the matrix variables rather than of *x*.
Examples include Trace(*X*) where *X* is a
symmetric matrix variable, or *u*^{T}*Xu* where *u* is a given
vector.

The function `defcx`

facilitates the derivation
of the *c* vector when the objective is an affine function of the
*matrix variables*. For the sake of illustration, consider
the linear objective

$$\text{Trace}\left(X\right)+{x}_{0}^{T}P{x}_{0}$$

where *X* and *P* are two symmetric
variables and *x*_{0} is a given vector. If
`lmsisys`

is the internal representation of the LMI system and
if *x*_{0}, *X,*
*P* have been declared by

x0 = [1;1] setlmis([]) X = lmivar(1,[3 0]) P = lmivar(1,[2 1]) : : lmisys = getlmis

the *c* vector such that $${c}^{T}x=\text{Trace}\left(X\right)+{x}_{0}^{T}P{x}_{0}$$ can be computed as follows:

n = decnbr(lmisys) c = zeros(n,1) for j=1:n, [Xj,Pj] = defcx(lmisys,j,X,P) c(j) = trace(Xj) + x0'*Pj*x0 end

The first command returns the number of decision variables in the problem and the
second command dimensions *c* accordingly. Then the
`for`

loop performs the following operations:

Evaluate the matrix variables

*X*and*P*when all entries of the decision vector*x*are set to zero except*x**j*: = 1. This operation is performed by the function`defcx`

. Apart from`lmisys`

and`j`

, the inputs of`defcx`

are the identifiers`X`

and`P`

of the variables involved in the objective, and the outputs`Xj`

and`Pj`

are the corresponding values.Evaluate the objective expression for

`X:= Xj`

and`P:= Pj`

. This yields the*j*-th entry of`c`

by definition.

In our example the result is

c = 3 1 2 1

Other objectives are handled similarly by editing the following generic skeleton:

n = decnbr(LMI system) c = zeros(n,1) for j=1:n, [matrix values] = defcx( LMI system,j,matrix identifiers) c(j) =objective(matrix values) end

### Feasibility Radius

When solving LMI problems with `feasp`

, `mincx`

, or `gevp`

, it is possible to constrain
the solution *x* to lie in the ball

*x*^{T}*x*
< *R*^{2}

where *R* > 0 is called the *feasibility
radius*. This specifies a maximum (Euclidean norm) magnitude for
*x* and avoids getting solutions of very large norm. This may
also speed up computations and improve numerical stability. Finally, the feasibility
radius bound regularizes problems with redundant variable sets. In rough terms, the
set of scalar variables is redundant when an equivalent problem could be formulated
with a smaller number of variables.

The feasibility radius *R* is set by the third entry of the
options vector of the LMI solvers. Its default value is *R* =
109. Setting *R* to a negative value means “no rigid
bound,” in which case the feasibility radius is increased during the
optimization if necessary. This “flexible bound” mode may yield
solutions of large norms.

### Well-Posedness Issues

The LMI solvers used in the LMI Lab are based on interior-point optimization
techniques. To compute feasible solutions, such techniques require that the system
of LMI constraints be strictly feasible, that is, the feasible set has a nonempty
interior. As a result, these solvers may encounter difficulty when the LMI
constraints are feasible but not *strictly feasible*, that is,
when the LMI

*L*(*x*) ≤ 0

has solutions while

*L*(*x*) < 0

has no solution.

For feasibility problems, this difficulty is automatically circumvented by
`feasp`

, which reformulates the
problem:

Find *x* such that

L(x) ≤ 0 | (2) |

as:

Minimize *t* subject to

*Lx* < *t* × *I*.

In this modified problem, the LMI constraint is always strictly feasible in
*x*, *t* and the original LMI Equation 2 is feasible if and only if the global minimum
*t*_{min} of Equation 2 satisfies

*t*_{min} ≤ 0

For feasible but not strictly feasible problems, however, the computational effort
is typically higher as `feasp`

strives to approach the
global optimum *t*_{min} = 0 to a high
accuracy.

For the LMI problems addressed by `mincx`

and `gevp`

, nonstrict feasibility
generally causes the solvers to fail and to return an “infeasibility”
diagnosis. Although there is no universal remedy for this difficulty, it is
sometimes possible to eliminate underlying algebraic constraints to obtain a
strictly feasible problem with fewer variables.

Another issue has to do with homogeneous feasibility problems such as

*A*^{T}*P + P
A* < 0, *P* > 0.

While this problem is technically well-posed, the LMI optimization is likely to
produce solutions close to zero (the trivial solution of the nonstrict problem). To
compute a nontrivial Lyapunov matrix and easily differentiate between feasibility
and infeasibility, replace the constraint *P* >
0-by-*P* > α*I* with α > 0. Note
that this does not alter the problem due to its homogeneous nature.

### Semi-Definite B(x) in gevp Problems

Consider the generalized eigenvalue minimization problem

Minimize *λ* subject to

A(x) <
λB(x),
B(x) > 0,
C(x) <0. | (3) |

Technically, the positivity of *B*(*x*) for
some *x* ∊
*R** ^{n}* is
required for the well-posedness of the problem and the applicability of
polynomial-time interior-point methods. Hence problems where

$$B\left(x\right)=\left(\begin{array}{cc}{B}_{1}\left(x\right)& 0\\ 0& 0\end{array}\right)$$

with *B*_{1}(*x*) > 0 strictly feasible, cannot be directly solved with `gevp`

. A simple remedy consists of
replacing the constraints

*A*(*x*) < *B(x*),
*B*(*x*) > 0

by

$$A\left(x\right)<\left(\begin{array}{cc}Y& 0\\ 0& 0\end{array}\right),\text{\hspace{1em}}Y<\lambda {B}_{1}\left(x\right),\text{\hspace{1em}}{B}_{1}\left(x\right)>0$$

where *Y* is an additional symmetric variable of proper
dimensions. The resulting problem is equivalent to Equation 3 and can be solved directly with `gevp`

.

### Efficiency and Complexity Issues

As explained in Tools for Specifying and Solving LMIs, the term-oriented description of LMIs used in the LMI Lab typically leads to higher efficiency than the canonical representation

A_{0} + x_{1}A_{1} + ... + x < 0._{N}A_{N} | (4) |

This is no longer true, however, when the number of variable terms is nearly equal
to or greater than the number *N* of decision variables in the
problem. If your LMI problem has few free scalar variables but many terms in each
LMI, it is therefore preferable to rewrite it as Equation 4 and to specify it in this form. Each scalar
variable *x** _{j}* is then
declared independently and the LMI terms are of the form

*x*

_{j}*A*

*.*

_{j}If *M* denotes the total row size of the LMI system and
*N* the total number of scalar decision variables, the flop
count per iteration for the `feasp`

and `mincx`

solvers is proportional
to

*N*^{3}when the least-squares problem is solved via Cholesky factorization of the Hessian matrix (default) [2].*M*-by-*N*^{2}when numerical instabilities warrant the use of QR factorization instead.

While the theory guarantees a worst-case iteration count proportional to
*M*, the number of iterations actually performed grows slowly
with *M* in most problems. Finally, while `feasp`

and `mincx`

are comparable in
complexity, `gevp`

typically demands more
computational effort. Make sure that your LMI problem cannot be solved with
`mincx`

before using `gevp`

.

### Solving M + P^{T}XQ + Q^{T}X^{T}P < 0

In many output-feedback synthesis problems, the design can be performed in two steps:

Compute a closed-loop Lyapunov function via LMI optimization.

Given this Lyapunov function, derive the controller state-space matrices by solving an LMI of the form

*M*+*P*+^{T}XQ*Q*< 0^{T}X^{T}P**(5)**where

*M*,*P*,*Q*are given matrices and*X*is an unstructured*m*-by-*n*matrix variable.

It turns out that a particular solution *X*c of Equation 5 can be computed via simple linear algebra
manipulations [1]. Typically,
*X*c corresponds to the center of the ellipsoid of matrices
defined by Equation 5.

The function` basiclmi`

returns the “explicit”
solution *X*c:

Xc = basiclmi(M,P,Q)

Since this central solution sometimes has large norm, `basiclmi`

also offers the option of computing an approximate least-norm solution of Equation 5. This is done by

X = basiclmi(M,P,Q,'Xmin')

and involves LMI optimization to minimize ∥*X* ∥.