Ilya Zakharevich on Fri, 12 Jan 2024 04:09:58 +0100


[Date Prev] [Date Next] [Thread Prev] [Thread Next] [Date Index] [Thread Index]

Implementation of forprime() and unextprime()


I looked through the code (in 2.13.4 I have at hand), and see a few
strange things.

  A) Why this bruhaha with special-casing primes near 2^63 in
     u_forprime_next() (in language/forprime.c)?  This should
     slow-down the general case — and I cannot immediately see any bug
     this would work around…

     Not even mentioning how (prime) T->p can be equal to 1<<63 (=HIGHBIT)…

>	 switch(T->p)
>	 {
>   #define retp(x) return T->p = (HIGHBIT+x <= T->b)? HIGHBIT+x: 0
>	   case HIGHBIT: retp(29);
>	   case HIGHBIT + 29: retp(99);
>	   case HIGHBIT + 99: retp(123);
>	   case HIGHBIT +123: retp(131);
>	   case HIGHBIT +131: retp(155);

  B) unextprime() (the ulong version of nextprime() from
     basemath/ifactor1.c) has a rather convoluted code to sieve mod
     210.  Why not:

      • just store tables with a bitmap of a few (say 8) next odd
      	residues mod N which are mutually prime with N.  (Indexed by
      	“an odd” residue mod N.)
      • For a given odd number, extract the corresponding bitmaps from
      	the tables for a few values of N (such as 210, 11*19, 13*17, 23 and 29).
      • “bit-AND” these bitmaps.
      • Use the count of trailing/leading 0 bits (should be cheap) to
      	jump ahead.

     I would not be surprised if this can speedup the code by 20%
     (maybe even a bit more for numbers closer to 2^64).

Thanks,
Ilya