Liren Large Software Subsidy 10/10.iso/l/l455/2.ddi/MUTOOLS1.DI$/MU.M

Text File | 1993-03-11 | 5.0 KB | 175 lines

% function [bnds,dvec,sens,pvec] = mu(matin,blk,opt)
%
%   Calculates structured singular value (mu)
%   of a VARYING matrix
%    inputs:
%    MATIN - input matrix, CONSTANT/VARYING
%    BLK   - block structure
%    OPT   - character string describing the
%        desired computation.
%                 'l' lower bound, power iteration
%                 't' increase iterations (to 300) in lower bound
%                 'r' start with D=I for EACH independent variable
%                 'c' upper bound, subgradient method
%                 'o' upper bound, modified Osborne's method
%                 'p' upper bound, Perron vector method
%                 'u' upper bound, depending on the block
%             structure, Osborne or Perron method
%          's' suppress printing information about calculation
%          'L' only compute lower bound
%          'R' start lower bound with RANDOM vectors
%        DEFAULT is 'lu'
%
%   See also: BLKNORM, DYPERT, H2NORM, HINFNORM, NORM, RANDEL,
%             UNWRAPD, UNWRAPP, VNORM, VRHO, and VSVD.

function [bnds,rvdl,sens,rvpert,pinfos] = mu(matin,blk,opt)
  if nargin < 1
    disp('usage: [bnds,dvec,sens,pvec] = mu(matin,blk,opt)')
    return
  end

  [mtype,mrows,mcols,mnum] = minfo(matin);

  if nargin < 2
    if mrows==mcols
        opt = 'lu';
        blk = ones(mrows,2);
    else
        disp('usage: [bnds,dvec,sens,pvec] = mu(matin,blk,opt)')
        disp('Nonsquare MATIN must have BLK specified')
        return
    end % if mrows
   elseif nargin == 2
    opt = 'lu';
   else
    opt = [opt 'lu'];
   end %if nargin
   initreset = 0;
   if any(opt)=='r'
     initreset = 1;
   end

%
  if mtype == 'syst'
    disp('MU for SYSTEM matrices is not implemented yet')
    return
  end
% set up pointers for block structure
  [nblk,dum] = size(blk);
  for i=1:nblk
    if blk(i,2) == 0 & blk(i,1) == 1
      blk(i,2) = 1;
    end
  end
  [blkp,repp,fulp,wd,wpert,nrd,ncd,rnblk,rblk,rblkp] = ptrs(blk);
% see PTRS.M for definitions of these variables
  if nrd ~= mcols | ncd ~= mrows
    error('MATIN dimensions incompatible with BLK dimensions')
    return
  end
  if strcmp(mtype,'cons')
%   CONSTANT matrix
    [bnds,dleft,drighti,rvpert,sens,pinfos] = ...
        mu_comp(matin,nblk,blk,blkp,wpert,repp,fulp,rnblk,rblk,rblkp,opt);
    [rvdl,condd] = dscl_vec(dleft,nblk,blk,blkp);
    rvdl(1:wd) = rvdl(1:wd)/rvdl(wd);
  elseif strcmp(mtype,'vary')
%   VARYING matrix - initialize output matrices
    initdl = eye(mrows);
    initdri = eye(mcols);
    ff = (mnum+1)*mrows;
    pp = 1:mrows:ff;
    ppm1 = pp(2:mnum+1) - 1;
%
    rvdl = zeros(mnum+1,wd+1);
    rvdl(1:mnum,wd+1) = matin(1:mnum,mcols+1);
    rvdl(mnum+1,wd+1) = inf;
    rvdl(mnum+1,wd) =mnum;
%
    bnds = zeros(mnum+1,3);
    bnds(1:mnum,3) = matin(1:mnum,mcols+1);
    bnds(mnum+1,3) = inf;
    bnds(mnum+1,2) = mnum;
%
    sens = zeros(mnum+1,nblk+1);
    sens(1:mnum,nblk+1) = matin(1:mnum,mcols+1);
    sens(mnum+1,nblk+1) = inf;
    sens(mnum+1,nblk) = mnum;
%
%   if length(opt) > 1
    if any(opt=='l') | any(opt=='t')
      rvpert = zeros(mnum+1,wpert+1);
      rvpert(1:mnum,wpert+1) = matin(1:mnum,mcols+1);
      rvpert(mnum+1,wpert+1) = inf;
      rvpert(mnum+1,wpert) = mnum;
      pinfos = bnds;
    end
%   main loop on varying matrix
    if nblk == 2 & length(repp) == 0
      idd = 0;
    end
    count = 0;
    if ~any(opt=='s')
      fprintf('points completed....\n')
    end
    for i = 1:mnum
      tmp = matin(pp(i):ppm1(i),1:mcols);
      if nblk == 2 & length(repp) == 0
        [up,dl,dri,sen,rowp,lb,err] = twobl(tmp,blk,idd);
        if err == 0
          idd = log(dl(1,1));
          bnd = [up up];
          pinfo = [err err];
        else
          idd = 0;
          bnd = [up lb];
          pinfo = [err err];
        end
      else
        [bnd,dl,dri,rowp,sen,pinfo] = ...
             mu_comp(tmp,nblk,blk,blkp,wpert,repp,fulp,...
             rnblk,rblk,rblkp,opt,initdl,initdri);
      end
      count = count + 1;
      if count < 18
        if ~any(opt=='s')
          fprintf([int2str(i) '.'])
        end
      else
        if ~any(opt=='s')
          fprintf('\n')
          fprintf([int2str(i) '.'])
        end
        count = 0;
      end
      if any(opt=='l') | any(opt=='t')
        rvpert(i,1:wpert) = rowp;
        pinfos(i,1:2) = pinfo;
      end
      [rvdl(i,1:wd),condd] = dscl_vec(dl,nblk,blk,blkp);
      if condd > 1e6
        disp(['WARNING: D scaling matrices are poorly conditioned...'])
        disp(['check upper and lower bound results carefully'])
      end
      bnds(i,1:2) = bnd;
      if initreset == 1
        initdl = eye(mrows);
        initdri = eye(mcols);
      else
        initdl = dl/rvdl(i,wd);
        initdri = dri*rvdl(i,wd);
      end
      sens(i,1:nblk) = sen;
    end % for i
    if ~any(opt=='s')
      fprintf('\n')
    end
  end

  if ~any(opt=='l') & ~any(opt=='t') & ~any(opt=='L')
    bnds = sel(bnds,1,1);
  end
%
% Copyright MUSYN INC 1991,  All Rights Reserved