Re: [Maxima-discuss] How does Maxima compute determinants ?

[Richard F. <fa...@gm...>]

 In my experience it is often critical to find out if det(M) is zero.  If
you
If you are not actually computing results by multiplying, adding (etc)
polynomials,
then you might not find out.  If you are just constructing a large
unsimplified
tree in order to save computation time, you can avoid even that work if you
just leave the result as det(M), which might be as useful, and would
certainly be
more compact!
It's an unnecessary complication to define M as a matrix of gensyms. You
can just fill it
in with  g[1,1], ... g[n,n]  by   newdet(genmatrix(g,10,10))...
On my computer the 10x10 determinant took. 2.9 seconds.  Using ratexpand(%)
to
find out how many terms there were in it ran out of memory on my system.
However, the 9x9 determinant has 362,880 terms.  9! is 362,880. I think we
can
predict the size of the determinant of the 10x10 !

This illustrates the unreality of this calculation except as a
benchmark.  It is unlikely that a human would want to see this displayed,
or that this is
in some way a valuable result to feed into another phase of computation.

 The first term of the 9x9 determinant is
g[1,1]*g[2,2]*g[3,3]*g[4,4]*g[5,5]*g[6,6]*g[7,7]*g[8,8]*g[9,9].
The last term, in the order given by ratexpand, is
-(g[1,8]*g[2,9]*g[3,7]*g[4,6]*g[5,5]*g[6,4]*g[7,3]*g[8,2]*g[9,1])

All rather predictable.  In fact, you can write a program generating, in
sequence,all the
terms in a determinant using combinatorial arrangements, and then, only if
needed, plug
in the actual matrix elements.  Thus we could write, instead of det(M), in
just those 6
characters, something like this:
   'sum(cc (M,j)  , j,1, n!)

where cc(M, j)  takes a matrix M of n^2 items and produces the jth term in
the
expanded determinant, where that term is a product of an appropriate
selection of n terms in M.
This would not answer the question of whether the determinant was zero or
not.
It could possibly be useful in that one could filter out terms that are a
product including
a matrix entry g[p,q] which is zero.  So writing the program cc() might be
of some
interest if you really want to try this out.  Of course you would not want
to expand that
sum(), first.  You would generate the terms one by one in a loop, and skip
over all terms with
a zero in the product.  This makes the cc program a bit more complicated,
but
knowing just one entry in a given matrix is zero substantially reduces the
number of terms.
Just experimenting with a general 6x6 determinant with 720 terms, one zero
at M[3,3] reduces it to
600 terms.  Still a lot of terms, but if that's what you are asking for,
there it is.
RJF




On Sun, Jan 7, 2024 at 7:45 AM Oscar Benjamin <osc...@gm...>
wrote:

