home *** CD-ROM | disk | FTP | other *** search
/ Liren Large Software Subsidy 10 / 10.iso / l / l455 / 2.ddi / MUTOOLS1.DI$ / HINFNORM.M < prev    next >
Encoding:
Text File  |  1993-03-11  |  2.8 KB  |  113 lines
  1. % function out = hinfnorm(sys,ttol)
  2. %
  3. %   Calculates the H_infinity norm of stable, SYSTEM
  4. %   matrices, using the Hamiltonian method, or of VARYING
  5. %   matrices, using PKVNORM. The second argument is used
  6. %   for SYSTEM matrices, and is the relative tolerance
  7. %   between the upper and lower bounds for the infinity
  8. %   norm when the search is terminated. There is an optional
  9. %   third argument, for a initial frequency guess, if desired.
  10. %   for SYSTEM matrices, OUT is a 1x3 row vector, with the
  11. %   upper bound, lower bound, and the frequency where the peak
  12. %   occurs. The default value for TTOL is 0.001.
  13. %
  14. %   See also: H2SYN, H2NORM, HINFCHK, HINFSYN, and HINFFI.
  15.  
  16. function out = hinfnorm(sys,ttol,iiloc)
  17.  if nargin == 0
  18.    disp('usage: out = hinfnorm(sys,ttol)')
  19.    return
  20.  end
  21.  if nargin == 1
  22.    ttol = 0.001;
  23.    iiloc = 5*abs(rand(1,1));
  24.  end
  25.  if nargin == 2
  26.    iiloc = 5*abs(rand(1,1));
  27.  end
  28.  [mtype,mrows,mcols,mnum] = minfo(sys);
  29.  if mtype == 'vary'
  30.    out = pkvnorm(sys);
  31.  elseif mtype == 'cons'
  32.    if nargout == 1
  33.      out = norm(sys);
  34.    else
  35.      disp(['SYSTEM is a constant, norm is ' num2str(norm(sys))])
  36.    end
  37.  else
  38.    [a,b,c,d] = unpck(sys);
  39.    if max(real(eig(a))) >= 0
  40.      error('SYSTEM has closed-right-half plane poles')
  41.      return
  42.    end
  43.    idn = eye(mnum);
  44.    [p,hesa] = hess(a);
  45.    hesb = p' * b;
  46.    hesc = c * p;
  47.    g = hesa \ hesb;
  48.    tfzero = d - hesc*g;
  49.    ntfz = norm(tfzero);
  50.    tfiiloc = (sqrt(-1)*iiloc*idn - hesa) \ hesb;
  51.    nnn = norm(hesc*tfiiloc + d);
  52.    nd = norm(d);
  53.    if nd >= ntfz
  54.      if nd >= nnn
  55.        nrma = nd;
  56.        loc = inf;
  57.      else
  58.        nrma = nnn;
  59.        loc = iiloc;
  60.      end
  61.    else
  62.      if ntfz >= nnn
  63.        nrma = ntfz;
  64.        loc = 0;
  65.      else
  66.        nrma = nnn;
  67.        loc = iiloc;
  68.      end
  69.    end
  70.    if nrma == 0
  71.      nrmt = 1e-8;
  72.    else
  73.      nrmt = (1+ttol)*nrma;
  74.    end
  75.    code = hinfchk(sys,nrmt,1e-9);
  76.    j = 0;
  77.  
  78.    maxit = 100;
  79.    while j < maxit
  80.      j = j+1;
  81.      if code == -1
  82.        if nargout == 0
  83.          disp(['norm between ' num2str(nrma) ' and ' num2str(nrmt)]);
  84.          disp(['achieved near ' num2str(loc)]);
  85.        else
  86.          out = [nrma nrmt loc];
  87.        end
  88.        j = maxit + 1;
  89.      else
  90.        if length(code) == 1
  91.          probfreq = code;
  92.        else
  93. %        probfreq = code(1:length(code)-1) + 0.5*diff(code);
  94.          probfreq = code(1:length(code)-1) + 0.5*diff(code);
  95.        end
  96.        tmp = 0;
  97.        for jj=1:length(probfreq)
  98.          g = ( (sqrt(-1)*probfreq(jj)) * idn - hesa) \ hesb;
  99.          try = norm(d + hesc*g);
  100.          if try > tmp
  101.            tmp = try;
  102.            loc = probfreq(jj);
  103.          end
  104.        end
  105.        nrma = tmp;
  106.        nrmt = (1+ttol)*nrma;
  107.        code = hinfchk(sys,nrmt,1e-9);
  108.       end
  109.     end
  110.   end
  111. %
  112. % Copyright MUSYN INC 1991,  All Rights Reserved
  113.