Opened 8 years ago
Closed 7 years ago
#10281 closed defect (fixed)
Multimodular echelon form over cyclotomic fields fails
Reported by: | mraum | Owned by: | jason, was |
---|---|---|---|
Priority: | critical | Milestone: | sage-5.0 |
Component: | linear algebra | Keywords: | cyclotomic, echelon |
Cc: | Merged in: | sage-5.0.beta14 | |
Authors: | William Stein | Reviewers: | Martin Raum |
Report Upstream: | N/A | Work issues: | |
Branch: | Commit: | ||
Dependencies: | Stopgaps: |
Description (last modified by )
The following example comes from Eisenstein series. It leads to previous_prime being called with invalid arguments.
sage: K.<rho> = CyclotomicField(106) sage: coeffs = [(18603/107*rho^51 - 11583/107*rho^50 - 19907/107*rho^49 - 13588/107*rho^48 - 8722/107*rho^47 + 2857/107*rho^46 - 19279/107*rho^45 - 16666/107*rho^44 - 11327/107*rho^43 + 3802/107*rho^42 + 18998/107*rho^41 - 10798/107*rho^40 + 16210/107*rho^39 - 13768/107*rho^38 + 15063/107*rho^37 - 14433/107*rho^36 - 19434/107*rho^35 - 12606/107*rho^34 + 3786/107*rho^33 - 17996/107*rho^32 + 12341/107*rho^31 - 15656/107*rho^30 - 19092/107*rho^29 + 8382/107*rho^28 - 18147/107*rho^27 + 14024/107*rho^26 + 18751/107*rho^25 - 8301/107*rho^24 - 20112/107*rho^23 - 14483/107*rho^22 + 4715/107*rho^21 + 20065/107*rho^20 + 15293/107*rho^19 + 10072/107*rho^18 + 4775/107*rho^17 - 953/107*rho^16 - 19782/107*rho^15 - 16020/107*rho^14 + 5633/107*rho^13 - 17618/107*rho^12 - 18187/107*rho^11 + 7492/107*rho^10 + 19165/107*rho^9 - 9988/107*rho^8 - 20042/107*rho^7 + 10109/107*rho^6 - 17677/107*rho^5 - 17723/107*rho^4 - 12489/107*rho^3 - 6321/107*rho^2 - 4082/107*rho - 1378/107, 1, 4*rho + 1), (0, 1, rho + 4)] sage: m = matrix(2, coeffs) sage: v = m.echelon_form() --------------------------------------------------------------------------- ValueError Traceback (most recent call last) /home/martin/<ipython console> in <module>() /opt/sage/local/lib/python2.6/site-packages/sage/matrix/matrix_cyclo_dense.so in sage.matrix.matrix_cyclo_dense.Matrix_cyclo_dense.echelon_form (sage/matrix/matrix_cyclo_dense.cpp:13469)() /opt/sage/local/lib/python2.6/site-packages/sage/matrix/matrix_cyclo_dense.so in sage.matrix.matrix_cyclo_dense.Matrix_cyclo_dense._echelon_form_multimodular (sage/matrix/matrix_cyclo_dense.cpp:14366)() /opt/sage/local/lib/python2.6/site-packages/sage/rings/arith.pyc in previous_prime(n) 1044 n = ZZ(n)-1 1045 if n <= 1: -> 1046 raise ValueError, "no previous prime" 1047 if n <= 3: 1048 return ZZ(n) ValueError: no previous prime
Apply:
Install:
Attachments (1)
Change History (17)
comment:1 Changed 7 years ago by
comment:2 Changed 7 years ago by
This turned out to be easier to resolve than it might have otherwise been, since Burcin Erocal and Martin Albrecht had improved matrices modn so they can now handle a modulus up to 2^{23. However, they hadn't adapted the chinese remainder theorem infrastructure to take this into account. }
I've written code that makes it so for dense linear algebra we have all the primes up to 2^{23 at our disposal, instead of just the primes up to 46341. There are more, bigger, and better primes up to 2}23:
sage: n = prod(prime_range(46431)) sage: len(str(n)) 20025 sage: time n = prod(prime_range(2^23)) Time: CPU 0.90 s, Wall: 0.90 s sage: len(str(n)) 3641670 sage: prime_pi(2^23) 564163 sage: prime_pi(46431) 4797 sage: prime_pi(2^23)//prime_pi(46431) 117
The upshot is that in various contexts we can deal with answers with up to 3.5 million digits, instead of only 20,025 digits, which was pretty feeble. Moreover, since the linear algebra is just as fast or faster, and we start with the biggest prime we can, certain linear algebra computations should be much faster with this new code.
I did *not* fix anything for sparse matrices, except avoiding introducing a major bug (caught by doctests). See #12679 for sparse matrices, as this should be another ticket.
Also, obviously, this fix is only a partial fix, and we'll someday have to deal with running out of primes on huge problems. But it'll take around 150 times as long to ever reach such problems, for what it is worth.
Note that it would probably be trivial to change the *old* matrix_modn_dense code to use 64-bit ints, hence support a modulus up to 2^{31. However the new Erocal-Albrecht code can't go beyond 2}23. So that code -- which they deprecated -- should probably be resurrected. That's a project for later.
comment:3 Changed 7 years ago by
- Status changed from new to needs_review
comment:4 Changed 7 years ago by
Note - the patch is bigger than expected since I fixed some of the formatting of docstrings in a file, just because.
comment:5 follow-up: ↓ 7 Changed 7 years ago by
Apply trac_10281-sage-5.0-beta9.patch
(for the patchbot)
I want the latest version to be tested on a variety of machines.
In principle, the patch looks good. It would be create if someone had a 32-bit system to test this, though.
comment:6 Changed 7 years ago by
Looks good to me. The only comment I have is that the changes adding
cdef extern from "stdint.h": ctypedef unsigned long uint_fast64_t cdef extern from "multi_modular.h": ctypedef uint_fast64_t mod_int
are unnecessary, declaring it as an unsigned long in Cython is sufficient.
Changed 7 years ago by
This patch is for beta9. I've made the one change that robertwb suggests to this patch (the formal cython change)
comment:7 in reply to: ↑ 5 Changed 7 years ago by
Replying to mraum:
Apply trac_10281-sage-5.0-beta9.patch
(for the patchbot)
I want the latest version to be tested on a variety of machines.
In principle, the patch looks good. It would be create if someone had a 32-bit system to test this, though.
I've started a 32-bit Ubuntu 10.4 virtual machine on boxen, and I'm running tests on it. I'll report back.
comment:8 Changed 7 years ago by
I tested on 32-bit and it works, but *only* after making a 1-line change to the Linbox spkg (which isn't needed on 64-bit, evidently).
I've posted the new linbox spkg here:
So to apply this ticket:
- Build new linbox: http://wstein.org/home/wstein/patches/linbox-1.1.6.p8.spkg
- Apply the patch to the sage repo: trac_10281-sage-5.0-beta9.patch
William
comment:9 Changed 7 years ago by
- Priority changed from major to critical
comment:10 follow-up: ↓ 11 Changed 7 years ago by
- Description modified (diff)
I have tested the spkg, and it is completely ok. The commit message could be a bit more informative, but since it refers to the ticket, it's ok also.
William, do you want to make the change suggested by Robert? Personally, I would keep the extern statements for clearness (they obviously don't do any harm).
If you don't intend to, then I would set this to positive review, so that everything can be run on the build farm.
comment:11 in reply to: ↑ 10 Changed 7 years ago by
Replying to mraum:
I have tested the spkg, and it is completely ok. The commit message could be a bit more informative, but since it refers to the ticket, it's ok also.
William, do you want to make the change suggested by Robert? Personally, I would keep the extern statements for clearness (they obviously don't do any harm).
I did make the changes and updated the second patch. I didn't bother with the beta2 patch, since that one will never get applied by the release manager.
If you don't intend to, then I would set this to positive review, so that everything can be run on the build farm.
Please do. Thanks!
I have to point out that I was running some big modular symbols computations over cyclotomic fields and found some "not enough primes" issues even with this patch. That's because we only get to use split primes. So it will be very important to revive the old matrix_modn_dense code, but using a 64-bit int instead of 32-bit, so we get modulo up to 2^32
(instead of 2^23
). (That is already #4968.) I plan to do that this week or next week. The current patch here is very good foundation for doing #4968 right.
comment:12 Changed 7 years ago by
- Status changed from needs_review to positive_review
Great! Let me know when you have made progress on #4968. I hope to have time to review it.
comment:13 Changed 7 years ago by
- Reviewers set to Martin Raum
comment:14 Changed 7 years ago by
- Description modified (diff)
comment:15 Changed 7 years ago by
- Milestone set to sage-5.0
comment:16 Changed 7 years ago by
- Merged in set to sage-5.0.beta14
- Resolution set to fixed
- Status changed from positive_review to closed
This is not just an issue with cyclotomic linear algebra. I'm doing some computations of modular symbols spaces of weight 4 with trivial character, and just hit "no previous prime". We need real 64-bit mod-p linear algebra to fix this issue. We need more primes!