%
% C. Tomlin Example 2 from Lecture Note 6 (originally from M. Johannson, 
% "Piecewise linear control systems", PhD thesis, Lund University, 1998).
%
% EECE 571M/491M
% Spring 2008
% M. Oishi
%


>> A1 = [-5 -4; -1 -2];
>> A2 =[-2 -4; 20 -2];
>> setlmis([]);
>> P1=lmivar(1,[2 1]);
>> P2=lmivar(1,[2 1]);
>> lmiterm([-1 1 1 P1],1,1);                       % LMI #1: P1
>> lmiterm([-2 1 1 P2],1,1);                       % LMI #2: P2
>> lmiterm([-3 1 1 P1],A1',1,'s');                 % LMI #3: A1'*P1+P1*A1
>> lmiterm([-3 1 1 P2],A2',1,'s');                 % LMI #3: A2'*P2+P2*A2
>> tmpsys=getlmis;
>> [tmin,pfeas]=feasp(tmpsys);

 Solver for LMI feasibility problems L(x) < R(x)
    This solver minimizes  t  subject to  L(x) < R(x) + t*I
    The best value of t should be negative for feasibility

 Iteration   :    Best value of t so far 
 
     1                        0.122000
     2                       -0.022603

 Result:  best value of t:    -0.022603
          f-radius saturation:  0.000% of R =  1.00e+09
 
>> %%%%% This is a good result since best value of t < 0
>> Pfinal1 = dec2mat(tmpsys,pfeas,P1)

Pfinal1 =

    0.9332   -0.7543
   -0.7543    0.7265

>> eig(Pfinal1)

ans =

    0.0685
    1.5912

>> Pfinal2 = dec2mat(tmpsys,pfeas,P2);

Pfinal2 =

    0.6477    0.2742
    0.2742    0.1442

>> eig(Pfinal2)

ans =

    0.0237
    0.7681

>> R = A1'*Pfinal1+Pfinal1*A1+A2'*Pfinal2+Pfinal2*A2

R =

    0.5526    0.0170
    0.0170    0.3586

>> eig(R)

ans =

    0.3571
    0.5541