#16880 closed enhancement (fixed)
previous_prime, previous_prime_power, next_prime_power
Reported by:  vdelecroix  Owned by:  

Priority:  major  Milestone:  sage6.7 
Component:  number theory  Keywords:  
Cc:  ncohen, jdemeyer  Merged in:  
Authors:  Vincent Delecroix  Reviewers:  Jeroen Demeyer 
Report Upstream:  N/A  Work issues:  
Branch:  35c0211 (Commits)  Commit:  
Dependencies:  Stopgaps: 
Description (last modified by )
In this ticket we add three methods to integers: previous_prime
(calling pari precprime
), next_prime_power
and previous_prime_power
. We also make changes in sage.rings.arith
to call these (faster) methods.
Timings on the current version (sage6.7.beta2)
sage: timeit("for i in range(200): next_prime_power(i)", number=5000, repeat=10) 5000 loops, best of 10: 311 µs per loop sage: a = 33622489 sage: timeit("next_prime_power(a)", number=5000, repeat=10) 5000 loops, best of 10: 53.6 µs per loop
Timings with the branch applied
sage: timeit("for i in range(200): next_prime_power(i)", number=5000, repeat=10) 5000 loops, best of 10: 184 µs per loop sage: a = 33622489 sage: timeit("next_prime_power(a)", number=5000, repeat=10) 5000 loops, best of 10: 16.9 µs per loop
Change History (30)
comment:1 Changed 6 years ago by
comment:2 Changed 6 years ago by
 Branch set to u/vdelecroix/16880
 Description modified (diff)
 Status changed from new to needs_review
comment:3 Changed 6 years ago by
 Commit set to 6c6cbe25b41dece8513ad3527c68415517328710
Branch pushed to git repo; I updated commit sha1. Last 10 new commits:
1d103da  Trac 15767: fix doctest in Ser()

62ca821  Trac 15767: explain application of Sturm's theorem

bf435f8  Merge tag '6.4.beta1' into ticket/15767

6e9f686  trac #16878: faster is_prime

d7dc389  trac #16878: review

175d6b5  trac #16878: review 2

d55dc67  trac #16878: bikeshedding

f73bd87  trac #16869: fast prime_powers

8259a7d  trac #16869: fix prime_powers(0,1)

6c6cbe2  trac #16880: implement prime_power_range

comment:4 Changed 6 years ago by
 Commit changed from 6c6cbe25b41dece8513ad3527c68415517328710 to e7b94588ee7ecb82c00afa48f905622f5c374cdd
Branch pushed to git repo; I updated commit sha1. This was a forced push. New commits:
e7b9458  trac #16880: prime_power_range

comment:5 Changed 6 years ago by
 Dependencies #16869, #16878 deleted
 Description modified (diff)
 Summary changed from even faster prime_powers to prime_power_range
Rebased on 6.4.beta2, no dependency anymore...
Vincent
comment:6 Changed 6 years ago by
 Commit changed from e7b94588ee7ecb82c00afa48f905622f5c374cdd to a3195f7bb4e16cc4eb3f3ae4ef58a4b5d0341ad2
Branch pushed to git repo; I updated commit sha1. This was a forced push. New commits:
6e9f686  trac #16878: faster is_prime

d7dc389  trac #16878: review

175d6b5  trac #16878: review 2

d55dc67  trac #16878: bikeshedding

7f61321  Reorganize is_prime_power() and next_prime() logic

5937af5  Sloane considers 1 to be a prime power

f0d716f  trac #16878: certificate flag in is_prime_power

bfd69d3  is_prime_power() reviewer fixes

9e04c7a  trac #16880: merge #16878 (faster is_prime)

a3195f7  trac #16880: prime_power_range

comment:7 Changed 6 years ago by
 Dependencies set to #16878
 Description modified (diff)
