mpz_prevprime

[Seth Troisi]

On Mon, Mar 9, 2020 at 5:03 AM Marco Bodrato <bodrato at mail.dm.unipi.it>
wrote:

> Ciao,
>
> Il 2020-03-09 11:59 Seth Troisi ha scritto:
> > On Mon, Mar 9, 2020 at 2:05 AM Marco Bodrato
> > <bodrato at mail.dm.unipi.it> wrote:
> >
> >> Ciao,
> >>
> >> Il 2020-03-09 02:56 Seth Troisi ha scritto:
>
> > It's highly possible I misunderstand how these macros are supposed to
> > be used.
> > I also dislike the boiler plate of the macros but I didn't like
>
> When I say "that code is messy", I mean my code.
> And you agree :-) those macros are just easy to misunderstand :-/
>
> > If you have suggestions or code pointers for a better pattern I'd
> > appreciate that.
>
> >> Did you try to use the gmp_primesieve function directly?
> >
> > I'm not sure what you're referencing here, I tried with
>
> I'm proposing to get rid of primegap, and use something like the
> following.
> This is just a copy-paste from your code, not a really working sequence!
>
> +      __sieve = TMP_ALLOC_LIMBS (primesieve_size (sieve_limit));
> +      prime_limit_pre = gmp_primesieve(__sieve, sieve_limit);
>
> +  next_mult = TMP_ALLOC_TYPE (prime_limit, unsigned int);
> [...]
> +  /* Invert p for modulo */
> +  mpz_neg(p, p);
>
> [..handle 3 outside the loop]
>
> +    LOOP_ON_SIEVE_BEGIN (prime, n_to_bit (5), n_to_bit (sieve_limit),
> 0, __sieve);
> +      m = mpz_fdiv_ui(p, prime);
> +      /* Only care about odd multiplies (even indexes) */
> +      if (m & 1)
> +        m += prime;
> +
> +      /* mark off any composites in sieve */
> +      for (; m < INCR_LIMIT; m += 2*prime)
> +        sieve[m] = 1;
> +      next_mult[i] = m;
> +    LOOP_ON_SIEVE_END;
> +  mpz_neg(p, p);
>
> [...]
>
> +      /* Sieve next segment */
> [..handle 3]
> +      memset(sieve, 0, INCR_LIMIT);
> +      LOOP_ON_SIEVE_BEGIN (prime, n_to_bit (5), n_to_bit (sieve_limit),
> 0, __sieve);
> +          m = next_mult[i] - INCR_LIMIT;
> +          for (; m < INCR_LIMIT; m += prime * 2)
> +            sieve[m] = 1;
> +          next_mult[i] = m;
> +      LOOP_ON_SIEVE_END;
>
When nbits is small (IMHO the most common case) we'd like to be able to
fallback to a constant array which I don't think exist for primesieve
(because gmp_limb_bits isn't constant).
So I lowered the threshold to 150M (50MB) in the future I have a trick that
will reduce memory usage and we can try to increase this limit.


> > Is there a way to unalloc a TMP_ALLOC_LIMBS before TMP_FREE?
>
> No, if you want to unalloc, you must use __GMP_ALLOCATE_FUNC_TYPE,
> __GMP_REALLOCATE_FUNC_TYPE, and __GMP_FREE_FUNC_TYPE.
>
> > Sieve next segment is rare, the average gap is ln(n), if INCR_LIMIT
> > was changed to be a variable it could be set
> > so /*sieve next segment*/ happened on average less than once, and then
>
> Why not using a variable for INCR_LIMIT?
>
It will be another change with more moving parts. I'm hoping to do it after
the large part of this change goes in.

>
> >> In the array char *sieve, you use only the even entries. Maybe you
> >> can reduce its size by half, and remove some "*2" around the code?
> >
> > it's only 4000 entries so I wasn't bothering at this time, but it
>
> I agree.
>
Fixed now.

>
> PS: did you consider mpz_Cdiv_ui, instead of _neg,_Fdiv_ui,_neg ?
>
Nope, this is much cleaner. Thanks.

>
> Ĝis,
> m
>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: nextprime_sieve_marco.patch
Type: text/x-patch
Size: 8691 bytes
Desc: not available
URL: <https://gmplib.org/list-archives/gmp-devel/attachments/20200309/bea30c2e/attachment.bin>