home *** CD-ROM | disk | FTP | other *** search
/ Liren Large Software Subsidy 10 / 10.iso / l / l455 / 2.ddi / MUTOOLS1.DI$ / HINFSYN.M < prev    next >
Encoding:
Text File  |  1993-03-11  |  8.4 KB  |  305 lines
  1. % function [k,g,gfin,ax,ay,hamx,hamy] = 
  2. %                 hinfsyn(p,nmeas,ncon,gmin,gmax,tol,ricmethd,epr,epp) 
  3. %
  4. %  This function computes the H-infinity (sub) optimal n-state
  5. %  controller, using Glover's and Doyle's 1988 result, for a system P.
  6. %  the system P is partitioned:   
  7. %                          | a   b1   b2   | 
  8. %              p    =      | c1  d11  d12  |
  9. %                          | c2  d21  d22  | 
  10. %   where b2 has column size of the number of control inputs (NCON)
  11. %   and c2 has row size of the number of measurements (NMEAS) being
  12. %   provided to the controller.
  13. % 
  14. %   inputs:
  15. %      P        -   interconnection matrix for control design
  16. %      NMEAS    -   # controller inputs (np2)
  17. %      NCON     -   # controller outputs (nm2)
  18. %      GMIN     -   lower bound on gamma
  19. %      GMAX     -   upper bound on gamma
  20. %      TOL      -   relative difference between final gamma values
  21. %      RICMETHD -   Riccati solution via
  22. %                     1 - eigenvalue reduction
  23. %                     2 - real schur decomposition  (default)
  24. %      EPR      -   measure of when real(eig) of Hamiltonian is zero
  25. %                    (default EPR = 1e-10)
  26. %      EPP      -  positive definite determination of xinf solution
  27. %                    (default EPP = 1e-6)
  28. %    outputs:
  29. %      K       -   H-infinity controller
  30. %      G       -   closed-loop system
  31. %      GFIN    -   final gamma value used in control design
  32. %      AX      -   X-Riccati solution as a VARYING matrix with 
  33. %                     independent variable gamma (optional)
  34. %      AY      -   Y-Riccati solution as a VARYING matrix with 
  35. %                      independent variable gamma (optional)
  36. %      HAMX    -   X-Hamiltonian matrix as a VARYING matrix with 
  37. %                      independent variable gamma (optional)
  38. %      HAMY    -   Y-Hamiltonian matrix as a VARYING matrix with 
  39. %                      independent variable gamma (optional)
  40. %
  41. %   See also: H2SYN, H2NORM, HINFFI, HINFNORM, RIC_EIG, and RIC_SCHR
  42.  
  43. function [k,g,gfin,ax,ay,hamx,hamy] = ...
  44.           hinfsyn(p,nmeas,ncon,gmin,gmax,tol,imethd,epr,epp) 
  45. %
  46. %   the following assumptions are made:
  47. %    (i)   (a,b2,c2) is stabilizable and detectable
  48. %    (ii)  d12 and d21 have full rank
  49. %     on return, the eigenvalues in egs must all be > 0
  50. %     for the closed-loop system to be stable and to
  51. %     have inf-norm less than GAM.
  52. %     this program provides a gamma iteration using a
  53. %     bisection method, it iterates on the value
  54. %     of GAM to approach the H-inf optimal controller.
  55. %                                    np  nm1  nm2
  56. %   The system p is partitioned:   | a   b1   b2   | np
  57. %                                  | c1  d11  d12  | np1
  58. %                                  | c2  d21  d22  | np2
  59. % 
  60. % Reference Paper:
  61. %       'State-space formulae for all stabilizing controllers
  62. %         that satisfy and h-infinity-norm bound and relations
  63. %         to risk sensitivity,' by keith glover and john doyle,
  64. %         systems & control letters 11, oct, 1988, pp 167-172.
  65.  
  66. storxinf = [];
  67. storyinf = [];
  68. storhx = [];
  69. storhy = [];
  70. gammavecxi = [];
  71. gammavecyi = [];
  72. gammavechx = [];
  73. gammavechy = [];
  74.  
  75. if nargin == 0
  76.  disp('usage: [k,g,gfin] = hinfsyn(p,nmeas,ncon,gmin,gmax,tol,ricmethd,epr,epp) ')
  77.  return
  78. end
  79. if nargin <= 5
  80.  disp('usage: [k,g,gfin] = hinfsyn(p,nmeas,ncon,gmin,gmax,tol,ricmethd,epr,epp) ')
  81.  error('incorrect number of input arguments')
  82.  return
  83. end
  84. % setup default cases
  85. if nargin == 6
  86.   imethd = 2;
  87.   epr = 1e-10;
  88.   epp = 1e-6;
  89. elseif nargin == 7
  90.   epr = 1e-10;
  91.   epp = 1e-6;
  92. elseif nargin == 8
  93.   epp = 1e-6;
  94. end
  95.   
  96. gam2=gmax;
  97. gam0=gmax;
  98. if gmin > gmax
  99.   error(' gamma min is greater than gamma max')
  100.   return
  101. end
  102. % save the input system to form the closed-loop system 
  103. pnom = p;
  104. niter = 1;
  105. npass = 0;
  106. first_it = 1;
  107. totaliter = 0;
  108. % fprintf('checking rank conditions: ');
  109.   [p,r12,r21,fail] = hinf_st(p,nmeas,ncon);
  110.   if fail == 1,
  111.     fprintf('\n')
  112.     return
  113.   end
  114.  
  115. %            gamma iteration
  116.  
  117. fprintf('Test bounds:  %10.4f <  gamma  <=  %10.4f\n\n',gmin,gmax)
  118. fprintf('  gamma    hamx_eig  xinf_eig  hamy_eig   yinf_eig   nrho_xy   p/f\n')
  119.  
  120.  while niter==1
  121.   totaliter = totaliter + 1;
  122.   if first_it ==1        % first iteration checks gamma max
  123.     gam = gmax;
  124.   else
  125.     gam=((gam2+gam0)/2.);
  126.     gamall=[gam0 gam gam2];
  127.   end
  128.   if gam < 1000
  129.     fprintf('%9.3f  ',gam)
  130.   else
  131.     fprintf('%9.3e  ',gam)
  132.   end
  133.   [xinf,yinf,f,h,fail,hmx,hmy,hxmin,hymin] = ... 
  134.                         hinf_gam(p,nmeas,ncon,epr,gam,imethd);
  135.   if nargout >= 4
  136.     if isempty(xinf)
  137.       % don't update the matrix
  138.     else
  139.       gammavecxi = [gammavecxi ; gam];
  140.       storxinf = [storxinf ; xinf];
  141.     end
  142.     if nargout >= 5
  143.       if isempty(yinf)
  144.         % don't update the matrix
  145.       else
  146.         gammavecyi = [gammavecyi ; gam];
  147.         storyinf = [storyinf ; yinf];
  148.       end
  149.       if nargout >= 6
  150.         if isempty(hmx)
  151.           % don't update the matrix
  152.         else
  153.           gammavechx = [gammavechx ; gam];
  154.           storhx = [storhx ; hmx];
  155.         end
  156.         if nargout >= 7
  157.           if isempty(hmy)
  158.             % don't update the matrix
  159.           else
  160.            gammavechy = [gammavechy ; gam];
  161.            storhy = [storhy ; hmy];
  162.           end
  163.         end
  164.       end
  165.     end
  166.   end
  167. %
  168. % check the conditions for a solution to the h-infinity control problem
  169. %
  170.   if fail == 0
  171.     xe = eig(xinf);
  172.     ye = eig(yinf);
  173.     xye = eig(xinf*yinf);
  174.     xeig = min(real(xe));
  175.     yeig = min(real(ye));
  176.     xyinfe = max(abs(xye));
  177.     index = find(real(xe)==xeig);
  178.     mprintf(hxmin,'%8.1e',' ');
  179.     mprintf(real(xe(index(1))),' %8.1e',' ');
  180.     index = find(real(ye)==yeig);
  181.     mprintf(hymin,' %8.1e ',' ');
  182.     mprintf(real(ye(index(1))),' %8.1e',' ');
  183.  
  184. % scaled rho(xy)
  185.     fprintf('  %7.4f   ',xyinfe/gam/gam)
  186.     if xeig < -epp | yeig < -epp | xeig == nan | yeig == nan;
  187.       npass=1;
  188.     end
  189.     if xyinfe/gam/gam > 1
  190.       npass=1;
  191.     end
  192.   else
  193.     if isempty(xinf)
  194.       mprintf(hxmin,'%8.1e',' ');
  195.       fprintf('  ******* ')
  196.       if isempty(yinf)
  197.         mprintf(hymin,' %8.1e ',' ');
  198.         fprintf('  ******* ')
  199.         fprintf('   ******   ')
  200.       else
  201.         ye = eig(yinf);
  202.         mprintf(hymin,' %8.1e ',' ');
  203.         mprintf(min(real(ye)),' %8.1e',' ');
  204.         fprintf('   ******   ')
  205.       end
  206.     else
  207.       xe = eig(xinf);
  208.       mprintf(hxmin,'%8.1e',' ');
  209.       mprintf(min(real(xe)),' %8.1e ',' ');
  210.       mprintf(hymin,' %8.1e ',' ');
  211.       fprintf(' ******* ')
  212.       fprintf('   ******   ')
  213.     end
  214.     npass = 1;
  215.   end
  216.  
  217.   if npass == 1
  218.     fprintf(' f \n')
  219.   else
  220.     fprintf(' p \n')
  221.   end
  222. %
  223. % check the gamma iteration
  224. %
  225.  if first_it == 1;
  226.    if npass == 0;
  227.      gam0 = gmin; 
  228.      gam2 = gmax; 
  229.      xinffin = xinf; 
  230.      yinffin = yinf;
  231.      fsave = f;
  232.      hsave = h; 
  233.      gfin = gmax;
  234.      first_it = 0;
  235.      if gmin == gmax
  236.        niter = 0;           % stop since gmax = gmin
  237.      end
  238.   else
  239.      niter = 0;
  240.      fprintf('Gamma max, %10.4f, is too small !!\n',gmax);
  241. % save output data
  242.      if nargout >= 4
  243.        ax = vpck(storxinf,gammavecxi);
  244.        if nargout >= 5
  245.          ay = vpck(storyinf,gammavecyi);
  246.          if nargout >= 6
  247.            hamx = vpck(storhx,gammavechx);
  248.            if nargout >= 7
  249.              hamy = vpck(storhy,gammavechy);
  250.         end
  251.          end
  252.        end
  253.      end
  254.      return
  255.    end    
  256.  else
  257.    if abs(gam - gam0) < tol
  258.      niter = 0;
  259.      if npass == 0
  260.        gfin = gam;
  261.        xinffin = xinf;
  262.        yinffin = yinf;
  263.      fsave = f;
  264.      hsave = h; 
  265.      else
  266.        gfin = gam2; 
  267.      end
  268.    else
  269.      if npass == 0;
  270.        gam2 = gam;
  271.        xinffin = xinf;
  272.        yinffin = yinf;
  273.      fsave = f;
  274.      hsave = h; 
  275.      else
  276.        gam0 = gam;
  277.      end
  278.    end
  279.    npass = 0;
  280.  end
  281. end
  282. xinf = xinffin; 
  283. yinf = yinffin; 
  284. %  finished with gamma iteration
  285. fprintf('\n Gamma value achieved: %10.4f\n',gfin)
  286.  
  287. [k] = hinf_c(p,nmeas,ncon,xinf,yinf,fsave,hsave,gfin,r12,r21);
  288. [g] = starp(pnom,k,nmeas,ncon);
  289.  
  290. % save output data
  291.   if nargout >= 4
  292.     ax = vpck(storxinf,gammavecxi);
  293.     if nargout >= 5
  294.       ay = vpck(storyinf,gammavecyi);
  295.       if nargout >= 6
  296.         hamx = vpck(storhx,gammavechx);
  297.         if nargout >= 7
  298.           hamy = vpck(storhy,gammavechy);
  299.         end
  300.       end
  301.     end
  302.   end
  303. %
  304. % Copyright MUSYN INC 1991,  All Rights Reserved
  305.