And with timings
using prime_powers
(EDIT: only for commit at a3195f7 and not at e3be1ce)
sage: timeit("prime_powers(1000)", number=4000) 4000 loops, best of 3: 86.7 µs per loop sage: timeit("prime_powers(5000)", number=2000) 2000 loops, best of 3: 314 µs per loop
using prime_power_range
instead
sage: timeit("prime_power_range(1000)", number=4000) 4000 loops, best of 3: 17.6 µs per loop sage: timeit("prime_power_range(5000)", number=2000) 2000 loops, best of 3: 72.3 µs per loop
using prime_power_range
with py_ints=True
sage: timeit("prime_power_range(1000,py_ints=True)", number=4000) 4000 loops, best of 3: 14.6 µs per loop sage: timeit("prime_power_range(5000,py_ints=True)", number=2000) 2000 loops, best of 3: 33.9 µs per loop
comment:8 Changed 6 years ago by
 Description modified (diff)
 Summary changed from prime_power_range to prime_powers and prime_power_range
comment:9 Changed 6 years ago by
 Commit changed from a3195f7bb4e16cc4eb3f3ae4ef58a4b5d0341ad2 to e3be1ce5b13d007baf6982e1b94a2299c68744c1
Branch pushed to git repo; I updated commit sha1. This was a forced push. New commits:
e3be1ce  trac #16880: prime_powers and prime_power_range

comment:10 Changed 6 years ago by
It looks more natural to me to follow what primes
and prime_range
do... this is what I did in e3be1ce.
Vincent
comment:11 Changed 6 years ago by
 Milestone changed from sage6.4 to sage6.6
 Status changed from needs_review to needs_work
Needs rebase.
comment:12 Changed 6 years ago by
 Dependencies changed from #16878 to #16878, #17852
comment:13 Changed 6 years ago by
 Commit changed from e3be1ce5b13d007baf6982e1b94a2299c68744c1 to 25176c520e538f0a98334c7ffdfe69c14207b32d
Branch pushed to git repo; I updated commit sha1. This was a forced push. Last 10 new commits:
b365b9e  trac #17852: fix binomial(n,0) for negative n

faed0c0  trac #17852: remove lcm on rationals (category takes care)

cbbb06f  trac #17852: make a reference to the symbolic binomial

1655367  trac #17852: replace powermod_ui by hash in the doc

d401269  trac #17852: fix output type of arith.binomial

8fc6b3d  trac #17852: type > parent

61b94c2  trac #17852: fix curious failing doctest in multi_polynomial_ideal

6697e68  trac #17852: be nicer on nonintegral domains

2260680  trac #17852: merge beta6.6

25176c5  trac #16880: cleanup for prime and prime powers

comment:14 Changed 6 years ago by
 Description modified (diff)
 Status changed from needs_work to needs_review
 Summary changed from prime_powers and prime_power_range to further cleanup in arith.py/integer.pyx
comment:15 Changed 6 years ago by
 Commit changed from 25176c520e538f0a98334c7ffdfe69c14207b32d to b3cfd3cce04b352ee86b7b84c3880f0c5d5829e9
Branch pushed to git repo; I updated commit sha1. New commits:
b3cfd3c  trac 16880: fix doctest related to random_primes

comment:16 Changed 6 years ago by
 Dependencies #16878, #17852 deleted
 Description modified (diff)
 Milestone changed from sage6.6 to sage6.7
 Summary changed from further cleanup in arith.py/integer.pyx to previous_prime, previous_prime_power, next_prime_power
comment:17 Changed 6 years ago by
 Commit changed from b3cfd3cce04b352ee86b7b84c3880f0c5d5829e9 to db1e869640dbb20d39e565970aaa8c5d58fa8953
comment:18 Changed 6 years ago by
 Description modified (diff)
comment:19 Changed 6 years ago by
 Description modified (diff)
comment:20 followup: ↓ 23 Changed 6 years ago by
 I don't know if
lesser or equal than
is proper English. How aboutless than or equal to
?
PARI
instead ofPari
.
 In
next_prime_power
andprevious_prime_power
, a small optimization is possible: you can replacempz_set(n.value, self.value) mpz_add_ui(n.value, n.value, 1)
by
mpz_add_ui(n.value, self.value, 1)
and in next_prime_power
, replace
if mpz_cmp_ui(n.value, 2) < 0: mpz_set_ui(n.value, 2) return n
by
if mpz_cmp_ui(self.value, 2) < 0: return smallInteger(2)
 Another algorithmic comment for
