This is what I covered in week 3 (this is almost an introductory
course in MIMO control).
Next week (today) I will talk about chapter 4 (system theory).
-Sigurd Skogestad
------------------------------------------------------------------
Chapter 3, continue
------------------------------------------------------------------
3.4 Control of MIMO plants
---------------------------
A few "simple-minded" approaches to MIMO control:
1. Decentralized control; K is diagonal
o Does nothing about interaction - poor performance (even nominally)
o But normally robust
o May often tune on-line (little modelling effort)
o Important issue: Select pairing (to achieve G diagonally dominant)
2. Compensator approaches and decoupling
o First make something simple which takes care of interactions
o Then design diagonal controller for performance etc.
o Decoupler: E,g. inverse -based controller
K = l G^{-1}
Often sensitive to uncertainty, Also: Decoupling
may not be required, e.g. for disturbance rejection.
o SVD-controller: More flexible. Usually based on SVD of G
around crossover frequency.
3.5 Introduction to MIMO RHP-zeros
----------------------------------
Square MIMO plant: Poles and zeros are in most cases poles and zeros of
det G(s).
More generally: Zeris are where G(s) looses rank.
NOTE: Poles are essentially poles of elements of G, but zeros are NOT
usually zeros of elements of G !!
Right half plane (RHP) zero: Fundamental limitation on achievable control
(as for SISO). Why? G(s)^-1 is unstable.
MIMO pole and zeros have directions.
Poles: Look at SVD of G(p): Pole direction is direction where G(p) is
infinite.
Zeros: Look at SVD of G(z). Zero direction is where gain of G(z) is
zero (i.e. it goes to zero for s=z).
Numerical: Use state-space form, see ch. 4
See 2x2 example with H-infinity design: Note that we can MOVE effect
of RHP-zero to particular outputs (unless zero is pinned).
Matlab-file is given below.
3.6 Introduction to robustness for MIMO plants
-----------------------------------------------
Example 1: Spinning satellite
*****************************
Usual stability margins (GM, PM) one loop at a time: Looks OK
BUT singular value of S and T has peak around 10.
This ALWAYS signals a robustness problem (and performance is alo
poor).
Matlab file given below.
Example 2: Distillation process
********************************
Here usual margins are OK + there is NO large peak on S and T.
Still poor robustness with input uncertainty, see simulation.
Whats wrong? The input uncertainty moves large input signals over
to direction where plant has large gain.
NOTE: We use an inverse-based controller. This problem CANNOT
happen with a diagonal controller with no peak on S or T.
Matlab file given below.
3.7 General control problem formulation
----------------------------------------
Important section!
Note that we can handle non-square plants, we need not measure all
the outputs, etc. etc.
Note that the general formulation (P) has already been used in the files for
H-infinty design (S/KS).
May use sysic in the MATLAB mu-toolbox to generate P.
Exercise: Find P for two-degrees of freedom controller with noise and
disturbance. z=[y-r u]^T. Also find N from w to z with two methods.
*** End lecture week 3 *********
________________________________________________________________________
THREE MATLAB FILES ARE GIVEN BELOW (RHP-zero, Satellite, distillation)
-----------------------------------------------------------------------
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Simple 2x2 plant with RHP zero.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
g11 = nd2sys(1, [0.2 1.2 1]); g12=g11; g22=mmult(g11,2);
g21 = nd2sys([2 1],[0.2 1.2 1]);
G = abv(sbs(g11,g12),sbs(g21,g22));
spoles(G) % 8 poles
G = sysbal(G); spoles(G) % minimal realization only 4 poles
format short e; szeros(G) % Mu-tools command yields extra zeros at infinity
[A,B,C,D]=unpck(G);
tzero(A,B,C,D) % Try Robust control toolbox instead, zero at z=0.5
% Analyze zero direction
G05 = frsp(G,-0.5i)
[u,s,v]=vsvd(G05) % look at last column in u and v
%Open-loop simulation
y1=trsp(G,[1;0],10,0.1);
y2=trsp(G,[0;1],10,0.1);
y3=trsp(G,[1;-1],10,0.1);
subplot(231);vplot('iv,d',y1)
subplot(232);vplot('iv,d',y2)
subplot(233);vplot('iv,d',y3)
title('OPEN-LOOP SIMULATION')
% H-INFINITY S-KS design
% Require wb, slope 1 and peak less than M
M=1.5; wb=0.5*z; wp1 = nd2sys([1/M wb], [1 1.e-6]); wp2=wp1; %weight1
Wp=daug(wp1,wp2); Wu=1.0*eye(2);
systemnames = 'G Wp Wu'
inputvar = '[ r(2); u(2)]';
outputvar = '[Wp; Wu; r-G]';
input_to_G = '[u]';
input_to_Wp = '[r-G]';
input_to_Wu = '[u]';
sysoutname = 'P';
cleanupsysic = 'yes';
sysic;
nmeas=2; nu=2; gmn=0.667; gmx=20; tol=0.001;
[khinf,ghinf,gopt] = hinfsyn(P,nmeas,nu,gmn,gmx,tol); K=khinf;
L = mmult(G,K); i = eye(2); S = minv(madd(i,L)); T=msub(i,S);
omega=logspace(-3,2,121);
Sf=frsp(S,omega); [u,Sfssav,v]=vsvd(Sf);
% Simulation of bad (y1) and good (y2) direction
y1sav=trsp(T,[1;-1],10,0.1);
% new performance weight
M=1.5; wb=50*z; wp2 = nd2sys([1/M wb], [1 1.e-6]); %weight2
% when wb2=10 we require much better performance for output 2
Wp=daug(wp1,wp2);
systemnames = 'G Wp Wu'
inputvar = '[ r(2); u(2)]';
outputvar = '[Wp; Wu; r-G]';
input_to_G = '[u]';
input_to_Wp = '[r-G]';
input_to_Wu = '[u]';
sysoutname = 'P';
cleanupsysic = 'yes';
sysic;
nmeas=2; nu=2; gmn=0.667; gmx=20; tol=0.001;
[khinf,ghinf,gopt] = hinfsyn(P,nmeas,nu,gmn,gmx,tol); K=khinf;
L = mmult(G,K); S = minv(madd(i,L)); T=msub(i,S);
Sf=frsp(S,omega); [u,Sfs,v]=vsvd(Sf);
y1=trsp(T,[1;-1],10,0.1);
% Plot
figure(1); subplot(211);vplot('liv,lm',Sfssav,Sfs,'--',1,':')
title('SINGULAR VALUES OF S:');
text(4,0.02,'DESIGN 1: SOLID LINE');text(4,0.005,'DESIGN 2: DASHED LINE');
axis([0.01,100,0.001,10]);
subplot(212);vplot('iv,d',y1sav,y1,'--',1,':',-1,':')
title('RESPONSE TO CHANGE IN REFERENCE, r = [1 -1]:');xlabel('Time');
text(4,-1.4,'OUTPUT 2');text(4,1.4,'OUTPUT 1');
axis([0,5,-2,2]);
drawnow;
------------------------------------------------------------------------
%%%%%%%%%%%%%%%%%%%%%
% Spinning satellite
%%%%%%%%%%%%%%%%%%%%%
alpha = 10;
a = [0 alpha ; -alpha 0];
b = [1 0; 0 1];
c = [1 alpha; -alpha 1];
d = [ 0 0; 0 0];
G = pck(a,b,c,d);
I2 = eye(2);
K=I2;
% KG and GK are both equal to G
L = mmult(G,K);
S = minv(madd(I2,L));
T = msub(I2,S);
hinfnorm(S) % H-infinity norm of S is 10.05 (peak of singular value of S).
hinfnorm(T) % H-infinity norm of T is 10.05 (peak of singular value of T).
--------------------------------------------------------------------------
%%%%%%%%%%%%%%%%%%%%%%%%
% DISTILLATION PROCESS
%%%%%%%%%%%%%%%%%%%%%%%%
%
% The example process is from:
% Skogestad, Morari and Doyle,``Robust Control of Ill-Conditioned
% Plants: High-Purity Distillation'', IEEE Atomat. Control 33,
% 1092-1105 (1988).
%
G0 = [87.8 -86.4; 108.2 -109.6];
[U,sval,V]=svd(G0);
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));
I2=eye(2); K = Kinv;
% Nominal
GK = mmult(G,K);
%% May perform minimal realization to get rid of the extra states
GK = sysbal(GK);
S = minv(madd(I2,GK)); T = msub(I2,S);
% With 20% INPUT uncertainty
Unc = [1.2 0; 0 0.8];
GKu = mmult(G,Unc,K);
Su = minv(madd(I2,GKu)); Tu = msub(I2,Su);
hinfnorm(S), hinfnorm(Su) % Peak of S; 1 nominally and 14.2 with uncertainty
hinfnorm(T), hinfnorm(Tu) % Peak of T; 1 nominally and 14.2 with uncertainty
% TIME simulation
kr=nd2sys(1,[5 1]); % 5 min filter on reference change
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);
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,'--');
axis([0 100 -0.2 2.7]);text(77,1.2,'y1'),text(77,0.2,'y2');title('OUTPUTS');
subplot(212);vplot(u,uu,'--'),title('INPUTS');figure(1);xlabel('TIME (min)');
% For EXERCISE.
% Now with relative uncertainty del BOTH at the input and output
del=0.11; % 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);
yu = trsp(Tru,[1;0],100,.1);
subplot(111);vplot(y,yu,'--');
title('11 percent UNC. ON INPUTS AND OUTPUTS:');
axis([0 100 -0.2 2.7]); text(77,1.2,'y1'), text(77,0.2,'y2');
xlabel('TIME (min)');
%----------------------------------------------------------------------