Accelerating the pace of engineering and science

# mussvextract

Extract muinfo structure returned by mussv

## Synopsis

[VDelta,VSigma,VLmi] = mussvextract(muinfo)

## Description

A structured singular value computation of the form

```[bounds,muinfo] = mussv(M,BlockStructure)
```

returns detailed information in the structure muinfo. mussvextract is used to extract the compressed information within muinfo into a readable form.

The most general call to mussvextract extracts three usable quantities: VDelta, VSigma, and VLmi. VDelta is used to verify the lower bound. VSigma is used to verify the Newlin/Young upper bound and has fields DLeft, DRight, GLeft, GMiddle, and GRight. VLmi is used to verify the LMI upper bound and has fields Dr, Dc, Grc, and Gcr. The relation/interpretation of these quantities with the numerical results in bounds is described below.

### Upper Bound Information

The upper bound is based on a proof that det(I - M*Delta) is nonzero for all block-structured matrices Delta with norm smaller than 1/bounds(1). The Newlin/Young method consists of finding a scalar β and matrices D and G, consistent with BlockStructure, such that

Here DL, DR, GL, GM and GR correspond to the DLeft, DRight, GLeft, GMiddle and GRight fields respectively.

Because some uncertainty blocks and M need not be square, the matrices D and G have a few different manifestations. In fact, in the formula above, there are a left and right D and G, as well as a middle G. Any such β is an upper bound of mussv(M,BlockStructure).

It is true that if BlockStructure consists only of complex blocks, then all G matrices will be zero, and the expression above simplifies to

$\overline{\sigma }\left({D}_{L}M{D}_{R}^{-1}\right)\le \beta .$

The LMI method consists of finding a scalar β and matrices D and G, consistent with BlockStructure, such that

${M}^{\prime }{D}_{r}M-{\beta }^{2}{D}_{c}+j\left({G}_{cr}M-{M}^{\prime }{G}_{rc}\right)\le 0$

is negative semidefinite. Again, D and G have a few different manifestations to match the row and column dimensions of M. Any such β is an upper bound of mussv(M,BlockStructure). If BlockStructure consists only of complex blocks, then all G matrices will be zero, and negative semidefiniteness of M´Dr M2Dc is sufficient to derive an upper bound.

### Lower Bound Information

The lower bound of mussv(M,BlockStructure) is based on finding a "small" (hopefully the smallest) block-structured matrix VDelta that causes det(I - M*VDelta) to equal 0. Equivalently, the matrix M*VDelta has an eigenvalue equal to 1. It will always be true that the lower bound (bounds(2)) will be the reciprocal of norm(VDelta).

## Examples

Suppose M is a 4-by-4 complex matrix. Take the block structure to be two 1-by-1 complex blocks and one 2-by-2 complex block.

```rng(0,'twister')
M = randn(4,4) + sqrt(-1)*randn(4,4);
BlockStructure = [1 1;1 1;2 2];
```

You can calculate bounds on the structured singular value using the mussv command and extract the scaling matrices using mussvextract.

```[bounds,muinfo] = mussv(M,BlockStructure);
[VDelta,VSigma,VLmi] = mussvextract(muinfo);
```

You can first verify the Newlin/Young upper bound with the information extracted from muinfo. The corresponding scalings are Dl and Dr.

```Dl = VSigma.DLeft
```
```Dl =

1.0000         0         0         0
0    0.7437         0         0
0         0    1.0393         0
0         0         0    1.0393

```
`Dr = VSigma.DRight `
```Dr =

1.0000         0         0         0
0    0.7437         0         0
0         0    1.0393         0
0         0         0    1.0393

```
`[norm(Dl*M/Dr) bounds(1)] `
```ans =

6.2950    6.2950
```

You can first verify the LMI upper bound with the information extracted from muinfo. The corresponding scalings are Dr and Dc.

```Dr = VLmi.Dr;
Dc = VLmi.Dc;
eig(M'*Dr*M - bounds(1)^2*Dc)
```
```ans =

-0.0000 - 0.0000i
-17.7242 - 0.0000i
-33.8550 + 0.0000i
-41.2013 - 0.0000i
```

Note that VDelta matches the structure defined by BlockStructure, and the norm of VDelta agrees with the lower bound,

```VDelta
```
```VDelta =

0.1301 - 0.0922i        0                  0                  0
0            -0.0121 - 0.1590i        0                  0
0                  0            -0.0496 - 0.0708i   0.1272 - 0.0075i
0                  0             0.0166 - 0.0163i   0.0076 + 0.0334i
```
`[norm(VDelta) 1/bounds(2)] `
```ans =

0.1595    0.1595
```

and that M*VDelta has an eigenvalue exactly at 1.

```eig(M*VDelta)
```
```ans =

1.0000 - 0.0000i
-0.2501 - 0.1109i
0.0000 + 0.0000i
-0.3022 + 0.2535i
```

Keep the matrix the same, but change BlockStructure to be a 2-by-2 repeated, real scalar block and two complex 1-by-1 blocks. Run mussv with the 'C' option to tighten the upper bound.

```BlockStructure2 = [-2 0; 1 0; 1 0];
[bounds2,muinfo2] = mussv(M,BlockStructure2,'C');
```

You can compare the computed bounds. Note that bounds2 should be smaller than bounds, because the uncertainty set defined by BlockStructure2 is a proper subset of that defined by BlockStructure.

```[bounds; bounds2]
```
```ans =

6.2950    6.2704
5.1840    5.1750
```

You can extract the D, G and Delta from muinfo2 using mussvextract.

```[VDelta2,VSigma2,VLmi2] = mussvextract(muinfo2);
```

As before, you can first verify the Newlin/Young upper bound with the information extracted from muinfo. The corresponding scalings are Dl, Dr, Gl, Gm and Gr.

```Dl = VSigma2.DLeft;
Dr = VSigma2.DRight;
Gl = VSigma2.GLeft;
Gm = VSigma2.GMiddle;
Gr = VSigma2.GRight;
dmd = Dl*M/Dr/bounds2(1) - sqrt(-1)*Gm;
SL = (eye(4)+Gl*Gl)^-0.25;
SR = (eye(4)+Gr*Gr)^-0.25;
norm(SL*dmd*SR)
```
```ans =

1.0000
```

You can first verify the LMI upper bound with the information extracted from muinfo. The corresponding scalings are Dr, Dc, Grc and Gcr.

```Dr = VLmi2.Dr;
Dc = VLmi2.Dc;
Grc = VLmi2.Grc;
Gcr = VLmi2.Gcr;
eig(M'*Dr*M - bounds(1)^2 *Dc + j*(Gcr*M-M'*Grc))
```
```ans =

-69.9757 + 0.0000i
-11.2139 - 0.0000i
-19.2766 - 0.0000i
-40.2869 - 0.0000i
```

VDelta2 matches the structure defined by BlockStructure, and the norm of VDelta2 agrees with the lower bound,

```VDelta2
```
```VDelta2 =

0.1932                  0                  0                  0
0             0.1932                  0                  0
0                  0            -0.1781 - 0.0750i        0
0                  0                  0             0.0941 + 0.1688i
```
`[norm(VDelta2) 1/bounds2(2)] `
```ans =

0.1932    0.1932
```

and that M*VDelta2 has an eigenvalue exactly at 1.

```eig(M*VDelta2)
```
```  ans =

1.0000 + 0.0000i
-0.4328 + 0.1586i
0.1220 - 0.2648i
-0.3688 - 0.3219i

```