----Fred J. wrote: -----
>I was expecting LU decomposition, ...
(%i2) load("linalg");
Warning - you are redefining the MACSYMA function RANK
(%o2) c:/maxima/linearalgebra/linalg.mac
(%i3) m : matrix([1,2,3],[4,5,6],[7,8,9])$
(%i4) lu_factor(m);
(%o4) [matrix([1,2,3],[4,-3,-6],[7,2,0]),[1,2,3]] <-- a packed form
(%i5) get_lu_factors(%);
(%o5)
[matrix([1,0,0],[0,1,0],[0,0,1]),matrix([1,0,0],[4,1,0],[7,2,1]),matrix([1,2,3],[0,-3,-6],[0,0,0])]
(%i6) %[1] . %[2] . %[3]; <--- let's check
(%o6) matrix([1,2,3],[4,5,6],[7,8,9])
(%i8) lu_factor(m,floatfield)$ <--- use double floats
(%i9) get_lu_factors(%);
(%o9)
[matrix([0,1,0],[0,0,1],[1,0,0]),matrix([1,0,0],[0.14285714285714285,1,0],
[0.5714285714285714,0.50000000000000022,1]),matrix([7.0,8.0,9.0],
[0,0.85714285714285721,1.7142857142857144],[0,0,0.0])]
(%i10) %[1] . %[2] . %[3];
(%o10) matrix([1.0,2.0,3.0],[4.0,5.0,6.0],[7.0,8.0,9.0])
(%i13) lu_factor(m,rationalfield)$ <--- factor using Common Lisp rational
numbers
(%i14) get_lu_factors(%);
(%o14)
[matrix([0,1,0],[0,0,1],[1,0,0]),matrix([1,0,0],[1/7,1,0],[4/7,1/2,1]),matrix([7,8,9],[0,6
/7,12/7],[0,0,0])]
If your Maxima doesn't have linalg, get it from
http://www.unk.edu/facstaff/profiles/willisb/.
Barton