Opened 12 months ago
Last modified 101 minutes ago
#33360 needs_review defect
is_prime for ideals uses factorization, can be VERY slow
| Reported by: | Gonzalo Tornaría | Owned by: | |
|---|---|---|---|
| Priority: | critical | Milestone: | |
| Component: | number fields | Keywords: | |
| Cc: | Lorenz Panny, Frédéric Chapoton | Merged in: | |
| Authors: | Gonzalo Tornaría, Lorenz Panny | Reviewers: | Lorenz Panny |
| Report Upstream: | N/A | Work issues: | |
| Branch: | public/33360 (Commits, GitHub, GitLab) | Commit: | b31f59bbae2d851431154103adf2e995ddda76d0 |
| Dependencies: | Stopgaps: |
Status badges
Description (last modified by )
The implementation of is_prime() for ideals uses pari idealfactor(). This is a bad idea if the norm of the ideal is hard to factor.
Pari has ismaximalideal(), but it is broken in pari 2.13.3. As a workaround, since a prime ideal has prime power norm, we avoid factorization except on this (easy) case.
Was: VERY slow doctest in sage/rings/number_field/number_field.py
With sage 9.5:
$ time sage -t --long --warn-long 60.0 --random-seed=40834229707572602474454926153679777311 src/sage/rings/number_field/number_field.py
Running doctests with ID 2022-02-16-15-58-06-58d0ee96.
Using --optional=build,pip,sage,sage_spkg,void
Features to be detected: 4ti2,benzene,bliss,buckygen,conway_polynomials,csdp,database_cremona_ellcurve,database_cremona_mini_ellcurve,database_jones_numfield,database_knotinfo,dvipng,graphviz,imagemagick,jupymake,kenzo,latte_int,lrslib,mcqd,meataxe,pandoc,pdf2svg,plantri,pynormaliz,rubiks,sage.combinat,sage.geometry.polyhedron,sage.graphs,sage.plot,sage.rings.number_field,sage.rings.real_double,sage.symbolic,sage_numerical_backends_coin,sagemath_doc_html,sphinx,tdlib
Doctesting 1 file.
sage -t --long --warn-long 60.0 --random-seed=40834229707572602474454926153679777311 src/sage/rings/number_field/number_field.py
**********************************************************************
File "src/sage/rings/number_field/number_field.py", line 2180, in sage.rings.number_field.number_field.NumberField_generic.?
Warning, slow doctest:
while not K.random_element().is_prime(): # long time
pass
Test ran for 1319.62 s
[2265 tests, 1367.95 s]
----------------------------------------------------------------------
All tests passed!
----------------------------------------------------------------------
Total time for all tests: 1368.3 seconds
cpu time: 1364.7 seconds
cumulative wall time: 1367.9 seconds
Features detected for doctesting:
real 22m51.641s
user 22m47.153s
sys 0m1.272s
Change History (24)
comment:1 Changed 12 months ago by
comment:2 Changed 11 months ago by
Ctrl-C works for me (I instantaneously get KeyboardInterrupt) with 9.6b1 on MacOS 11.5.2.
comment:3 Changed 11 months ago by
For the record:
- I have a working patch for
is_prime()that I'll polish and submit later - About the issue with Ctrl-C, it turns out this is caused by a recent upgrade of prompt_toolkit (a system package) from 3.0.24 to 3.0.26. I didn't try with 3.0.25 (the bundled version in sage is 3.0.22 and that one works ok)
comment:4 Changed 11 months ago by
| Cc: | Lorenz Panny added |
|---|
comment:5 Changed 11 months ago by
| Cc: | Lorenz Panny removed |
|---|
comment:6 Changed 11 months ago by
| Cc: | Lorenz Panny added |
|---|
comment:7 Changed 11 months ago by
| Authors: | → Gonzalo Tornaría |
|---|---|
| Branch: | → u/tornaria/33360-ideal.is_prime |
| Commit: | → 6330cee1a69229b2f5d8ac4a09a73d8152f5e14d |
| Status: | new → needs_review |
To fix the issue, we only factor the ideal if the norm is a prime power.
Note that a simpler way would be to use pari ismaximalideal() as in:
K = self.number_field().pari_nf()
I = self.pari_hnf()
self._pari_prime = K.idealismaximal(I) or None
return self._pari_prime is not None
However it turns out that ismaximalideal() is broken in pari 2.13.3. For a quick example in gp:
? K = nfinit(y^2 + 1); P = idealhnf(K, 3); idealismaximal(K, P) %1 = 0
This is fixed in current development version of pari by this commit.
New commits:
| 6330cee | trac 33360: avoid factoring in is_prime()
|
comment:8 Changed 10 months ago by
| Component: | doctest framework → number fields |
|---|---|
| Description: | modified (diff) |
| Priority: | major → critical |
| Summary: | VERY slow doctest in sage/rings/number_field/number_field.py → is_prime for ideals uses factorization, can be VERY slow |
comment:9 follow-up: 10 Changed 10 months ago by
| Reviewers: | → Lorenz Panny |
|---|
The patch looks good as far as I can tell.
Just one question: How much faster is idealismaximal() when using a non-broken version of PARI? Could it be worthwhile to detect the PARI version at runtime to switch between idealismaximal() and the workaround?
comment:10 Changed 10 months ago by
Replying to lorenz:
The patch looks good as far as I can tell.
Just one question: How much faster is
idealismaximal()when using a non-broken version of PARI? Could it be worthwhile to detect the PARI version at runtime to switch betweenidealismaximal()and the workaround?
I doubt the difference is significative. For the example in the doctest:
sage: K.<a,b,c> = NumberField([x^2-2,x^2-3,x^2-5]) sage: t = (((-2611940*c + 1925290/7653)*b - 1537130/7653*c ....: + 10130950)*a + (1343014/7653*c - 8349770)*b ....: + 6477058*c - 2801449990/4002519) sage: I = K.ideal(t).pari_hnf() sage: sage: sage: K.<a,b,c> = NumberField([x^2-2,x^2-3,x^2-5]) sage: t = (((-2611940*c + 1925290/7653)*b - 1537130/7653*c ....: + 10130950)*a + (1343014/7653*c - 8349770)*b ....: + 6477058*c - 2801449990/4002519) sage: I = K.ideal(t) sage: Kp = K.pari_nf() sage: Ip = I.pari_hnf() sage: %timeit Kp.idealismaximal(Ip) 14.4 µs ± 125 ns per loop (mean ± std. dev. of 7 runs, 100,000 loops each) sage: %timeit Kp.idealnorm(Ip) 2.26 µs ± 11.1 ns per loop (mean ± std. dev. of 7 runs, 100,000 loops each)
which probably just means pari is optimizing in a different way but that will depend on the ideal. In the long run, I think using idealismaximal() is better since pari is in a better position to optimize (and any effort to optimize is better spent in pari), but not worth complicating things.
Using idealfactor() on that ideal would take 20-30 minutes so a microsecond here and there is just insignificant in contrast.
For a prime ideal, I guess the time would be dominated by isprime() anyway...
comment:11 Changed 10 months ago by
The idealnorm() is indeed negligible, but idealfactor() is much slower than idealismaximal() even in the prime-power case.
Example:
R.<x> = QQ[] K.<a,b,c> = NumberField([x^2-2,x^2-3,x^2-5]) Kpari = K.pari_nf() for bits in (50, 100, 200, 300, 500, 800): p = random_prime(2^bits, proof=False) I = choice([f for f,_ in K.ideal(p).factor()]) # a prime above p Ipari = I.pari_hnf() print(f'{bits} bits') %timeit Kpari.idealismaximal(Ipari) %timeit assert Kpari.idealnorm(Ipari).isprimepower(); Kpari.idealfactor(Ipari)
Result:
50 bits 790 µs ± 2.97 µs per loop (mean ± std. dev. of 7 runs, 1,000 loops each) 799 µs ± 5 µs per loop (mean ± std. dev. of 7 runs, 1,000 loops each) 100 bits 3.35 ms ± 69.5 µs per loop (mean ± std. dev. of 7 runs, 100 loops each) 17.4 ms ± 617 µs per loop (mean ± std. dev. of 7 runs, 100 loops each) 200 bits 7.52 ms ± 394 µs per loop (mean ± std. dev. of 7 runs, 100 loops each) 52.9 ms ± 733 µs per loop (mean ± std. dev. of 7 runs, 10 loops each) 300 bits 11.9 ms ± 389 µs per loop (mean ± std. dev. of 7 runs, 100 loops each) 178 ms ± 6.06 ms per loop (mean ± std. dev. of 7 runs, 10 loops each) 500 bits 64 ms ± 1.3 ms per loop (mean ± std. dev. of 7 runs, 10 loops each) 973 ms ± 11 ms per loop (mean ± std. dev. of 7 runs, 1 loop each) 800 bits 61.4 ms ± 5 ms per loop (mean ± std. dev. of 7 runs, 10 loops each) 4.45 s ± 131 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)
comment:12 Changed 10 months ago by
We have to be careful with idealismaximal() though because (according to the PARI documentation) it uses a probabilistic primality test. So if Sage's flag for strict primality proving is set, we should still run isprimepower() prior to calling idealismaximal(), or ensure primality of the relevant entry of the structure returned by idealismaximal() afterwards.
Also, a small optimization: Calling self.norm() instead of idealnorm() will cache the norm instead of computing it from scratch every time.
comment:13 Changed 9 months ago by
| Milestone: | sage-9.6 → sage-9.7 |
|---|
comment:14 Changed 4 months ago by
| Milestone: | sage-9.7 → sage-9.8 |
|---|
comment:15 Changed 4 months ago by
| Authors: | Gonzalo Tornaría → Gonzalo Tornaría, Lorenz Panny |
|---|---|
| Branch: | u/tornaria/33360-ideal.is_prime → public/33360 |
| Commit: | 6330cee1a69229b2f5d8ac4a09a73d8152f5e14d → b938f7b81edb954aa3639043f0e09f15b2b98d1a |
I implemented my suggestions from comment:9 and comment:12. Someone else will have to review the new changes.
New commits:
| 6330cee | trac 33360: avoid factoring in is_prime()
|
| 7f19cea | Merge branch 'u/tornaria/33360-ideal.is_prime' into public/33360
|
| b938f7b | opportunistically use idealismaximal() when PARI version is recent enough
|
comment:16 Changed 4 months ago by
| Commit: | b938f7b81edb954aa3639043f0e09f15b2b98d1a → e3af4e3a1e11be18a340946dd97909bd911af1e5 |
|---|
comment:17 Changed 4 months ago by
| Dependencies: | → #34537 |
|---|
comment:18 Changed 4 months ago by
| Commit: | e3af4e3a1e11be18a340946dd97909bd911af1e5 → 58511ea3015159e1c1533cda12281aa9289e36c0 |
|---|
Branch pushed to git repo; I updated commit sha1. New commits:
| 58511ea | remove unused import
|
comment:19 Changed 3 months ago by
| Commit: | 58511ea3015159e1c1533cda12281aa9289e36c0 → 7eee1c4317586bc31dcded5e82837148c4bf8d54 |
|---|
comment:20 Changed 3 months ago by
| Cc: | Frédéric Chapoton added |
|---|---|
| Dependencies: | #34537 |
Reverted earlier changes from comment:16 as #34537 appears to be stuck. Can we get this in?
comment:21 Changed 3 months ago by
| Commit: | 7eee1c4317586bc31dcded5e82837148c4bf8d54 → 689eb880afd39b6bb8916be32c1b92df5c425bf3 |
|---|
Branch pushed to git repo; I updated commit sha1. New commits:
| 689eb88 | fix docstring markup
|
comment:22 Changed 3 months ago by
| Commit: | 689eb880afd39b6bb8916be32c1b92df5c425bf3 → 4ac2c9ba19ccdb80b227d02920d128fe360101e5 |
|---|
Branch pushed to git repo; I updated commit sha1. New commits:
| 4ac2c9b | fix docstring markup some more
|
comment:23 Changed 7 hours ago by
| Milestone: | sage-9.8 |
|---|
comment:24 Changed 101 minutes ago by
| Commit: | 4ac2c9ba19ccdb80b227d02920d128fe360101e5 → b31f59bbae2d851431154103adf2e995ddda76d0 |
|---|
Tracked this down to
Behind the scenes this is handled by pari as follows:
It is possible to use pari to test for primality much faster without factoring as follows:
I'll try to cook up a patch for
is_prime()to useidealismaximal()instead ofidealfactor().Another issue: during the call to
idealfactor(), the one that takes 27 minutes, hitting control C doesn't work for me, but I'm using all system packages (in particular cypari, cysignals, python 3.10, etc) can someone confirm that Ctrl-C is broken on a standard build of sage?