[sage-devel] Re: p-adic precision loss in charpoly
Instead I worked around this by computing the determinants mod 2 and mod 13 and using CRT (if the determinants were both units). The time was then almost trivial. Suppose I replace this problem over ZZ/25ZZ or ZZ/256ZZ. I would still hope that the problem would NOT be lifted to ZZ for computation, since this would certainly not terminate in reasonable time for a dense matrix. In general doing this via the CRT requires factoring, so some cutoff must be made, which will become outdated... Nick --~--~-~--~~~---~--~~ To post to this group, send email to sage-devel@googlegroups.com To unsubscribe from this group, send email to sage-devel-unsubscr...@googlegroups.com For more options, visit this group at http://groups.google.com/group/sage-devel URLs: http://www.sagemath.org -~--~~~~--~~--~--~---
[sage-devel] Re: p-adic precision loss in charpoly
On Fri, Mar 20, 2009 at 10:18 AM, Nick Alexander ncalexan...@gmail.com wrote: Instead I worked around this by computing the determinants mod 2 and mod 13 and using CRT (if the determinants were both units). The time was then almost trivial. Suppose I replace this problem over ZZ/25ZZ or ZZ/256ZZ. I would still hope that the problem would NOT be lifted to ZZ for computation, since this would certainly not terminate in reasonable time for a dense matrix. In general doing this via the CRT requires factoring, so some cutoff must be made, which will become outdated... Writing fast code that works for a range of inputs almost always involves making one more explicit cutoffs, and every single one that is ever made will eventually be non-optimal in the future, or even on differently configured computers now. I just spent some time making the linear algebra over ZZ and QQ 5-250 times faster in lots of cases by special casing small dimensions to much more quickly use PARI (etc.), then switching to some asymptotically fast native sage algorithm for big input. This literally speeds up a lot of code by a factor of 5-250 that uses small linear algebra. I just hardcoded some defaults for now -- e.g., use PARI for dets of matrices with less than 50 rows. Of course, it will be much better to factor out all this tuning info in a global object that can be optimized for a range of hardware and be updated over time (sort of like a Sage version of ATLAS). That said, even with hardcoded values, doing this special casing is *dramatically* -- I mean *dramatically* -- better than not doing anything. It also makes it much easier to factor out the tuning object, when we do that later. Bill Hart -- you might want to make some remarks at this point, since you are constantly dealing with such tuning issues in MPIR and FLINT. -- William --~--~-~--~~~---~--~~ To post to this group, send email to sage-devel@googlegroups.com To unsubscribe from this group, send email to sage-devel-unsubscr...@googlegroups.com For more options, visit this group at http://groups.google.com/group/sage-devel URLs: http://www.sagemath.org -~--~~~~--~~--~--~---
[sage-devel] Re: p-adic precision loss in charpoly
Tuning is always difficult, especially when the search space is multidimensional. FLINT 2.0 will come with a tuning facility which will tune all linear cutoffs at build time. I'm afraid I don't have anything to offer for multidimensional cutoffs. Usually the search space is just too large to have an efficient tuning implementation. Something I do notice is that multiprecision arithmetic tends to need less adjustments to its tuning than stuff over Z/pZ for word sized p for example, or for single precision stuff. Trying to find a single algorithm that doesn't require some tuning cutoffs is almost certainly impossible, if you want maximum performance. Bill. 2009/3/20 William Stein wst...@gmail.com: On Fri, Mar 20, 2009 at 10:18 AM, Nick Alexander ncalexan...@gmail.com wrote: Instead I worked around this by computing the determinants mod 2 and mod 13 and using CRT (if the determinants were both units). The time was then almost trivial. Suppose I replace this problem over ZZ/25ZZ or ZZ/256ZZ. I would still hope that the problem would NOT be lifted to ZZ for computation, since this would certainly not terminate in reasonable time for a dense matrix. In general doing this via the CRT requires factoring, so some cutoff must be made, which will become outdated... Writing fast code that works for a range of inputs almost always involves making one more explicit cutoffs, and every single one that is ever made will eventually be non-optimal in the future, or even on differently configured computers now. I just spent some time making the linear algebra over ZZ and QQ 5-250 times faster in lots of cases by special casing small dimensions to much more quickly use PARI (etc.), then switching to some asymptotically fast native sage algorithm for big input. This literally speeds up a lot of code by a factor of 5-250 that uses small linear algebra. I just hardcoded some defaults for now -- e.g., use PARI for dets of matrices with less than 50 rows. Of course, it will be much better to factor out all this tuning info in a global object that can be optimized for a range of hardware and be updated over time (sort of like a Sage version of ATLAS). That said, even with hardcoded values, doing this special casing is *dramatically* -- I mean *dramatically* -- better than not doing anything. It also makes it much easier to factor out the tuning object, when we do that later. Bill Hart -- you might want to make some remarks at this point, since you are constantly dealing with such tuning issues in MPIR and FLINT. -- William --~--~-~--~~~---~--~~ To post to this group, send email to sage-devel@googlegroups.com To unsubscribe from this group, send email to sage-devel-unsubscr...@googlegroups.com For more options, visit this group at http://groups.google.com/group/sage-devel URLs: http://www.sagemath.org -~--~~~~--~~--~--~---
[sage-devel] Re: p-adic precision loss in charpoly
William, Thanks for pointing out that I was wrong about the lifting problem and original of this slow behavior. Rather than just putting in place the solution by lifting to ZZ, which you show to be much faster for reasonable size, I hope someone can profile arithmetic or linear algebra for ZZ/nZZ and figure out why this determinant is so dog slow. Are there too many GCD's or what is going on? I hope the original p-adic problem also gets a track rather than being lost by my side comment, since there the problem was that a wrong answer was being returned -- a much more serious issue! Cheers, David --~--~-~--~~~---~--~~ To post to this group, send email to sage-devel@googlegroups.com To unsubscribe from this group, send email to sage-devel-unsubscr...@googlegroups.com For more options, visit this group at http://groups.google.com/group/sage-devel URLs: http://www.sagemath.org -~--~~~~--~~--~--~---
[sage-devel] Re: p-adic precision loss in charpoly
On Fri, Mar 20, 2009 at 2:22 PM, David Kohel drko...@gmail.com wrote: William, Thanks for pointing out that I was wrong about the lifting problem and original of this slow behavior. Rather than just putting in place the solution by lifting to ZZ, which you show to be much faster for reasonable size, I hope someone can profile arithmetic or linear algebra for ZZ/nZZ and figure out why this determinant is so dog slow. Are there too many GCD's or what is going on? It's slow becase it's using the O(n!) expansion by minors algorithm, since one isn't over a domain and hence Sage doesn't use Hessenberg form. William --~--~-~--~~~---~--~~ To post to this group, send email to sage-devel@googlegroups.com To unsubscribe from this group, send email to sage-devel-unsubscr...@googlegroups.com For more options, visit this group at http://groups.google.com/group/sage-devel URLs: http://www.sagemath.org -~--~~~~--~~--~--~---
[sage-devel] Re: p-adic precision loss in charpoly
Hi, The algorithm used doesn't need to be specific to the p-adics, but shouldn't apply an algorithm for fields or integral domains. Rather than the matrix code, it looks like this is the source of the problem: sage: M.parent().base_ring() 101-adic Ring with capped relative precision 2 sage: R.is_integral_domain() True With such a result, the wrong linear algebra is applied. One could possibly identify which algorithms to apply for real complex fields, p-adics, and power series rings from this: sage: R.is_exact() False Inexact rings are only approximations to mathematical objects. However, in this case the intended ring is the quotient ring Z/p^nZ (= Z_p/p^n Z_p), which is 'exact' (i.e. a well-defined ring) but is not an integral domain. I would hope that the matrix code would then pick up a valid algorithm to apply, even if not the most efficient. --David --~--~-~--~~~---~--~~ To post to this group, send email to sage-devel@googlegroups.com To unsubscribe from this group, send email to sage-devel-unsubscr...@googlegroups.com For more options, visit this group at http://groups.google.com/group/sage-devel URLs: http://www.sagemath.org -~--~~~~--~~--~--~---
[sage-devel] Re: p-adic precision loss in charpoly
Hi guys, Thanks for looking into this. I ended up working around the problem by lifting to Z and doing the charpoly there (I know in advance the output precision, and the complexity is not that important for my application, so it's no big deal). I've put up a patch for review here: http://sagetrac.org/sage_trac/ticket/793 (wink wink nudge nudge) david On Mar 19, 10:00 am, David Kohel drko...@gmail.com wrote: Hi, The algorithm used doesn't need to be specific to the p-adics, but shouldn't apply an algorithm for fields or integral domains. Rather than the matrix code, it looks like this is the source of the problem: sage: M.parent().base_ring() 101-adic Ring with capped relative precision 2 sage: R.is_integral_domain() True With such a result, the wrong linear algebra is applied. One could possibly identify which algorithms to apply for real complex fields, p-adics, and power series rings from this: sage: R.is_exact() False Inexact rings are only approximations to mathematical objects. However, in this case the intended ring is the quotient ring Z/p^nZ (= Z_p/p^n Z_p), which is 'exact' (i.e. a well-defined ring) but is not an integral domain. I would hope that the matrix code would then pick up a valid algorithm to apply, even if not the most efficient. --David --~--~-~--~~~---~--~~ To post to this group, send email to sage-devel@googlegroups.com To unsubscribe from this group, send email to sage-devel-unsubscr...@googlegroups.com For more options, visit this group at http://groups.google.com/group/sage-devel URLs: http://www.sagemath.org -~--~~~~--~~--~--~---
[sage-devel] Re: p-adic precision loss in charpoly
Hi, One can always work around a problem, but it would be better to discuss the best approach and put this in place. A related problem I had recently was in finding a random element of GL_n(ZZ/26ZZ) where n was 20-30. It was failing to terminate in the determinant computation. My guess is that a determinant over ZZ was being computed and reduced but that the resulting determinant was too big. I didn't verify this, but invite someone to check. Instead I worked around this by computing the determinants mod 2 and mod 13 and using CRT (if the determinants were both units). The time was then almost trivial. Suppose I replace this problem over ZZ/25ZZ or ZZ/256ZZ. I would still hope that the problem would NOT be lifted to ZZ for computation, since this would certainly not terminate in reasonable time for a dense matrix. Taking the approach of lifting to ZZ for charpolys of matrices of non-trivial size will undoubtably lead to exactly the same coefficient explosion which is the presumed source of the problem with determinants over ZZ/26ZZ. --David --~--~-~--~~~---~--~~ To post to this group, send email to sage-devel@googlegroups.com To unsubscribe from this group, send email to sage-devel-unsubscr...@googlegroups.com For more options, visit this group at http://groups.google.com/group/sage-devel URLs: http://www.sagemath.org -~--~~~~--~~--~--~---
[sage-devel] Re: p-adic precision loss in charpoly
On Thu, Mar 19, 2009 at 4:00 PM, David Kohel drko...@gmail.com wrote: Hi, One can always work around a problem, but it would be better to discuss the best approach and put this in place. A related problem I had recently was in finding a random element of GL_n(ZZ/26ZZ) where n was 20-30. It was failing to terminate in the determinant computation. My guess is that a determinant over ZZ was being computed and reduced but that the resulting determinant was too big. I didn't verify this, but invite someone to check. It is trivial to compute the determinant of an nxn matrix over ZZ when n = 30 and the entries of the matrix have 2 digits. That would be the case in your example. Sage wasn't computing your det over ZZ; if it had, then doing that computation would have worked fine. For example: -- | Sage Version 3.4, Release Date: 2009-03-11 | | Type notebook() for the GUI, and license() for information.| -- sage: a = random_matrix(Integers(26), 7) sage: time a.det() CPU times: user 0.03 s, sys: 0.00 s, total: 0.03 s Wall time: 0.05 s 6 sage: a = random_matrix(Integers(26), 8) sage: time a.det() CPU times: user 0.15 s, sys: 0.00 s, total: 0.15 s Wall time: 0.16 s 9 sage: a = random_matrix(Integers(26), 9) sage: time a.det() CPU times: user 1.37 s, sys: 0.01 s, total: 1.37 s Wall time: 1.38 s 23 sage: time Integers(26)(a.lift().det()) CPU times: user 0.02 s, sys: 0.01 s, total: 0.03 s Wall time: 0.39 s 23 sage: time Integers(26)(a.lift().det()) CPU times: user 0.00 s, sys: 0.00 s, total: 0.00 s Wall time: 0.00 s 23 sage: a = random_matrix(Integers(26), 9) sage: time Integers(26)(a.lift().det()) CPU times: user 0.00 s, sys: 0.00 s, total: 0.00 s Wall time: 0.00 s 10 sage: a = random_matrix(Integers(26), 30) sage: time Integers(26)(a.lift().det()) CPU times: user 0.02 s, sys: 0.00 s, total: 0.02 s Wall time: 0.02 s 20 sage: a = random_matrix(Integers(26), 200) sage: time Integers(26)(a.lift().det()) CPU times: user 0.30 s, sys: 0.04 s, total: 0.33 s Wall time: 0.34 s 15 It would thus be far better for now if Sage were to lift to ZZ, do the det, then reduce again. For square-free modulus, a multimodular algorithm would of course be even better. This is now trac #5570: http://trac.sagemath.org/sage_trac/ticket/5570 Instead I worked around this by computing the determinants mod 2 and mod 13 and using CRT (if the determinants were both units). The time was then almost trivial. Suppose I replace this problem over ZZ/25ZZ or ZZ/256ZZ. I would still hope that the problem would NOT be lifted to ZZ for computation, since this would certainly not terminate in reasonable time for a dense matrix. Yes it would. Even for a dense 200x200 matrix over ZZ/256ZZ it only take 0.35 seconds total time to lift *and* compute the determinant, then reduce the result. sage: a = random_matrix(Integers(256), 30) sage: time Integers(256)(a.lift().det()) CPU times: user 0.01 s, sys: 0.00 s, total: 0.01 s Wall time: 0.01 s 222 sage: a = random_matrix(Integers(256), 200) sage: time Integers(256)(a.lift().det()) CPU times: user 0.32 s, sys: 0.04 s, total: 0.35 s Just out of curiosity, why didn't you try any of the above before writing this email? :-) Taking the approach of lifting to ZZ for charpolys of matrices of non-trivial size will undoubtably lead to exactly the same coefficient explosion which is the presumed source of the problem with determinants over ZZ/26ZZ. --David Except it isn't anywhere nearly so bad as you think: sage: a = random_matrix(Integers(256), 500) sage: time Integers(256)(a.lift().det()) CPU times: user 3.70 s, sys: 0.49 s, total: 4.19 s Wall time: 4.04 s 188 sage: a = random_matrix(Integers(2^20), 500) sage: time Integers(256)(a.lift().det()) CPU times: user 7.23 s, sys: 0.81 s, total: 8.04 s Wall time: 7.94 s 208 All timings above on my 2GB of RAM Macbook laptop. William --~--~-~--~~~---~--~~ To post to this group, send email to sage-devel@googlegroups.com To unsubscribe from this group, send email to sage-devel-unsubscr...@googlegroups.com For more options, visit this group at http://groups.google.com/group/sage-devel URLs: http://www.sagemath.org -~--~~~~--~~--~--~---
