mpz_gcd_ext(NULL, ...)

[Marco Bodrato]

Ciao,

Il Ven, 25 Novembre 2016 2:11 pm, Marc Glisse ha scritto:
> On Fri, 25 Nov 2016, Niels Möller wrote:
>
>> Marc Glisse <marc.glisse at inria.fr> writes:
>>> a user was asking if we could support calling mpz_gcd_ext with a NULL
>>> first argument (the gcd), since they are only interested in the
>>> coefficients s and t and would like to save the unnecessary
>>> allocation. I doubt it would save that much, but it seems trivial to
>>> add a check if(g!=NULL) similar to the tests for s and t. Does it make
>>> sense to you?
>>
>> Sounds reasonable to me.
>
> Done.

The modified manual now says that "If @var{s}, @var{t} or @var{g} is
@code{NULL} then that value is not computed.", but g is computed anyway,
even if it is not returned...

While we are looking at gcdext, I propose a small change to save an
allocation when t must be computed. Can you review it, to check if the
size assumptions are safe?

*** /tmp/extdiff.sTxw3v/gmp-repo.e121607cb19e/mpz/gcdext.c	2016-12-03
02:53:21.823688086 +0100
***************
*** 77,94 ****
  	}
        return;
      }

    TMP_MARK;

!   TMP_ALLOC_LIMBS_2 (tmp_ap, asize, tmp_bp, bsize);
    MPN_COPY (tmp_ap, PTR (a), asize);
    MPN_COPY (tmp_bp, PTR (b), bsize);

-   TMP_ALLOC_LIMBS_2 (tmp_gp, bsize, tmp_sp, bsize + 1);
-
    gsize = mpn_gcdext (tmp_gp, tmp_sp, &tmp_ssize, tmp_ap, asize, tmp_bp,
bsize);

    ssize = ABS (tmp_ssize);
    tmp_ssize = SIZ (a) >= 0 ? tmp_ssize : -tmp_ssize;

    if (t != NULL)
--- 77,95 ----
  	}
        return;
      }

    TMP_MARK;

!   /* Reserve 2 more limbs for tmp_{ap,bp}, to reuse this space if t !=
NULL */
!   TMP_ALLOC_LIMBS_3 (tmp_ap, asize + bsize + 2,
! 		     tmp_gp, bsize, tmp_sp, bsize + 1);
!   tmp_bp = tmp_ap + asize;
    MPN_COPY (tmp_ap, PTR (a), asize);
    MPN_COPY (tmp_bp, PTR (b), bsize);

    gsize = mpn_gcdext (tmp_gp, tmp_sp, &tmp_ssize, tmp_ap, asize, tmp_bp,
bsize);

    ssize = ABS (tmp_ssize);
    tmp_ssize = SIZ (a) >= 0 ? tmp_ssize : -tmp_ssize;

    if (t != NULL)
***************
*** 99,111 ****
        PTR (&gtmp) = tmp_gp;
        SIZ (&gtmp) = gsize;

        PTR (&stmp) = tmp_sp;
        SIZ (&stmp) = tmp_ssize;

!       MPZ_TMP_INIT (x, ssize + asize + 1);
        mpz_mul (x, &stmp, a);
        mpz_sub (x, &gtmp, x);
        mpz_divexact (t, x, b);
      }

    if (s != NULL)
--- 100,115 ----
        PTR (&gtmp) = tmp_gp;
        SIZ (&gtmp) = gsize;

        PTR (&stmp) = tmp_sp;
        SIZ (&stmp) = tmp_ssize;

!       ASSERT (asize + bsize + 2 >= ssize + asize + 1);
!       PTR (x) = tmp_ap; /* asize + bsize + 2 */
!       ALLOC (x) = ssize + asize + 1;
!
        mpz_mul (x, &stmp, a);
        mpz_sub (x, &gtmp, x);
        mpz_divexact (t, x, b);
      }

    if (s != NULL)


Next step to clean up this function should be to change the type of gtmp
and stmp, from __mpz_struct to mpz_t. Agreed?

Best regards,
m

-- 
http://bodrato.it/papers/