home *** CD-ROM | disk | FTP | other *** search
/ The Datafile PD-CD 5 / DATAFILE_PDCD5.iso / utilities / r / rlab / CTB / lyapcheck < prev    next >
Encoding:
Text File  |  1995-11-15  |  1.6 KB  |  69 lines

  1.  
  2. rfile lyap2
  3. rfile lyap
  4.  
  5. A=rand(4,4);
  6.  
  7. E=eig(A);
  8. Eval=[-4.0,-3.0,-2.0,-5.0];
  9. Anew=conj(E.vec')*diag(Eval)*E.vec;
  10. A=real(Anew);
  11. Aold=A;
  12.  
  13. B=rand(4,4);
  14. E=eig(B);
  15. Eval=[-12,-3.4,-5.8,-1.3];
  16. Bnew=conj(E.vec')*diag(Eval)*E.vec;
  17. B=real(Bnew);
  18. Bold=B;
  19.  
  20. C=rand(4,4);
  21. E=eig(C);
  22. Eval=[-0.4,-3.4,-9.5,-6.3];
  23. Cnew=conj(E.vec')*diag(Eval)*E.vec;
  24. C=real(Cnew);
  25. Cold=C;
  26.  
  27. #####################################################################
  28. printf ("\nStart A,B checks\n");
  29. Xab2=lyap2(A,B);
  30. printf ("The following checks ought to be TRUE (1)\n");
  31. Bcheck = all (all (Bold == B))    # these ought to be identical (TRUE)
  32. Acheck = all( all (Aold == A))
  33.  
  34. printf ("lyap2 solution self-check (this matrix ought to be SMALL\n");
  35. (A*Xab2 + Xab2*A') + B
  36.  
  37. printf ("The following checks ought to be TRUE (1)\n");
  38. Xab=lyap(A,B);
  39. Bcheck = all (all (Bold == B))
  40. Acheck = all (all (Aold == A))
  41.  
  42. printf ("lyap solution self-check (this matrix ought to be SMALL\n");
  43. (A*Xab + Xab*A') + B
  44.  
  45. printf("%s","More A,B Solution checks.\n");
  46. Xnorm=norm(Xab-Xab2)
  47. printf ("\nCheck lyap against lyap2 (X - X2)\n");
  48. Xab-Xab2
  49.  
  50. #####################################################################
  51. printf ("\nStart A,B,C checks\n");
  52. Xabc2=lyap2(A,B,C);
  53. printf ("lyap2 solution self-check (this matrix ought to be SMALL\n");
  54. (A*Xabc2 + Xabc2*B) + C
  55.  
  56. Xabc=lyap(A,B,C);
  57. printf ("\nlyap solution self-check (this matrix ought to be SMALL\n");
  58. (A*Xabc + Xabc*B) + C
  59.  
  60. printf("%s","\nMore A,B,C Solution checks.\n");
  61. Xnorm=norm(Xabc-Xabc2)
  62. printf ("\nCheck lyap against lyap2 (X - X2)\n");
  63. Xabc-Xabc2
  64.  
  65. Bcheck = all( all (Bold == B))    # these ought to be identical
  66. Acheck = all (all (Aold == A))
  67. Ccheck = all (all (Cold == C))
  68.  
  69.