mpz_prevprime

Tue Sep 3 06:59:08 UTC 2019

Trying to sum up the proposals for how this api should work,

1. add mpz_prevprime(rop, op) or mpz_precprime(rop, op)
2. change mpz_nextprime(-10) => -7 (instead of 2)
3. add mpz_nextprime_seq(a, b) => next prime a + k*b with k >= 1

Considerations
1. on negative number return X?
Should prevprime return "input had correct sign"? should nextprime?
2. this is an incompatible change for the interface (but probably not one
many people depenedend on)
3. a bigger change than other two proposals.

---

Personally I'm a fan of number 1. because I can clearly see the steps
needed.
I'm wary of time it would take for me to understand and code number 3. to
the
quality that it would be submitted.

On Mon, Sep 2, 2019 at 11:41 PM Seth Troisi <braintwo at gmail.com> wrote:

> I've attached an updated diff, thanks for the comments Marco.
>
> There doesn't seem to be a way to set s->reps without setting a high
> precision so
> I reused SPEED_BLOCK_SIZE from speed_gcd and ran the inner loop that many
> times
>
> Would folks rather I time
> 1. s->reps next_primes chained (next_prime(p, p); next_prime(p, p);)
> 2. s->reps of next_prime on random s->size bit numbers?
>
> 1. will be slightly slower (as it will have to search the full gap each
> time) and more likely
> for intermediate results to overflow s->size bits (which is seems fine)
> than 2.
>
> Currently it takes about 20 seconds to test -s 16-128, and 15 minutes to
> test -s 16-1028
> with -t 1.04 the resulting graph is much smoother than the old graph.
>
>
> On Sat, Aug 24, 2019 at 8:33 AM Marco Bodrato <bodrato at mail.dm.unipi.it>
> wrote:
>
>> Ciao,
>>
>> Il Gio, 22 Agosto 2019 8:48 am, Marco Bodrato ha scritto:
>> > Il Gio, 22 Agosto 2019 7:24 am, Niels Möller ha scritto:
>> >> Maybe it would make sense to generalize a bit further, and write a
>> >> function to find the first prime (if any) in an arithmetic sequence.
>> >> I.e., find the smallest k >= such that A + kB is prime (if any). A
>>
>> ... here too we should decide what to do with negative numbers.
>> which prime should we return if:
>>  1) A = -123; B= 42; (I'd suggest A+3*B = +3)
>>  2) A = -123; B= 30; (I'd suggest A+4*B = -3)
>>  3) A = -123; B= 26; (I'd suggest A+1*B = -97)
>>
>> > I know that my proposal is an incompatible change for the interface...
>> but
>> > I'd suggest to extend the range for mpz_nextprime. Can we consider both
>> > positive and negative primes?
>>
>> By the way... our documentation does not mention "positive" primes :-) and
>> I'm quite sure that there aren't programs using the fact that the current
>> implementation of nextprime returns 2 for negative numbers...
>>
>> Il Gio, 22 Agosto 2019 10:32 am, Niels Möller ha scritto:
>> > I think that's a bit unusual. We can compare with the corresponding
>> > functions in pari/gp. These differ from gmp by nextprime(p) == p if p is
>> > prime. And they don't consider negative numbers to be primes.
>>
>> Il Mer, 21 Agosto 2019 10:02 am, Seth Troisi ha scritto:
>> > Finally for mpz_prevprime(rop, op <= 2) it's not clear what rop should
>> be
>> > set to, so I use the smallest prime (2).
>>
>> Well, in pari/gp there is a function
>> precprime(x): largest pseudoprime <= x, 0 if x<=1.
>>
>> I'm not sure I like that interface...
>>
>> My proposal to extend the range of nextprime to negative numbers:
>>
>> diff -r 643c931da9bd mpz/nextprime.c
>> --- a/mpz/nextprime.c   Thu Aug 22 15:14:09 2019 +0200
>> +++ b/mpz/nextprime.c   Sat Aug 24 16:52:21 2019 +0200
>> @@ -56,21 +56,24 @@
>>    mp_size_t pn;
>>    mp_bitcnt_t nbits;
>>    unsigned incr;
>> +  mpz_t tp;
>>    TMP_SDECL;
>>
>>    /* First handle tiny numbers */
>> -  if (mpz_cmp_ui (n, 2) < 0)
>> +  if (mpz_cmp_ui (n, 2) < 0 && mpz_cmpabs_ui (n, 4) < 0)
>>      {
>>        mpz_set_ui (p, 2);
>> +      if (mpz_cmpabs_ui (n, 3) == 0)
>> +       mpz_neg (p, p);
>>        return;
>>      }
>>    mpz_add_ui (p, n, 1);
>>    mpz_setbit (p, 0);
>>
>> -  if (mpz_cmp_ui (p, 7) <= 0)
>> +  if (mpz_cmpabs_ui (p, 7) <= 0)
>>      return;
>>
>> -  pn = SIZ(p);
>> +  pn = ABSIZ(p);
>>    MPN_SIZEINBASE_2EXP(nbits, PTR(p), pn, 1);
>>    if (nbits / 2 >= NUMBER_OF_PRIMES)
>>      prime_limit = NUMBER_OF_PRIMES - 1;
>> @@ -88,7 +91,7 @@
>>        prime = 3;
>>        for (i = 0; i < prime_limit; i++)
>>         {
>> -         moduli[i] = mpz_tdiv_ui (p, prime);
>> +         moduli[i] = mpz_fdiv_ui (p, prime);
>>           prime += primegap[i];
>>         }
>>
>> @@ -115,7 +118,7 @@
>>           difference = 0;
>>
>>           /* Miller-Rabin test */
>> -         if (mpz_millerrabin (p, 25))
>> +         if (mpz_millerrabin (mpz_roinit_n (tp, PTR (p), ABSIZ (p)), 25))
>>             goto done;
>>         next:;
>>           incr += 2;
>>
>>
>> After those changes, the _precprime function can be implemented with:
>>
>> void
>> mpz_precprime (mpz_ptr p, mpz_srcptr n)
>> {
>>   mpz_t tn;
>>
>>   mpz_nextprime (p, mpz_roinit_n (tn, PTR (n), -SIZ (n)));
>>   mpz_neg (p, p);
>> }
>>
>>
>> Ĝis,
>> m
>>
>> --
>> http://bodrato.it/papers/
>>
>>