share/matrix/nchrpl.mac

   1
   2 /* MATTRACE computes the trace [sum of the elements on the main diagonal] of
   3 a square matrix.  It is used by NCHARPOLY, an alternative to MACSYMA's
   4 CHARPOLY.  NCHARPOLY works by computing traces of powers the given matrix,
   5 which are known to be equal to sums of powers of the roots of the
   6 characteristic polynomial.  From these quantities the symmetric functions of
   7 the roots can be calculated, which are nothing more than the coefficients of
   8 the characteristic polynomial.  CHARPOLY works by forming the determinant of
   9 VAR * IDENT [N] - A.  Thus NCHARPOLY wins, for example, in the case of large
  10 dense matrices filled with integers, since it avoids polynomial arithmetic
  11 altogether.
  12
  13 Programmed by DRB@MC. */
  14
  15 mattrace(a) := block([acc : 0, i],
  16   if matrixp(a) then (
  17     i : lmin(matrix_size(a)),
  18     while i > 0 do (
  19       acc : acc + inpart(a,i,i),
  20       i : i - 1),
  21     acc)
  22   else funmake('mattrace, [a]));
  23
  24 ncharpoly(a,var):=
  25   block
  26     ([ak:a,trlist:[mattrace(a)],
  27       symlist:[1],
  28       k:0,p:0,
  29       maperror:false,mapprint:false],
  30      thru length(a)-1 do block([],ak:a.ak,trlist:cons(mattrace(ak),trlist)),
  31      trlist:reverse(trlist),
  32       /*
  33      MAP(LAMBDA([X],K:K+1,SYMLIST:
  34      CONS(APPLY("+",MAPLIST("*",SYMLIST,TRLIST))/-K,SYMLIST)),TRLIST),
  35     */
  36      for i in trlist do (k:k+1,symlist:
  37        cons(apply("+",maplist("*",symlist,trlist))/-k,symlist)),
  38      for i:0 unless symlist=[] do
  39        block([],p:p+first(symlist)*var^i,symlist:rest(symlist)),
  40      ratsimp(p,var));
  41
  42