> This is what I get with newdet. It crashes at 10x10:
>
> (%i1) g: lambda([i,j], gensym());
> (%o1)                      lambda([i, j], gensym())
> (%i2) showtime:true;
> Evaluation took 0.0000 seconds (0.0000 elapsed) using 0 bytes.
> (%o2)                                true
> (%i3) newdet(genmatrix(g, 6, 6))$
> Evaluation took 0.0022 seconds (0.0000 elapsed) using 447.844 KB.
> (%i4) newdet(genmatrix(g, 7, 7))$
> Evaluation took 0.0097 seconds (0.0080 elapsed) using 2.374 MB.
> (%i5) newdet(genmatrix(g, 8, 8))$
> Evaluation took 0.0715 seconds (0.0720 elapsed) using 17.447 MB.
> (%i6) newdet(genmatrix(g, 9, 9))$
> Evaluation took 0.4691 seconds (0.4680 elapsed) using 153.718 MB.
> (%i7) newdet(genmatrix(g, 10, 10))$
> Heap exhausted during garbage collection
>
> The newdet function computes the result in rat form and the rat form
> of this expression is necessarily factorially large so both runtime
> and memory consumption grow factorially.
>
> With determinant (and without copy-tree) we can go much larger:
>
> (%i3) determinant(genmatrix(g, 10, 10))$
> Evaluation took 0.0886 seconds (0.0880 elapsed) using 8.467 MB.
> (%i4) determinant(genmatrix(g, 11, 11))$
> Evaluation took 0.1761 seconds (0.1760 elapsed) using 19.808 MB.
> (%i5) determinant(genmatrix(g, 12, 12))$
> Evaluation took 0.4656 seconds (0.4640 elapsed) using 45.755 MB.
> (%i6) determinant(genmatrix(g, 13, 13))$
> Evaluation took 1.5604 seconds (1.5600 elapsed) using 114.221 MB.
> (%i7) determinant(genmatrix(g, 14, 14))$
> Evaluation took 4.9642 seconds (4.9641 elapsed) using 237.625 MB.
> (%i8) determinant(genmatrix(g, 15, 15))$
> Evaluation took 17.9706 seconds (17.9724 elapsed) using 535.027 MB.
> (%i9) determinant(genmatrix(g, 16, 16))$
> Evaluation took 71.0962 seconds (71.1095 elapsed) using 1196.405 MB.
> (%i10) determinant(genmatrix(g, 17, 17))$
> Evaluation took 316.7744 seconds (316.8670 elapsed) using 2659.191 MB.
>
> Here the growth in both runtime and memory consumption is exponential
> rather than factorial so we can go to larger matrices with more
> symbols. The key is that not expanding products like (a+b)*(c+d)
> allows a representation with fewer than factorially many *distinct*
> subexpressions.
>
> The Bird method can take this much larger. If products are not
> expanded then the method can in principle be O(n^4) in the sense that
> there are O(n^4) arithmetic operations and O(n^4) distinct
> subexpressions in the result (I think...). There are reasons why this
> is slower than it should be with SymPy expressions but in Maxima's
> case it would likely not be O(n^4) because of copy-tree.
>
> This is the Bird method with SymPy's current master branch (and one
> small change because I noticed something slowing it down
> unnecessarily):
>
> In [1]: mkmatrix = lambda n: Matrix(symbols(f'x:{n**2}')).reshape(n, n)
>
> In [2]: %time d = mkmatrix(15).det(method='bird')
> CPU times: user 3.14 s, sys: 4.3 ms, total: 3.14 s
> Wall time: 3.14 s
>
> In [3]: %time d = mkmatrix(16).det(method='bird')
> CPU times: user 4.06 s, sys: 12.4 ms, total: 4.07 s
> Wall time: 4.07 s
>
> In [4]: %time d = mkmatrix(17).det(method='bird')
> CPU times: user 5.1 s, sys: 4.35 ms, total: 5.1 s
> Wall time: 5.1 s
>
> In [5]: %time d = mkmatrix(18).det(method='bird')
> CPU times: user 6.39 s, sys: 8.19 ms, total: 6.4 s
> Wall time: 6.4 s
>
> In [6]: %time d = mkmatrix(19).det(method='bird')
> CPU times: user 8.05 s, sys: 4.16 ms, total: 8.05 s
> Wall time: 8.05 s
>
> In [7]: %time d = mkmatrix(20).det(method='bird')
> CPU times: user 9.77 s, sys: 4.06 ms, total: 9.78 s
> Wall time: 9.78 s
>
> In [8]: %time d = mkmatrix(21).det(method='bird')
> CPU times: user 12.1 s, sys: 4.11 ms, total: 12.1 s
> Wall time: 12.1 s
>
> I would expect that a Maxima implementation (without copy-tree) would
> be faster than this but note that the scaling looks at least
> approximately like n^4.
>
> --
> Oscar
>
> On Sun, 7 Jan 2024 at 15:06, Michel Talon <ta...@lp...> wrote:
> >
> > Note there exists another determinant program in maxima, called newdet
> >
> > which is notably faster than determinant.  It has been written by R.
> Fateman, and
> >
> > does some memoizing if i remember well. Of course all those determinant
> programs
> >
> > are only usable on small matrices.
> >
> >
> > Le 07/01/2024 à 14:13, Oscar Benjamin a écrit :
> >
> > I found that Maxima's built in determinant routine seemed to be faster
> > than your description here so I looked into it and it turns out that
> > it suffers from the problem with copy-tree that I mentioned in the
> > other thread "what is slow about large expressions". I was initially
> > testing with my local build that has copy-tree calls removed. If the
> > copy-tree calls are removed then it is faster and uses less memory
> > etc:
> >
> > --
> > Michel Talon
> >
> > _______________________________________________
> > Maxima-discuss mailing list
> > Max...@li...
> > https://lists.sourceforge.net/lists/listinfo/maxima-discuss
>
>
> _______________________________________________
> Maxima-discuss mailing list
> Max...@li...
> https://lists.sourceforge.net/lists/listinfo/maxima-discuss
>