More accurate calculation of sieve array size

[Adrien Prost-Boucle]

Hi,

Thank you for the clarifications.
I was wrong indeed for the oddfac implementation.

For primorial, when thinking on the code I considered that the call to
  mpz_prodlimbs (x, factors, j)
is the one responsible for (re)sizing x to whatever appropriate for the result.
(and I think I got confused with the reuse of name x in nested loop)
If sizing x correctly before the call mpz_prodlimbs() garantees that mpz_prodlimbs()
will not resize it, then the current implementation is efficient indeed.
This could have deserved a bit of comment ;-)

Thanks,
Adrien


On Thu, 2022-08-11 at 16:32 +0200, Marco Bodrato wrote:
> Ciao Adrien,
> 
> Il 2022-08-10 15:37 Adrien Prost-Boucle ha scritto:
> > Looking into oddfac and primorial computation code,
> > I think there is a suboptimality in the computation of the prime sieve 
> > size.
> > Could anyone confirm?
> 
> Maybe there are some too-conservative assumptions, but we must be 
> careful with
> the code.
> 
> > --- mpz/oddfac_1.old    2022-08-10 15:14:04.511298108 +0200
> > +++ mpz/oddfac_1.c      2022-08-10 15:31:34.721354647 +0200
> > @@ -381,8 +381,8 @@
> >           TMP_MARK;
> > 
> >           flag--;
> > -         size = n / GMP_NUMB_BITS + 4;
> > -         ASSERT (primesieve_size (n - 1) <= size - (size / 2 + 1));
> > +         size = (n / 3 / GMP_NUMB_BITS + 1) * 2;
> > +         ASSERT (primesieve_size (n - 1) <= size / 2);
> >           /* 2-multiswing(n) < 2^(n-1)*sqrt(n/pi) < 2^(n+GMP_NUMB_BITS);
> >              one more can be overwritten by mul, another for the sieve */
> >           MPZ_TMP_INIT (mswing, size);
> 
> Here, mswing is not just used for the sieve, it will contain the result
> of the computation of the "2-multiswing.factorial of n". A portion of
> the same allocated space is used for the sieve (the ASSERT checks, is
> the size needed for the sieve less than half the size we allocate? Yes,
> ok we can go on...).
> 
> The comment says that /* 2-multiswing(n) < 2^(n+GMP_NUMB_BITS) */, then
> it makes sense to allocate n / GMP_NUMB_BITS + 4 limbs.
> 
> Surely this is an over-estimate, and maybe we can estimate it better.
> But to reduce the memory usage, we have to estimate the maximal size
> of the value that will be stored in mswing, not the sieve only.
> 
> > --- mpz/primorial_ui.c.old      2022-08-10 14:44:04.437867901 +0200
> > +++ mpz/primorial_ui.c  2022-08-10 14:55:41.721238758 +0200
> > @@ -82,8 +82,7 @@
> >        mp_limb_t prod;
> >        TMP_DECL;
> > 
> > -      size = n / GMP_NUMB_BITS;
> > -      size = size + (size >> 1) + 1;
> > +      size = n / 3 / GMP_NUMB_BITS + 1;
> >        ASSERT (size >= primesieve_size (n));
> >        sieve = MPZ_NEWALLOC (x, size);
> >        size = (gmp_primesieve (sieve, n) + 1) / log_n_max (n) + 1;
> 
> The same here. sieve is stored at the same address as x, but we need
> to estimate the size of the value that x will store, not the sieve.
> 
> x will store the result, primorial (n). The question is: is the
> allocated size reasonable to store the result x=primorial(n)?
> The other question, "Can it also contain the sieve?", is secondary.
> 
> Ĝis,
> m
>