home *** CD-ROM | disk | FTP | other *** search
/ Liren Large Software Subsidy 10 / 10.iso / l / l455 / 7.ddi / ROBUST.DI$ / DRICCOND.M < prev    next >
Encoding:
Text File  |  1993-03-11  |  2.4 KB  |  78 lines
  1. function [tot] = driccond(a,b,q,r,p1,p2);
  2. %
  3. % [TOT] = DRICCOND(A,B,Q,R,P1,P2) produces a Riccati condition
  4. %    analysis using the concept of Arnold and Laub (1984).
  5. %
  6. %    INPUTS:
  7. %
  8. %      System: (A,B) matrices, and weighting matrices Q, R
  9. %      Riccati Solution: P1,P2 (P = P2/P1)
  10. %
  11. %    OUTPUT:
  12. %
  13. %         TOT = [norA norQ norR conr conP1 conArn conBey res]
  14. %     where                                           -1
  15. %        norA, norQ, norR ----- F-Norm of A, Q, and  BR B'
  16. %        conr ----- condition number of R
  17. %        conP1 ---- condition number of P1
  18. %        conArn --- Arnold and Laub's Riccati condition number
  19. %        conBey --- Byers's Riccati condition number
  20. %        res ------ residual of Riccati equation
  21. %
  22.  
  23. % R. Y. Chiang & M. G. Safonov 9/18/1988
  24. % Copyright (c) 1988 by the MathWorks, Inc.
  25. % All Rights Reserved.
  26. % -----------------------------------------------------------------
  27. [ns,ns] = size(a);
  28. %
  29. rc = b/r*b';
  30. p = p2/p1;
  31. perr = a'*p*a-p-a'*p*b/(rc+b'*p*b)*b'*p*a+q;
  32. %
  33. clc
  34. format short e
  35. disp(' -----------------------------------');
  36. disp('    Riccati Condition Analysis ');
  37. disp(' -----------------------------------');
  38. disp('                         -1')                   
  39. disp(' 1). F-Norm of A, Q and BR B`:');
  40. norA = norm(a,'fro');
  41. norQ = norm(q,'fro');
  42. norRc = norm(rc,'fro');
  43. norAQR = [norA norQ norRc],
  44. disp('                        hit <ret> to continue');
  45. pause
  46. %
  47. disp(' 2). Condition no. of R:');
  48. conr = cond(r),
  49. disp('                        hit <ret> to continue');
  50. pause
  51. %
  52. disp(' 3). Condition no. of P1:');
  53. conP1 = cond(p1),
  54. disp('                        hit <ret> to continue');
  55. pause
  56. %
  57. %disp(' 4). Arnold and Laub Riccati Condition no.:');
  58. %acl = inv(eye(ac)+rc*p)*ac;
  59. %sepArn = min(svd(kron(eye(acl'),acl')+kron(acl',eye(acl'))));
  60. %conArn = norm(qc,'fro')/norm(p,'fro')/sepArn,
  61. %disp('                        hit <ret> to continue');
  62. %pause
  63. %
  64. disp(' 4). Byers Riccati Condition no.:');
  65. acl = inv(eye(a)+rc*p)*a;
  66. sepBye = min(svd(kron(acl',acl')-eye(size(acl)*[1;0]*2)));
  67. conBye = (norm(q,'fro')/norm(p,'fro')+2*norm(a,'fro')^2...
  68.          +norm(a,'fro')^2*norm(rc,'fro')*norm(p,'fro'))/sepBye,
  69. disp('                        hit <ret> to continue');
  70. pause
  71. %
  72. disp(' 5). Residual = norm(perr,1)/norm(p,1)');
  73. res = norm(perr,1)/norm(p,1),
  74. %
  75. tot = [norA norQ norRc conr conP1 conBye res]';
  76. %
  77. % ------ End of RICCOND.M % RYC/MGS %
  78.