home *** CD-ROM | disk | FTP | other *** search
/ Liren Large Software Subsidy 10 / 10.iso / l / l455 / 10.ddi / CONTROL.DI$ / LYAP2.M < prev    next >
Encoding:
Text File  |  1993-03-11  |  1.2 KB  |  49 lines
  1. function X = lyap2(A, B, C)
  2. %LYAP2    Lyapunov equation solution using eigenvalue decomposition. 
  3. %    X = LYAP2(A,C) solves the special form of the Lyapunov matrix 
  4. %    equation:
  5. %
  6. %        A*X + X*A' = -C
  7. %
  8. %    X = LYAP2(A,B,C) solves the general form of the Lyapunov matrix
  9. %    equation:
  10. %
  11. %        A*X + X*B = -C
  12. %
  13. %    LYAP2 is faster and generally more accurate than LYAP except when 
  14. %    A or B have repeated roots. 
  15. %
  16. %    See also DLYAP.
  17.  
  18. %    A.C.W. Grace  10-25-89
  19. %    Copyright (c) 1986-93 by The MathWorks, Inc.
  20.  
  21. [ma,na] = size(A);
  22. [mb,nb] = size(B);
  23. if ((ma ~= na) | (mb ~= nb))
  24.      error('Dimensions of A and B do not agree.');
  25. elseif ma==0,
  26.     X = []; return, end
  27. end
  28.  
  29. % Perform eigenvalue  decomposition on A and B 
  30. [T,DA]=eig(A);
  31. if nargin==3
  32.     [mc,nc] = size(C);
  33.     if ((mc ~= ma) | (nc ~= mb)), error('Dimensions of C do not agree'); end
  34.     [U,DB]=eig(B);
  35.     X=-T*(T\C*U./(DA*ones(ma,mb)+ones(ma,mb)*DB))/U;
  36. else
  37. % Shortcut for AX+XA'+C=0
  38. % Note B is really C for two input arguments
  39.     IT=inv(T);  
  40.     DA=DA*ones(ma,ma);
  41.     X=-T*(IT*B*IT.'./(DA+DA.'))*T.';
  42.     C=0;
  43. end
  44.     
  45. % ignore complex part if real inputs (better be small)
  46.     if ~(any(any(imag(A))) | any(any(imag(B))) | any(any(imag(C))))
  47.         X = real(X);
  48.     end
  49.