On Feb 3, 10:52 pm, John H Palmieri <[email protected]> wrote:
Hi John,
> Almost two years ago, Linbox's implementation of Smith normal form was
> taken out of Sage because it was too buggy. After some work, I managed
> to reinstate it, hoping that the bugs might have been fixed. Here's a
> partial status report:
>
> 1. I haven't tested it very broadly, but it seems to work on all of
> the matrices that I feed it -- it agrees with Pari, anyway.
>
> 2. The Linbox implementation and the Pari implementation differ by a
> transpose; that is, if you really want the two to agree on some matrix
> m, you need to compare
>
> m.elementary_divisors(algorithm='pari')
>
> to
>
> m.transpose().elementary_divisors(algorithm='linbox').
>
> Actually, I just stuck the transpose into the elementary_divisors
> function...
Cool.
> 3. The timing is a bit odd. For a 105x455 sparse integer matrix m2,
> this is what happened on my iMac with 2GB RAM:
>
> sage: timeit('m2.dense_matrix().elementary_divisors
> (algorithm="pari")')
> 5 loops, best of 3: 431 ms per loop
> sage: timeit('m2.dense_matrix().elementary_divisors
> (algorithm="linbox")')
> 5 loops, best of 3: 365 ms per loop
>
> For a 719 x 1347 sparse integer matrix m6:
>
> sage: timeit("m6.dense_matrix().elementary_divisors
> (algorithm='pari')", number=1)
> 1 loops, best of 3: 23.9 s per loop
> sage: timeit("m6.dense_matrix().elementary_divisors
> (algorithm='linbox)", number=1)
> 1 loops, best of 3: 397 s per loop
>
> One comment here: although the first timing says "23.9 s per loop", it
> took about an hour for this to actually complete.
Well, it says "1 loops, best of 3: 23.9 s per loop", so I would used
time not timeit to measure performance. Since it takes many seconds at
least for those examples it is pointless anyway to run them in a loop.
We might cache the pari result, but not the LinBox one, so that might
explain the discrepancies.
> On sage.math, I got timings which are consistent with this (except
> that here, the timings accurately reflect the wall time for the
> command). Here, m7 is a 694 x1300 sparse integer matrix:
>
> sage: timeit("m7.dense_matrix().elementary_divisors
> (algorithm='pari')", number=1)
> 1 loops, best of 3: 14.7 s per loop
> sage: timeit("m7.dense_matrix().elementary_divisors
> (algorithm='linbox')", number=1)
> 1 loops, best of 3: 219 s per loop
>
> So Linbox is faster for smaller matrices and Pari might be faster on
> larger ones, depending on your computer.
I don't think that is possible, but you might be seeing some artifact.
How do you create those matrices? Maybe they have a special structure?
> As far as I can tell, neither
> has an implementation for sparse matrices which is accessible from
> Sage. Does anyone know if Pari or Linbox have code explicitly for
> Smith normal form for sparse matrices? It might be fun to see what
> such code could do with these matrices...
>
> Anyway, I'm not sure what conclusions to draw from this, except that I
> should do my computations on sage.math :)
:)
> Finally, I was unable to do this successfully on the old ubuntu box in
> my office: calling the linbox version of elementary_divisors on just a
> 2 x 3 matrix led to a RuntimeError. I have no idea how to fix this...
>
> John
Cheers,
Michael
--~--~---------~--~----~------------~-------~--~----~
To post to this group, send email to [email protected]
To unsubscribe from this group, send email to
[email protected]
For more options, visit this group at http://groups.google.com/group/sage-devel
URLs: http://www.sagemath.org
-~----------~----~----~----~------~----~------~--~---