checkpp.m

线性时变系统控制器设计的工具包

function vout = checkpp(vlpv,ny,nu,gam,xmat,ymat,fparm,gparm,gest)

% function vout = checkpp(vlpv,ny,nu,gam,xmat,ymat,fparm,gparm,gest)
% or
% function vout = checkpp(vlpv,ny,nu,gam,xmat,ymat)
%
% Form the matrix PHI for which the (frozen-point) closed-loop poles have 
% real part > -r (i.e. slow dynamics) if we modify the coupling LMI to be
% [X I;I Y] + 1/2r * (PHI + PHI') > 0
% 
% IMPORTANT: D11 is assumed to be zero. 
%
% Inputs: vlpv		varying matrix of the lpv plant at parm values
%	  ny,nu		number of measurements & controls, respectively
%	  gam		optimal gamma from vlpvsol
%	  xmat,ymat 	LMI solutions from vlpvsol in varying matrices
%	  fparm,gparm 	basis-function values at parm values
%			If X or Y is constant, enter [] for fparm or gparm.
%			If both are constant, exclude fparm and gparm.
%	  gest		estimate for gam
%
% Output: vout is PHI

if nargin == 0
 disp('usage: [vctrl,vclp] = checkpp(vlpv,ny,nu,gam,xmat,ymat,fparm,gparm)');
 disp('       or');
 disp('       [vctrl,vclp] = checkpp(vlpv,ny,nu,gam,xmat,ymat)'); 
 return;
end

if ~exist('gest')
 gest = gam;
end

% If vlpv, xmat, or ymat aren't varying matrices, vpck them
[typ,dum2,dum3,npts1] = minfo(vlpv);
if typ == 'syst'
  vlpv = vpck(vlpv,1);
  [typ,dum2,dum3,npts1] = minfo(vlpv);
end
viv = getiv(vlpv);
[typ,xrow,xcol,numx] = minfo(xmat);
if typ == 'cons'
  xmat = vpck(xmat,1);
  [typ,xrow,xcol,numx] = minfo(xmat);
end
[typ,yrow,ycol,numy] = minfo(ymat);
if typ == 'cons'
  ymat = vpck(ymat,1);
  [typ,yrow,ycol,numy] = minfo(ymat);
end
if (xrow ~= xcol) | (yrow ~= ycol)
  error('X and Y must be square');
end

% If constant X and/or Y
if ~exist('fparm') 
 fparm = vpck(ones(npts1,1),1:npts1);
end
if ~exist('gparm')
 gparm = vpck(ones(npts1,1),1:npts1);
end
if isempty(fparm)
 fparm = vpck(ones(npts1,1),1:npts1);
end
if isempty(gparm)
 gparm = vpck(ones(npts1,1),1:npts1);
end

% Check that the number of grid points and parameters is consistent
 [dum1,nbasisx,ckvx,npts2] = minfo(fparm);
 [dum1,nbasisy,ckvy,npts3] = minfo(gparm);
 if (ckvx ~= 1) | (ckvy ~= 1)
   error('Error: Basis data should be column vectors.');
 end
 if (nbasisx ~= numx) | (nbasisy ~= numy)
   error('Error: Basis dimension not consistent');
 end
 if (npts1 ~= npts2) | (npts2 ~= npts3) | (npts3 ~= npts1) 
   disp(['There are ' int2str(npts1) ' points in vlpv']);
   disp(['There are ' int2str(npts2) ' points in fparm']);
   disp(['There are ' int2str(npts3) ' points in gparm']);
   error('Error: number of grid points inconsistent');
 end  
grid_pts = npts1;

sys = xtracti(vlpv,1,1);
[systype,no,ni,nx] = minfo(sys);
ne = no-ny;
ne1 = ne-nu;
nd = ni-nu;
nd1 = nd-ny;
if (nx ~= xrow) | (nx ~= yrow)
 error('State dimension inconsistent with Lyapunov matrices');
end

if min(numx,numy) > 1
 disp('Warning: Both X and Y are varying. Output applies to LTI case only.')
end

vout = zeros(2*grid_pts*nx,2*nx);

for i=1:grid_pts

 % Form the lyapunov matrices for parm_i
 X = zeros(nx);
 Y = zeros(nx);
 fvec = xtracti(fparm,i,1);
 gvec = xtracti(gparm,i,1);
 for k=1:nbasisx
  X = X + fvec(k) * xtracti(xmat,k,1);
 end
 for k=1:nbasisy
  Y = Y + gvec(k) * xtracti(ymat,k,1);
 end
 if numy == 1
  M = I;
  N = eye(nx) - X*Y;
 elseif numx == 1
  M = eye(nx) - Y*X;
  N = eye(nx);
 else
  M = eye(nx) - Y*X;
  N = eye(nx);
 end
  
% Get state-space data for parm_i.
 sys = xtracti(vlpv,i,1);
 [a,b1,b2,c1,c2,d11,d12,d21,d22,r12inv,r21inv,q12,q21] = ...
	transf(sys,ny,nu);
 if any(any(d11))
  error('This m-file is invalid for nonzero D11');
  [trow,tcol]=size(d11);
  if min(eig(eye(tcol)-d11'*d11)) <= 0
    disp('I - D11*D11 < 0');
  end
 end
  
 b11 = b1(:,1:nd1);
 b12 = b1(:,nd1+1:nd);
 c11 = c1(1:ne1,:);
 c12 = c1(ne1+1:ne,:);
 d1111 = d11(1:ne1,1:nd1);
 d1112 = d11(1:ne1,nd1+1:nd);
 d1121 = d11(ne1+1:ne,1:nd1);
 d1122 = d11(ne1+1:ne,nd1+1:nd);

 F = -(gam*b2'/Y + c12);
 L = -(X\c2'*gam + b12);
 akh = -X*(b1*b1'+L*b12')/gam -(c1'*c1+c12'*F)*Y/gam;
 dk = zeros(nu,ny);
 ck = F*Y/M';
 bk = N\X*L;
 ak = N \ (akh - a' - X*(a+b2*F+L*c2)*Y) / M';

 if exist('gest')
  atld = a-b12*c2;
  ahat = a-b2*c12;
  out = ...
   [X*atld - gam*c2'*c2, c2'*b12'+ c12'*b2' - (X*b11*b11'+c11'*c11*Y)/gest-a'; 
   a, ahat*Y - gam*b2*b2'];
 else
  out = [X*(a+L*c2), akh; a, (a+b2*F)*Y];
 end

 vout((i-1)*2*nx+1:i*2*nx,:) = out;
   
end

if grid_pts > 1
 vout = vpck(vout,viv);
end