next_prime_power
(an analogous comment holds forprevious_prime_power
): instead of computing the next power of 2 and doingwhile mpz_cmp(m.value, n.value) == 1:
I would instead check if mpz_sizeinbase(n.value, 2)
increases during the loop. In that case, we know that we skipped a power of 2.
comment:21 Changed 6 years ago by
 Status changed from needs_review to needs_work
comment:22 Changed 6 years ago by
 Commit changed from db1e869640dbb20d39e565970aaa8c5d58fa8953 to 1e0956157ef8f49bc19eb145dc60589f8e59587a
Branch pushed to git repo; I updated commit sha1. New commits:
1e09561  Trac 16880: optimization + Pari > PARI

comment:23 in reply to: ↑ 20 Changed 6 years ago by
 Status changed from needs_work to needs_review
Replying to jdemeyer:
 Another algorithmic comment for
next_prime_power
(an analogous comment holds forprevious_prime_power
): instead of computing the next power of 2 and doingwhile mpz_cmp(m.value, n.value) == 1:I would instead check if
mpz_sizeinbase(n.value, 2)
increases during the loop. In that case, we know that we skipped a power of 2.
I used mpz_tstbit
instead. Should be the fastest.
Potential micro optimization that would need to use mpn_*
functions
 the call to
mpz_sizeinbase
can be replaced bympn_count_leading_zeros
and some arithmetic  we can avoid a useless bit shift in
mpz_tstbit
by only considering the&
with the most significant limb. Moreover, the allocation of the limbs will not change inside the loop. So we can already use the good pointer to the concerned limb.
I think this is too much for a naive algorithm...
comment:24 Changed 6 years ago by
 Reviewers set to Jeroen Demeyer
 Status changed from needs_review to needs_work
skept
> skipped
This should be moved after the comparison with 2 (and self
should be compared with 2):
cdef Integer n = PY_NEW(Integer)
Another simplification that I missed: you can skip the block
if mpz_even_p(n.value): if n.is_prime_power(proof=proof): return n mpz_add_ui(n.value, n.value, 1)
if you instead just add 2 if self
was odd:
mpz_add_ui(n.value, self.value, 1 if mpz_even_p(self.value) else 2)
Then the bit_index
should be based on the value of self
, such that the check works correctly if self
was 2^n  1
.
comment:25 Changed 6 years ago by
 Commit changed from 1e0956157ef8f49bc19eb145dc60589f8e59587a to 35c0211766537e6769cd060ed528ad9d77e60068
Branch pushed to git repo; I updated commit sha1. New commits:
35c0211  Trac 16880: optimization

comment:26 Changed 6 years ago by
 Status changed from needs_work to needs_review
comment:27 Changed 6 years ago by
 Status changed from needs_review to positive_review
comment:28 Changed 6 years ago by
 Branch changed from u/vdelecroix/16880 to 35c0211766537e6769cd060ed528ad9d77e60068
 Resolution set to fixed
 Status changed from positive_review to closed
comment:29 followup: ↓ 30 Changed 6 years ago by
 Commit 35c0211766537e6769cd060ed528ad9d77e60068 deleted
Some copy/paste accidents (could probably be corrected on the followups):
+ def previous_prime(self, proof=None): + r""" + Returns the previous prime before self. + + This method calls the PARI ``precprime`` function. + + INPUT: + +  ``proof``  if ``True`` ensure that the returned value is the next + prime power and if set to ``False`` uses probabilistic methods + (i.e. the result is not guaranteed). By default it uses global + configuration variables to determine which alternative to use (see + :mod:`proof.arithmetic` or :mod:`sage.structure.proof`).
+ def previous_prime_power(self, proof=None): + r""" + Return the previous prime power before self. + + INPUT: + +  ``proof``  if ``True`` ensure that the returned value is the next + prime power and if set to ``False`` uses probabilistic methods + (i.e. the result is not guaranteed). By default it uses global + configuration variables to determine which alternative to use (see + :mod:`proof.arithmetic` or :mod:`sage.structure.proof`).
comment:30 in reply to: ↑ 29 Changed 6 years ago by
Replying to leif:
Some copy/paste accidents (could probably be corrected on the followups):
Ha right.
+1 for "primes"