-Sigurd
%
% Chapter 3.
% This m-file contains simulations using various controllers
% for the distillation process.
% File name: ieee_more.m
%
% Includes:
% 1. Remark 3: S/KS-design
% 2. Ex. 3.6: Attempt to robustify inverse-based controller
% 3. Ex. 3.7: SVD-controllers
% 4. Ex. 3.8: Combined input and output uncertainty
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% S-KS DESIGN (see Remark 3)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
G0 = [87.8 -86.4; 108.2 -109.6];
dyn = nd2sys(1,[75 1]); Dyn = daug(dyn,dyn); G = mmult(Dyn,G0);
A= 1.e-4; M=2; wb=0.05;
wp = nd2sys([1/M wb], [1 wb*A]); Wp=daug(wp,wp);
wu = 1; Wu = daug(wu,wu);
systemnames = 'G Wp Wu'
inputvar = '[ r(2); u(2)]';
outputvar = '[Wp; Wu; r-G]'; % NEGATIVE feedback
input_to_G = '[u]';
input_to_Wp = '[r-G]';
input_to_Wu = '[u]';
sysoutname = 'P';
cleanupsysic = 'yes'; sysic;
nmeas=2; nu=2; gmn=1; gmx=100; tol=0.01;
[khinf,ghinf,gopt] = hinfsyn(P,nmeas,nu,gmn,gmx,tol);
K=khinf;
% TIME simulation
% Nominal
GK = mmult(G,K); I2=eye(2); S = minv(madd(I2,GK)); T = msub(I2,S);
kr=nd2sys(1,[5 1]); Kr=daug(kr,kr); Tr = mmult(T,Kr);
y = trsp(Tr,[1;0],100,.1);
u = trsp(mmult(K,S,Kr),[1;0],100,.1);
% With 20% uncertainty
Unc = [1.2 0; 0 0.8]; GKu = mmult(G,Unc,K);
Su = minv(madd(I2,GKu)); Tu = msub(I2,Su); Tru = mmult(Tu,Kr);
yu = trsp(Tru,[1;0],100,.1);
uu = trsp(mmult(K,Su,Kr),[1;0],100,.1);
subplot(211);vplot(y,yu,'--');title('OUTPUTS')
subplot(212);vplot(u,uu,'--');title('INPUTS');
xlabel('TIME (min)');
% Simulation shows that the controller is very sensitive
% to uncertainty. Large and Long peak in y1 and y2.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Ex. 3.6: Try to robustify inverse-based
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
G0 = [87.8 -86.4; 108.2 -109.6];
dyn = nd2sys(1,[75 1]); Dyn = daug(dyn,dyn); G = mmult(Dyn,G0);
% Inverse-based controller
dynk = nd2sys([75 1],[1 1.e-6],0.7); Dynk = daug(dynk,dynk);
Kinv = mmult(Dynk,minv(G0));
% Try to robustify with respect to coprime uncerainty
Gs = mmult(G,Kinv);
[a,b,c,d]=unpck(Gs);
gamrel=1.1;
[Ac,Bc,Cc,Dc,gammin]=coprimeunc(a,b,c,d,gamrel); % gammin=1.4142
Ks=pck(Ac,-Bc,Cc,-Dc); % Change from positive to negative feedback
K = mmult(Kinv,Ks);
% TIME simulation
% Nominal
GK = mmult(G,K); I2=eye(2); S = minv(madd(I2,GK)); T = msub(I2,S);
kr=nd2sys(1,[5 1]); Kr=daug(kr,kr); Tr = mmult(T,Kr);
y = trsp(Tr,[1;0],100,.1);
u = trsp(mmult(K,S,Kr),[1;0],100,.1);
% With 20% uncertainty
Unc = [1.2 0; 0 0.8]; GKu = mmult(G,Unc,K);
Su = minv(madd(I2,GKu)); Tu = msub(I2,Su); Tru = mmult(Tu,Kr);
yu = trsp(Tru,[1;0],100,.1);
uu = trsp(mmult(K,Su,Kr),[1;0],100,.1);
subplot(211);vplot(y,yu,'--');title('OUTPUTS')
subplot(212);vplot(u,uu,'--');title('INPUTS');
xlabel('TIME (min)');
% Not very successful - almost identical with the original Kinv
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Ex. 3.7 SVD designs
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
G0 = [87.8 -86.4; 108.2 -109.6];
dyn = nd2sys(1,[75 1]); Dyn = daug(dyn,dyn); G = mmult(Dyn,G0);
[U,sval,V]=svd(G0);
% You must try these one at a time:
c1=0.005; c2=0.005; % Decentralized: Robust, but poor response
%c1 = 0.005; c2=0.05; % Good response and also robust
%c1=0.7/197; c2=0.7/1.39;% Inverse-based: Excellent nominal, but not robust
k1 = nd2sys([75 1],[1 1.e-6],c1); k2 = nd2sys([75 1],[1 1.e-6],c2);
Kdiag = daug(k1,k2);
K = mmult(V, Kdiag,U');
% TIME simulation
% Nominal
GK = mmult(G,K); I2=eye(2); S = minv(madd(I2,GK)); T = msub(I2,S);
kr=nd2sys(1,[5 1]); Kr=daug(kr,kr); Tr = mmult(T,Kr);
y = trsp(Tr,[1;0],100,.1);
u = trsp(mmult(K,S,Kr),[1;0],100,.1);
% With 20% uncertainty
Unc = [1.2 0; 0 0.8]; GKu = mmult(G,Unc,K);
Su = minv(madd(I2,GKu)); Tu = msub(I2,Su); Tru = mmult(Tu,Kr);
yu = trsp(Tru,[1;0],100,.1);
uu = trsp(mmult(K,Su,Kr),[1;0],100,.1);
subplot(211);vplot(y,yu,'--');title('OUTPUTS')
subplot(212);vplot(u,uu,'--');title('INPUTS');
xlabel('TIME (min)');
% Time simulation with higher-order dynamics at the input
dyni1 = nd2sys(1,[0.2 1]); dyni = mmult(dyni1,dyni1,dyni1,dyni1,dyni1);
dyni2 = daug(dyni,dyni); GKu2 = mmult(G,dyni2,Unc,K);
Su2 = minv(madd(I2,GKu2)); Tu2 = msub(I2,Su2);
hinfnorm(Su2), hinfnorm(Tu2)
Tru2 = mmult(Tu2,Kr);
yu2 = trsp(Tru2,[1;0],100,.1);
uu2 = trsp(mmult(K,Su2,Kr),[1;0],100,.1);
subplot(211);vplot(y,yu,'--',yu2,'-.');title('OUTPUTS')
subplot(212);vplot(u,uu,'--',uu2,'-.');title('INPUTS');
% The inverse-based is the only one affected by this higher-order dynamics
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Ex. 3.8 Combined input and output unc. for inverse-based control
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
G0 = [87.8 -86.4; 108.2 -109.6];
dyn = nd2sys(1,[75 1]); Dyn = daug(dyn,dyn); G = mmult(Dyn,G0);
% Inverse-based controller
dynk = nd2sys([75 1],[1 1.e-6],0.7); Dynk = daug(dynk,dynk);
Kinv = mmult(Dynk,minv(G0)); K=Kinv;
% TIME simulation
% Nominal
GK = mmult(G,K); I2=eye(2); S = minv(madd(I2,GK)); T = msub(I2,S);
kr=nd2sys(1,[5 1]); Kr=daug(kr,kr); Tr = mmult(T,Kr);
y = trsp(Tr,[1;0],100,.1);
u = trsp(mmult(K,S,Kr),[1;0],100,.1);
% With relative uncertainty del both at the input and output
del=0.12; % Inverse-controller: Unstable for del larger than 0.12
Unci = [1+del 0; 0 1-del]; Unco = [1-del 0; 0 1+del];
GKu = mmult(Unco,G,Unci,K);
Su = minv(madd(I2,GKu)); Tu = msub(I2,Su); Tru = mmult(Tu,Kr);
Su=sysbal(Su);spoles(Su)
yu = trsp(Tru,[1;0],100,.1);
uu = trsp(mmult(K,Su,Kr),[1;0],100,.1);
subplot(211);vplot(y,yu,'--');title('OUTPUTS')
subplot(212);vplot(u,uu,'--');title('INPUTS');
xlabel('TIME (min)');