share/linearalgebra/rtest_eigens_by_jacobi.mac

   1 (listequalp(p,q) := block([listarith : true], every(lambda([x], is(abs(x) < ?flonum\-epsilon * 10000)), sort(p) - sort(q))),0);
   2 0$
   3
   4 (matrixequalp(p,q) := every(lambda([x], is(abs(x) < ?flonum\-epsilon * 10000)),p-q),0);
   5 0$
   6
   7 (m : matrix([1]),0);
   8 0$
   9
  10 eigens_by_jacobi(m,floatfield);
  11 [[1.0], matrix([1.0])];
  12
  13 eigens_by_jacobi(m,bigfloatfield);
  14 [[1.0b0], matrix([1.0b0])];
  15
  16 (m : matrix([0]),0);
  17 0$
  18
  19 eigens_by_jacobi(m,floatfield);
  20 [[0.0], matrix([1.0])];
  21
  22 eigens_by_jacobi(m,bigfloatfield);
  23 [[0.0b0], matrix([1.0b0])];
  24
  25 (m : matrix([%i]),0);
  26 0$
  27
  28 errcatch(eigens_by_jacobi(m,floatfield));
  29 []$
  30
  31 errcatch(eigens_by_jacobi(m,bigfloatfield));
  32 []$
  33
  34 (m : matrix([0,0],[0,0]),0);
  35 0$
  36
  37 eigens_by_jacobi(m,floatfield);
  38 [[0.0,0.0], matrix([1.0,0.0],[0.0,1.0])]$
  39
  40 eigens_by_jacobi(m,bigfloatfield);
  41 [[0.0b0,0.0b0], matrix([1.0b0,0.0b0],[0.0b0,1.0b0])]$
  42
  43
  44 /*-------*/
  45 (m : matrix([1,2],[2,1]),0);
  46 0$
  47
  48 (e : eigens_by_jacobi(m),0);
  49 0$
  50
  51 listequalp(first(e),[-1,3]);
  52 true$
  53
  54 matrixequalp(transpose(second(e)) . m . second(e), apply('diag_matrix, first(e)));
  55 true$
  56
  57 /*(e : eigens_by_jacobi(m,bigfloatfield),0);
  58 0$
  59
  60 listequalp(first(e),[-1,3]);
  61 true$
  62
  63 matrixequalp(transpose(second(e)) . m . second(e), apply('diag_matrix, first(e)));
  64 true$
  65 */
  66 /*-------*/
  67 (m : matrix([14,32,50],[32,77,122],[50,122,194]),0);
  68 0$
  69
  70 (e : eigens_by_jacobi(m),0);
  71 0$
  72
  73 listequalp(first(e), map('rhs, allroots(charpoly(m,z))));
  74 true$
  75
  76 matrixequalp(transpose(second(e)) . m . second(e), apply('diag_matrix, first(e)));
  77 true$
  78
  79 (remvalue(m,e),0);
  80 0$
  81
  82 (remfunction(matrixequalp, listequalp),0);
  83 0$
  84
  85
  86