%  Example 5.14.  Robust stability, nominal performance, and robust 
%    performance analysis using the structured singular value.
clear
%  Define the plant.
A=[0  1   0   0
   0  0   0 -1.6
   0  0   0   1
   0 1.8  0   0];
Bu=[ 0   0
    11   0
     0   0
     0  12];
Cm=eye(4);
Cy=[1 0 0 0
    0 0 1 0];
Kr=Cy';
%  Generate a state feedback controller for the plant.
p=[-5+10*j -5-10*j -6+12*j -6-12*j];
K=place(A,Bu,p)
%  Generate the nominal closed loop system Mydwd.
Mydwd=pck(A-Bu*K*Cm,Bu,-0.5*K*Cm,zeros(2));
%  Evaluate robust stability by computing the maximum singular value
%  of the transfer function Mydwd.  Note:  There is only one uncertainty block
%  so the singular value is used for robust stability analysis.
normmMydwd=hinfnorm(Mydwd)
%  Compute the principle gains for plotting.
w=logspace(-1,2);
frydwd = frsp(Mydwd,w);
pgydwd=vsvd(frydwd);
%  Evaluate nominal performance using the singular values of Myw.
%  Generate the nominal closed loop system without weighting.
Arod=pck(A-Bu*K*Cm,Bu*K*Kr,-Cy,eye(2));
%  Define the weighting.
Aw=[0 1;-1 -2]; Aw=[Aw zeros(size(Aw));zeros(size(Aw)) Aw];
Bw=[0;20];
Bw=[Bw zeros(size(Bw));zeros(size(Bw)) Bw];
Cw=[1 0];
Cw=[Cw zeros(size(Cw));zeros(size(Cw)) Cw];
Wi=pck(Aw,Bw,Cw,zeros(2));
%  Use sysic to combine the nominal plant with the weighting functions.
systemnames='Arod Wi';
inputvar='inputt(2)';
outputvar='[Wi]';
input_to_Wi='[Arod]';
input_to_Arod='[inputt]';
sysoutname='Myw';
sysic
%  Compute the maximum singular value of Myd.
normmyw=hinfnorm(Myw)
%  Compute the principle gains for plotting.
fryw = frsp(Myw,w);
pgyw=vsvd(fryw);
%  Evaluate robust performance using the structured singular values of M.
%  Generate the nominal closed loop system without weighting.
D=[zeros(2) 0.5*K*Kr;zeros(2) eye(2)];
MArod=pck(A-Bu*K*Cm,[Bu Bu*K*Kr],[-0.5*K*Cm;-Cy],D);
%  Use "sysic" to generate the nominal closed loop system M.
systemnames='MArod Wi';
inputvar='inputt(4)';
outputvar='[MArod(1:2); Wi]';
input_to_Wi='[MArod(3:4)]';
input_to_MArod='[inputt]';
sysoutname='M';
sysic
%  Compute the structured singular value.
fr = frsp(M,w);
blk=[2 2;2 2];
muM=mu(fr,blk);
pgydwdmax=veval('max',pgydwd);
pgywmax=veval('max',pgyw);
muMu=veval('max',muM);
figure(1)
set(0,'DefaultAxesFontName','times')
set(0,'DefaultAxesFontSize',16)
set(0,'DefaultTextFontName','times')
clf
vplot('liv,lm',sel(muM,1,1),'r-',pgydwdmax,'r-.',pgywmax,'r--')
xlabel('Frequency (rad/sec)')
ylabel('Magnitude')
grid
legend('Robust performance', 'Robust stability', 'Nominal performance')
%  Compute the maximum of the upper bound on the SSV.
maxmu=max(muM(:,1))
