mpz_prevprime

[Seth Troisi]

On Thu, Oct 15, 2020 at 12:13 AM Niels Möller <nisse at lysator.liu.se> wrote:

> Seth Troisi <braintwo at gmail.com> writes:
>
> > I modified the patch a tiny bit. Still hoping to get this in for an
> > upcoming prime-gap search project.
>
> Looks pretty good to me.
>
I updated and attached a new patch with the tests discussed on Oct 3rd

>
> > +static int
> > +findnext (mpz_ptr p,
> > +          unsigned long(*nextmod_func)(const mpz_t, unsigned long),
> > +          void(*nextseq_func)(mpz_t, const mpz_t, unsigned long))
>
> I'd name the function pointer arguments to be more similar to functions
> they refer to, and without the "next" and "_func" parts. Maybe mod_ui
> and add_ui/incr_ui/update_ui?
>
> And is there a good reason you need different mod functions? I haven't
> been following along closely, so I don't know why the current code uses
> mpz_cdiv_ui rather than mpz_fdiv_ui. It would make things a bit simpler
> if we could use mpz_fdiv_ui (the standard mathematical mod operation)
> always.
>

I use cdiv in next_prime to find the offset to the next multiple of a prime
<=> (-n) % p, for previous prime this is n % p (e.g. fdiv)

for nextprime sieving the interval 1,000,000 to 1,000,500
p is the start of the interval (1,000,000) in this case
let prime = 401
(-p) % prime = mpz_cdiv_ui(p, prime) = 94
the 94th number in the sieve (p + 94 = 100094) is the first number
divisible by 401 in the interval

for previous prime we can reuse almost all of the logic but with fdiv_ui
sieving the interval 1,000,000 to 999,500
set p to the "start" of the interval (1,000,000) in this case
p % prime = mpz_fdiv_ui(p, prime) = 307
the 307th number in the sieve (p - 307 = 999693) is the first number
divisible by 401 in the interval.

Technically I can restructure the code to make both use with mpz_fdiv_ui
My quick and dirty attempt

diff -r 1e638b839cbe mpz/nextprime.c
--- a/mpz/nextprime.c   Thu Oct 15 01:13:43 2020 -0700
+++ b/mpz/nextprime.c   Thu Oct 15 01:30:03 2020 -0700
@@ -156,7 +156,7 @@

 static int
 findnext (mpz_ptr p,
-          unsigned long(*nextmod_func)(const mpz_t, unsigned long),
+          char direction,
           void(*nextseq_func)(mpz_t, const mpz_t, unsigned long))
 {
   char *composite;
@@ -234,11 +234,17 @@
       unsigned long incr, prime;
       int primetest;

+      /* set p to the top of the interval */
+      if (direction == 1) {
+        mpz_add_ui(p, p, 2 * (odds_in_composite_sieve-1));
+      }
+
       memset (composite, 0, odds_in_composite_sieve);
       prime = 3;
       for (i = 0; i < prime_limit; i++)
         {
-          m = nextmod_func(p, prime);
+          /* m = p % prime which is the distance from top of sieve
interval to previous multiple of p */
+          m = mpz_fdiv_ui(p, prime);
           /* Only care about odd multiplies of prime. */
           if (m & 1)
             m += prime;
@@ -250,10 +256,15 @@
           prime += primegap[i];
         }

+      /* reset p to old value */
+      if (direction == 1) {
+        mpz_sub_ui(p, p, 2 * (odds_in_composite_sieve-1));
+      }
+
       difference = 0;
       for (incr = 0; incr < odds_in_composite_sieve; difference += 2, incr
+= 1)
         {
-          if (composite[incr])
+          if (direction < 0 ? composite[incr] :
composite[odds_in_composite_sieve - 1 - incr])
             continue;

           nextseq_func(p, p, difference);
@@ -288,7 +299,7 @@
   mpz_add_ui (p, n, 1);
   mpz_setbit (p, 0);

-  findnext(p, mpz_cdiv_ui, mpz_add_ui);
+  findnext(p, 1, mpz_add_ui);
 }

 int
@@ -309,7 +320,7 @@
   mpz_sub_ui (p, n, 2);
   mpz_setbit (p, 0);

-  return findnext(p, mpz_fdiv_ui, mpz_sub_ui);
+  return findnext(p, -1, mpz_sub_ui);
 }


> > +int
> > +mpz_prevprime (mpz_ptr p, mpz_srcptr n)
>
> Interface looks good to me. And if we later add a function to find the
> first prime in an arithmetic progression, I think that should fit well.
>
> Regards,
> /Niels
>
> --
> Niels Möller. PGP-encrypted email is preferred. Keyid 368C6677.
> Internet email is subject to wholesale government surveillance.
>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: prevprime_20201015.patch
Type: text/x-patch
Size: 16680 bytes
Desc: not available
URL: <https://gmplib.org/list-archives/gmp-devel/attachments/20201015/0203c2c5/attachment-0001.bin>