Re: faster for very big primes and small exponents than "isprimepower()"

Karim Belabas <[email protected]>
Newsgroups gmane.comp.mathematics.pari.user
Message-ID <[email protected]>
Hi Hermann,

  for a fair comparison between your procedure and isprimepower, then
you must include a primality proof (isprime) in the former.

If the goal is to try and produce a base n given N = n^k without checking
whether n is prime, then benchmark your procedure against ispower.

Cheers,

    K.B.
-- 
Pr. Karim Belabas, U. Bordeaux, Vice-président en charge du Numérique
Institut de Mathématiques de Bordeaux UMR 5251 - (+33) 05 40 00 29 77
http://www.math.u-bordeaux.fr/~kbelabas/


* [email protected] [2025-05-18 01:30]:
> On mersenne forum I learned about a method to determine prime of prime
> power:
> https://www.mersenneforum.org/node/1055690?p=1077586#post1077586
> 
> After I simplified modulus and exponent, it can be proven with Euler's
> theorem and induction that for odd prime powers N=p^n
> 2^N ≡ 2 (mod N)
> 
> This determines the prime for most odd prime powers:
> lift(gcd(Mod(2,N)^N-2,N))
> 
> The exceptions are the two Wiefrich primes below 18*10^18
> (1093 and 3511), for which the square of the prime is returned.
> 
> 
> PARI/GP "isprimepower(N,&r);return(r)" is nearly always much faster.
> I created script to demonstrate that the gcd is faster for very big primes
> and low exponents:
> https://www.mersenneforum.org/node/1055690?p=1077586#post1077586
> 
> hermann@7950x:~$ gp -q < <(curl --no-progress-meter https://gist.githubusercontent.com/Hermann-SW/e4c2ccccd918e59766dc82b01993d41b/raw/e7480819b8a57331d84cd30f6a9f89650a44388e/measure.gp)
> runtimes in milliseconds, with N=nextprime(base)^expo
> 
> isprimepower(N,&r);return(r)
>       |  2|  3|  4|  5|  6|  7|  8|  9| 10| 11| 12| 13| 14| 15|
>       |---|---|---|---|---|---|---|---|---|---|---|---|---|---|
> 10^ 10|  0|  0|  0|  0|  0|  0|  0|  0|  0|  0|  0|  0|  0|  0|
> 10^ 20|  0|  0|  0|  0|  0|  0|  0|  1|  0|  0|  0|  0|  0|  0|
> 10^ 30|  1|  1|  1|  1|  1|  1|  1|  1|  1|  1|  0|  0|  1|  1|
> 10^ 40|  2|  2|  2|  2|  2|  2|  1|  2|  2|  2|  2|  2|  1|  1|
> 10^ 50|  1|  0|  0|  1|  0|  0|  1|  0|  0|  1|  0|  0|  0|  0|
> 10^ 60|  5|  5|  5|  5|  5|  5|  5|  4|  4|  5|  5|  5|  5|  5|
> 10^ 70| 23| 23| 23| 23| 23| 23| 23| 23| 23| 23| 23| 23| 23| 23|
> 10^ 80| 13| 13| 13| 13| 12| 12| 12| 12| 12| 12| 12| 12| 12| 12|
> 10^ 90| 18| 18| 19| 19| 19| 19| 18| 18| 18| 17| 16| 16| 17| 16|
> 10^100| 23| 23| 24| 24| 24| 24| 23| 23| 23| 23| 23| 23| 23| 24|
> 10^110|  0|  1|  0|  1|  0|  1|  0|  1|  0|  1|  0|  1|  0|  1|
> 10^120| 51| 50| 51| 51| 51| 51| 51| 51| 51| 51| 51| 51| 51| 51|
> 10^130|109|109|108|108|109|109|108|109|108|109|108|108|109|109|
> 10^140|115|115|115|115|115|115|115|115|114|115|115|115|115|115|
> 10^150|136|136|136|136|136|136|136|136|136|136|136|136|136|136|
> 10^160|217|216|216|216|216|216|216|216|216|216|216|217|216|216|
> 
> lift(gcd(Mod(2,N)^N-2,N))
>       |  2|  3|  4|  5|  6|  7|  8|  9| 10| 11| 12| 13| 14| 15|
>       |---|---|---|---|---|---|---|---|---|---|---|---|---|---|
> 10^ 10|  0|  0|  0|  0|  0|  0|  0|  0|  0|  0|  0|  0|  0|  0|
> 10^ 20|  0|  0|  0|  0|  0|  0|  0|  0|  0|  0|  0|  0|  1|  0|
> 10^ 30|  0|  0|  0|  0|  0|  0|  1|  0|  0|  1|  0|  1|  0|  1|
> 10^ 40|  0|  0|  0|  0|  0|  0|  1|  0|  1|  1|  1|  1|  2|  1|
> 10^ 50|  0|  0|  0|  0|  1|  0|  1|  1|  1|  2|  1|  2|  3|  3|
> 10^ 60|  0|  0|  0|  1|  0|  1|  1|  1|  2|  3|  3|  3|  5|  5|
> 10^ 70|  0|  0|  0|  1|  1|  1|  1|  2|  2|  4|  4|  6|  6|  8|
> 10^ 80|  0|  0|  0|  0|  1|  1|  3|  3|  4|  5|  6|  7|  9| 11|
> 10^ 90|  0|  0|  0|  1|  2|  2|  3|  4|  5|  7|  8| 10| 12| 15|
> 10^100|  0|  0|  1|  1|  1|  3|  4|  5|  6|  9| 11| 13| 16| 20|
> 10^110|  0|  0|  0|  2|  2|  4|  5|  7|  8| 10| 14| 18| 21| 26|
> 10^120|  0|  0|  1|  2|  3|  4|  6|  8| 11| 14| 18| 22| 28| 32|
> 10^130|  0|  1|  1|  2|  4|  5|  8| 10| 13| 17| 22| 28| 33| 41|
> 10^140|  0|  1|  2|  3|  4|  7|  9| 13| 16| 21| 27| 33| 40| 49|
> 10^150|  1|  1|  2|  3|  5|  8| 11| 15| 20| 26| 32| 39| 48| 55|
> 10^160|  0|  1|  2|  4|  6|  9| 13| 17| 24| 31| 38| 47| 54| 58|
> hermann@7950x:~$
> 
> 
> Regards,
> 
> Hermann.
lmpx.com only provides a reader for public news (NNTP) servers. It is not affiliated with the servers or forums shown here and is not responsible for the content of articles, which is written by their respective authors.