Re: mini-gmp mpz_gcdext Bézout coefficients do not match documentation

2024-02-20 Thread Guido Vranken
I've confirmed with my fuzzer that the committed patch resolves the issues.
Thanks.

On Sun, Feb 18, 2024 at 8:27 PM Niels Möller  wrote:

> marco.bodr...@tutanota.com writes:
>
> > The only reason why I prefer my solution is: when cmp<0, there is no
> need to compute
> > mpz_sub (t1, t0, t1);
>
> That's certainly a micro optimization, but let's keep things simple.
>
> I've pushed this fix now, I think there were only comment and ChangeLog
> changes since the previous version of the patch.
>
> Thanks for the review,
> /Niels
>
> --
> Niels Möller. PGP key CB4962D070D77D7FCB8BA36271D8F1FF368C6677.
> Internet email is subject to wholesale government surveillance.
>
___
gmp-bugs mailing list
gmp-bugs@gmplib.org
https://gmplib.org/mailman/listinfo/gmp-bugs


Re: mini-gmp mpz_gcdext Bézout coefficients do not match documentation

2024-02-18 Thread Niels Möller
marco.bodr...@tutanota.com writes:

> The only reason why I prefer my solution is: when cmp<0, there is no need to 
> compute
> mpz_sub (t1, t0, t1);

That's certainly a micro optimization, but let's keep things simple.

I've pushed this fix now, I think there were only comment and ChangeLog
changes since the previous version of the patch.

Thanks for the review,
/Niels

-- 
Niels Möller. PGP key CB4962D070D77D7FCB8BA36271D8F1FF368C6677.
Internet email is subject to wholesale government surveillance.
___
gmp-bugs mailing list
gmp-bugs@gmplib.org
https://gmplib.org/mailman/listinfo/gmp-bugs


Re: mini-gmp mpz_gcdext Bézout coefficients do not match documentation

2024-02-18 Thread marco . bodrato
Ciao Niels,

18 feb 2024, 10:10 da ni...@lysator.liu.se:

> marco.bodr...@tutanota.com writes:
>
>> What about;
>>
>>
>> /* Arrange so that |s| < |u| / 2g */
>>   mpz_add (s1, s0, s1);
>>   {
>>     int cmp = mpz_cmpabs (s0, s1);
>>     if (cmp >= 0)
>>   {
>>     mpz_swap (s0, s1);
>>     mpz_sub (t1, t0, t1);
>>     if (cmp != 0 || mpz_cmpabs (t0, t1) > 0)
>>   mpz_swap (t0, t1);
>>   }
>>   }
>>
>
> I think swapping of s and t must go together? I've tried out this
> similar variant 
>

You are right, it's better with the two swaps together.


>  mpz_add (s1, s0, s1);
>  mpz_sub (t1, t0, t1);
>  cmp = mpz_cmpabs (s0, s1);
>  if (cmp > 0 || (cmp == 0 && mpz_cmpabs (t0, t1) > 0))
>  {
>  mpz_swap (s0, s1);
>  mpz_swap (t0, t1);
>  }
>

The only reason why I prefer my solution is: when cmp<0, there is no need to 
compute
mpz_sub (t1, t0, t1);


> Seems to work fine, and it's a rather clear way to express the "lexical"
> preference that we want the cofactor pair with the smallest of |s|, |s'|, but
> in case of a tie, we choose the pair with smallest |t|, |t'|.
>

I agree.

Ĝis,
m
___
gmp-bugs mailing list
gmp-bugs@gmplib.org
https://gmplib.org/mailman/listinfo/gmp-bugs


Re: mini-gmp mpz_gcdext Bézout coefficients do not match documentation

2024-02-18 Thread Niels Möller
marco.bodr...@tutanota.com writes:

> s==0 is not a special case, it's one of the cases with smaller |s|.

Right, s == 0 is not a special case in itself, what's special is the
case g = a = b, then we have too candidate cofactor pairs, s=1, t=0,
s'=0, t'= 1. And we have the asymmetric preference for the first.

> What about;
>
> /* Arrange so that |s| < |u| / 2g */
>   mpz_add (s1, s0, s1);
>   {
>     int cmp = mpz_cmpabs (s0, s1);
>     if (cmp >= 0)
>   {
>     mpz_swap (s0, s1);
>     mpz_sub (t1, t0, t1);
>     if (cmp != 0 || mpz_cmpabs (t0, t1) > 0)
>   mpz_swap (t0, t1);
>   }
>   }

I think swapping of s and t must go together? I've tried out this
similar variant 

  mpz_add (s1, s0, s1);
  mpz_sub (t1, t0, t1);
  cmp = mpz_cmpabs (s0, s1);
  if (cmp > 0 || (cmp == 0 && mpz_cmpabs (t0, t1) > 0))
{
  mpz_swap (s0, s1);
  mpz_swap (t0, t1);
}

Seems to work fine, and it's a rather clear way to express the "lexical"
preference that we want the cofactor pair with the smallest of |s|, |s'|, but
in case of a tie, we choose the pair with smallest |t|, |t'|.

Updated version of the patch appended below.

Regrads,
/Niels

diff -r 1c566525476a mini-gmp/mini-gmp.c
--- a/mini-gmp/mini-gmp.c   Wed Dec 27 19:35:36 2023 +0100
+++ b/mini-gmp/mini-gmp.c   Sun Feb 18 10:02:56 2024 +0100
@@ -2809,6 +2809,7 @@
   mpz_t tu, tv, s0, s1, t0, t1;
   mp_bitcnt_t uz, vz, gz;
   mp_bitcnt_t power;
+  int cmp;
 
   if (u->_mp_size == 0)
 {
@@ -2960,12 +2961,21 @@
   mpz_tdiv_q_2exp (t0, t0, 1);
 }
 
-  /* Arrange so that |s| < |u| / 2g */
+  /* Choose small cofactors (they should satify
+
+   |s| < |u| / 2g and |t| < |v| / 2g,
+
+ except for documented exceptions). Always choose the smallest s,
+ if there are two choices for s with same absolute value, choose
+ the one with smallest corresponding t (this asymmetryic condition
+ is needed to prefer s = 0, |t| = 1 when g = |a| = |b|). */
   mpz_add (s1, s0, s1);
-  if (mpz_cmpabs (s0, s1) > 0)
+  mpz_sub (t1, t0, t1);
+  cmp = mpz_cmpabs (s0, s1);
+  if (cmp > 0 || (cmp == 0 && mpz_cmpabs (t0, t1) > 0))
 {
   mpz_swap (s0, s1);
-  mpz_sub (t0, t0, t1);
+  mpz_swap (t0, t1);
 }
   if (u->_mp_size < 0)
 mpz_neg (s0, s0);
diff -r 1c566525476a mini-gmp/tests/t-gcd.c
--- a/mini-gmp/tests/t-gcd.cWed Dec 27 19:35:36 2023 +0100
+++ b/mini-gmp/tests/t-gcd.cSun Feb 18 10:02:56 2024 +0100
@@ -55,6 +55,10 @@
   if (mpz_sgn (g) <= 0)
 return 0;
 
+  /* Require that s==0 iff g==abs(b) */
+  if (!mpz_sgn (s) != !mpz_cmpabs (g, b))
+goto fail;
+
   mpz_init (ta);
   mpz_init (tb);
   mpz_init (r);
@@ -79,20 +83,20 @@
   if (mpz_sgn (r) != 0)
 goto fail;
 
-  /* Require that 2 |s| < |b/g|, or |s| == 1. */
-  if (mpz_cmpabs_ui (s, 1) > 0)
+  /* Require that 2 |s| < |b/g|, or s == sgn(a) */
+  if (mpz_cmp_si (s, mpz_sgn (a)) != 0)
 {
   mpz_mul_2exp (r, s, 1);
-  if (mpz_cmpabs (r, tb) > 0)
+  if (mpz_cmpabs (r, tb) >= 0)
goto fail;
 }
 
-  /* Require that 2 |t| < |a/g| or |t| == 1*/
-  if (mpz_cmpabs_ui (t, 1) > 0)
+  /* Require that 2 |t| < |a/g| or t == sgn(b) */
+  if (mpz_cmp_si (t, mpz_sgn (b)) != 0)
 {
   mpz_mul_2exp (r, t, 1);
-  if (mpz_cmpabs (r, ta) > 0)
-   return 0;
+  if (mpz_cmpabs (r, ta) >= 0)
+   goto fail;
 }
 
   mpz_clear (ta);
@@ -102,17 +106,53 @@
   return 1;
 }
 
+static void test_one(const mpz_t a, const mpz_t b)
+{
+  mpz_t g, s, t;
+
+  mpz_init (g);
+  mpz_init (s);
+  mpz_init (t);
+
+  mpz_gcdext (g, s, t, a, b);
+  if (!gcdext_valid_p (a, b, g, s, t))
+{
+  fprintf (stderr, "mpz_gcdext failed:\n");
+  dump ("a", a);
+  dump ("b", b);
+  dump ("g", g);
+  dump ("s", s);
+  dump ("t", t);
+  abort ();
+}
+
+  mpz_gcd (s, a, b);
+  if (mpz_cmp (g, s))
+{
+  fprintf (stderr, "mpz_gcd failed:\n");
+  dump ("a", a);
+  dump ("b", b);
+  dump ("r", g);
+  dump ("ref", s);
+  abort ();
+}
+
+  mpz_clear (g);
+  mpz_clear (s);
+  mpz_clear (t);
+}
+
 void
 testmain (int argc, char **argv)
 {
   unsigned i;
-  mpz_t a, b, g, s, t;
+  mpz_t a, b, g, s;
+  int ai, bi;
 
   mpz_init (a);
   mpz_init (b);
   mpz_init (g);
   mpz_init (s);
-  mpz_init (t);
 
   for (i = 0; i < COUNT; i++)
 {
@@ -129,6 +169,15 @@
}
 }
 
+  /* Exhaustive test of small inputs */
+  for (ai = -30; ai <= 30; ai++)
+for (bi = -30; bi <= 30; bi++)
+  {
+   mpz_set_si (a, ai);
+   mpz_set_si (b, bi);
+   test_one (a, b);
+  }
+
   for (i = 0; i < COUNT; i++)
 {
   unsigned flags;
@@ -147,32 +196,11 @@
   if (flags & 2)
mpz_neg (b, b);
 
-  mpz_gcdext (g, s, t, a, b);
-  if (!gcdext_valid_p (a, b, g, s, t))
-   {
- fprintf (stderr, "mpz_gcdext failed:\n");
- dump ("a", a);
- dump ("b", b);
- dump ("g", g);
- dump ("s", s);
- dump ("t", 

Re: mini-gmp mpz_gcdext Bézout coefficients do not match documentation

2024-02-17 Thread marco . bodrato
Ciao Niels,

17 feb 2024, 21:26 da ni...@lysator.liu.se:

> marco.bodr...@tutanota.com writes:
>
>> Maybe we should also add the test:
>>
>>   /* Require that s==0 iff g==abs(b) */
>>   if (!mpz_sgn (s) ^ !mpz_cmpabs (g, b))
>>   goto fail;
>>
>
> Good point, that's better than only checking the special case |a| ==
> |b|. (But maybe more readable with != instead of ^).
>

Yes, with != it's more readable!


> To get mini-gmp to conform, I find no simpler way than special casing s
> == 0, like
>
s==0 is not a special case, it's one of the cases with smaller |s|.

What about;

/* Arrange so that |s| < |u| / 2g */
  mpz_add (s1, s0, s1);
  {
    int cmp = mpz_cmpabs (s0, s1);
    if (cmp >= 0)
  {
    mpz_swap (s0, s1);
    mpz_sub (t1, t0, t1);
    if (cmp != 0 || mpz_cmpabs (t0, t1) > 0)
  mpz_swap (t0, t1);
  }
  }

Ĝis,
m
___
gmp-bugs mailing list
gmp-bugs@gmplib.org
https://gmplib.org/mailman/listinfo/gmp-bugs


Re: mini-gmp mpz_gcdext Bézout coefficients do not match documentation

2024-02-17 Thread Niels Möller
marco.bodr...@tutanota.com writes:

> Maybe we should also add the test:
>
>   /* Require that s==0 iff g==abs(b) */
>   if (!mpz_sgn (s) ^ !mpz_cmpabs (g, b))
>   goto fail;

Good point, that's better than only checking the special case |a| ==
|b|. (But maybe more readable with != instead of ^).

To get mini-gmp to conform, I find no simpler way than special casing s
== 0, like

  if (mpz_sgn (s0) != 0)
{
  /* Arrange so that |s| < |u| / 2g and |t| < |v| / 2g, if possible. */
  mpz_add (s1, s0, s1);
  mpz_sub (t1, t0, t1);
  if (mpz_cmpabs (s0, s1) > 0 || mpz_cmpabs (t0, t1) > 0)
{
  mpz_swap (s0, s1);
  mpz_swap (t0, t1);
}
}

And I think this is the only condition that isn't fully symmetrical, so
perhaps not surprising if it needs a special case.

Regards,
/Niels

-- 
Niels Möller. PGP key CB4962D070D77D7FCB8BA36271D8F1FF368C6677.
Internet email is subject to wholesale government surveillance.
___
gmp-bugs mailing list
gmp-bugs@gmplib.org
https://gmplib.org/mailman/listinfo/gmp-bugs


Re: mini-gmp mpz_gcdext Bézout coefficients do not match documentation

2024-02-17 Thread marco . bodrato
Ciao,

17 feb 2024, 19:15 da :

> Ciao Niels,
>
> 17 feb 2024, 17:53 da ni...@lysator.liu.se:
>
>> Niels Möller  writes:
>> The documented conditions say that gmp should return the same cofactors
>>
>
> The documentation also says "If abs(a) = abs(b), then s = 0, t = sgn(b)."
>
Maybe we should also add the test:

  /* Require that s==0 iff g==abs(b) */
  if (!mpz_sgn (s) ^ !mpz_cmpabs (g, b))
  goto fail;

Ĝis,
m
___
gmp-bugs mailing list
gmp-bugs@gmplib.org
https://gmplib.org/mailman/listinfo/gmp-bugs


Re: mini-gmp mpz_gcdext Bézout coefficients do not match documentation

2024-02-17 Thread Niels Möller
marco.bodr...@tutanota.com writes:

> And if I correctly patched and tested your proposed code. with equal numbers 
> I get t=0, instead of s=0.

Thanks for testing, I have to look into that case, then.

> Shat about simply changing the test from > to >= ?
>
>    /* Arrange so that |s| < |u| / 2g */
>    mpz_add (s1, s0, s1);
> -  if (mpz_cmpabs (s0, s1) > 0)
> +  if (mpz_cmpabs (s0, s1) >= 0)

I don't think that's right. My understanding is that if |s| = |s'|
(notation from my mail, the variables in the code are reassigned several
times, which is a bit confusing), then we ncessarily have |t| != |t'|,
and we must choose candidate to get the smaller one.

Regards,
/Niels

-- 
Niels Möller. PGP key CB4962D070D77D7FCB8BA36271D8F1FF368C6677.
Internet email is subject to wholesale government surveillance.
___
gmp-bugs mailing list
gmp-bugs@gmplib.org
https://gmplib.org/mailman/listinfo/gmp-bugs


Re: mini-gmp mpz_gcdext Bézout coefficients do not match documentation

2024-02-17 Thread marco . bodrato
Ciao Niels,

17 feb 2024, 17:53 da ni...@lysator.liu.se:

> Niels Möller  writes:
> The documented conditions say that gmp should return the same cofactors
>

The documentation also says "If abs(a) = abs(b), then s = 0, t = sgn(b)."


>
> canonicalization logic at the end wasn't right. Appending a patch with a
> fix + stricter tests.
>

And if I correctly patched and tested your proposed code. with equal numbers I 
get t=0, instead of s=0.


> Without the fix, we could make the wrong choice in case
>
>  |s'| = |s| = |b| / 2g
>
> Since it's a bit subtle, it would be nice with a review of this code
> before I commit it.
>

Shat about simply changing the test from > to >= ?

   /* Arrange so that |s| < |u| / 2g */
   mpz_add (s1, s0, s1);
-  if (mpz_cmpabs (s0, s1) > 0)
+  if (mpz_cmpabs (s0, s1) >= 0)
 {
   mpz_swap (s0, s1);
   mpz_sub (t0, t0, t1);

Ĝis,
m
___
gmp-bugs mailing list
gmp-bugs@gmplib.org
https://gmplib.org/mailman/listinfo/gmp-bugs


Re: mini-gmp mpz_gcdext Bézout coefficients do not match documentation

2024-02-17 Thread Niels Möller
Niels Möller  writes:

> Guido Vranken  writes:
>
>> Computing extended gcd using mpz_gcdext where a = 1, b = 2:
>>
>> libgmp: g: 1, s: 1, t: 0
>> mini-gmp: g: 1, s: -1, t: 1
>> This violates the docs: "s = sgn(a) if abs(b) = 2", e.g. s must be 1
>>
>> Computing extended gcd using mpz_gcdext where a = 6, b = 4:
>> libgmp: g: 2, s: 1, t: -1
>> mini-libgmp: g: 2, s: -1, t: 2
>
> I think mini-gmp is wrong here, I'll have to investigate. Thanks for
> reporting.

The documented conditions say that gmp should return the same cofactors
as one would get naturally with Euclid's algorithm. It appears mini-gmp
uses binary gcdext (I had forgotten that detail), and the
canonicalization logic at the end wasn't right. Appending a patch with a
fix + stricter tests.

At this place in the algorithm, we have a, b >= 0, and we have a
candidate s, t, with

  -b / g <= s <= 0, 0 <= t <= a / g

(not sure at which edges we might have equality), and then we construct
another candidate s', t' as

  s' = s + b / g >= 0, t' = t - a / g <= 0

And we then return s', t' if either |s'| < |s| or |t'| < |t|.

Without the fix, we could make the wrong choice in case

  |s'| = |s| = |b| / 2g

Since it's a bit subtle, it would be nice with a review of this code
before I commit it.

Regards,
/Niels

diff -r 1c566525476a mini-gmp/mini-gmp.c
--- a/mini-gmp/mini-gmp.c   Wed Dec 27 19:35:36 2023 +0100
+++ b/mini-gmp/mini-gmp.c   Sat Feb 17 17:17:21 2024 +0100
@@ -2960,12 +2960,13 @@
   mpz_tdiv_q_2exp (t0, t0, 1);
 }
 
-  /* Arrange so that |s| < |u| / 2g */
+  /* Arrange so that |s| < |u| / 2g and |t| < |v| / 2g, if possible. */
   mpz_add (s1, s0, s1);
-  if (mpz_cmpabs (s0, s1) > 0)
+  mpz_sub (t1, t0, t1);
+  if (mpz_cmpabs (s0, s1) > 0 || mpz_cmpabs (t0, t1) > 0)
 {
   mpz_swap (s0, s1);
-  mpz_sub (t0, t0, t1);
+  mpz_swap (t0, t1);
 }
   if (u->_mp_size < 0)
 mpz_neg (s0, s0);
diff -r 1c566525476a mini-gmp/tests/t-gcd.c
--- a/mini-gmp/tests/t-gcd.cWed Dec 27 19:35:36 2023 +0100
+++ b/mini-gmp/tests/t-gcd.cSat Feb 17 17:17:21 2024 +0100
@@ -79,20 +79,20 @@
   if (mpz_sgn (r) != 0)
 goto fail;
 
-  /* Require that 2 |s| < |b/g|, or |s| == 1. */
-  if (mpz_cmpabs_ui (s, 1) > 0)
+  /* Require that 2 |s| < |b/g|, or s == sgn(a) */
+  if (mpz_cmp_si (s, mpz_sgn (a)) != 0)
 {
   mpz_mul_2exp (r, s, 1);
-  if (mpz_cmpabs (r, tb) > 0)
+  if (mpz_cmpabs (r, tb) >= 0)
goto fail;
 }
 
-  /* Require that 2 |t| < |a/g| or |t| == 1*/
-  if (mpz_cmpabs_ui (t, 1) > 0)
+  /* Require that 2 |t| < |a/g| or t == sgn(b) */
+  if (mpz_cmp_si (t, mpz_sgn (b)) != 0)
 {
   mpz_mul_2exp (r, t, 1);
-  if (mpz_cmpabs (r, ta) > 0)
-   return 0;
+  if (mpz_cmpabs (r, ta) >= 0)
+   goto fail;
 }
 
   mpz_clear (ta);
@@ -102,17 +102,53 @@
   return 1;
 }
 
+static void test_one(const mpz_t a, const mpz_t b)
+{
+  mpz_t g, s, t;
+
+  mpz_init (g);
+  mpz_init (s);
+  mpz_init (t);
+
+  mpz_gcdext (g, s, t, a, b);
+  if (!gcdext_valid_p (a, b, g, s, t))
+{
+  fprintf (stderr, "mpz_gcdext failed:\n");
+  dump ("a", a);
+  dump ("b", b);
+  dump ("g", g);
+  dump ("s", s);
+  dump ("t", t);
+  abort ();
+}
+
+  mpz_gcd (s, a, b);
+  if (mpz_cmp (g, s))
+{
+  fprintf (stderr, "mpz_gcd failed:\n");
+  dump ("a", a);
+  dump ("b", b);
+  dump ("r", g);
+  dump ("ref", s);
+  abort ();
+}
+
+  mpz_clear (g);
+  mpz_clear (s);
+  mpz_clear (t);
+}
+
 void
 testmain (int argc, char **argv)
 {
   unsigned i;
-  mpz_t a, b, g, s, t;
+  mpz_t a, b, g, s;
+  int ai, bi;
 
   mpz_init (a);
   mpz_init (b);
   mpz_init (g);
   mpz_init (s);
-  mpz_init (t);
 
   for (i = 0; i < COUNT; i++)
 {
@@ -129,6 +165,15 @@
}
 }
 
+  /* Exhaustive test of small inputs */
+  for (ai = -30; ai <= 30; ai++)
+for (bi = -30; bi <= 30; bi++)
+  {
+   mpz_set_si (a, ai);
+   mpz_set_si (b, bi);
+   test_one (a, b);
+  }
+
   for (i = 0; i < COUNT; i++)
 {
   unsigned flags;
@@ -147,32 +192,11 @@
   if (flags & 2)
mpz_neg (b, b);
 
-  mpz_gcdext (g, s, t, a, b);
-  if (!gcdext_valid_p (a, b, g, s, t))
-   {
- fprintf (stderr, "mpz_gcdext failed:\n");
- dump ("a", a);
- dump ("b", b);
- dump ("g", g);
- dump ("s", s);
- dump ("t", t);
- abort ();
-   }
+  test_one (a, b);
+}
 
-  mpz_gcd (s, a, b);
-  if (mpz_cmp (g, s))
-   {
- fprintf (stderr, "mpz_gcd failed:\n");
- dump ("a", a);
- dump ("b", b);
- dump ("r", g);
- dump ("ref", s);
- abort ();
-   }
-}
   mpz_clear (a);
   mpz_clear (b);
   mpz_clear (g);
   mpz_clear (s);
-  mpz_clear (t);
 }
diff -r 1c566525476a tests/mpz/t-gcd.c
--- a/tests/mpz/t-gcd.c Wed Dec 27 19:35:36 2023 +0100
+++ b/tests/mpz/t-gcd.c Sat Feb 17 17:17:21 2024 +0100
@@ -440,8 +440,8 @@
   if (mpz_sgn 

Re: mini-gmp mpz_gcdext Bézout coefficients do not match documentation

2024-02-17 Thread Niels Möller
Guido Vranken  writes:

> Computing extended gcd using mpz_gcdext where a = 1, b = 2:
>
> libgmp: g: 1, s: 1, t: 0
> mini-gmp: g: 1, s: -1, t: 1
> This violates the docs: "s = sgn(a) if abs(b) = 2", e.g. s must be 1
>
> Computing extended gcd using mpz_gcdext where a = 6, b = 4:
> libgmp: g: 2, s: 1, t: -1
> mini-libgmp: g: 2, s: -1, t: 2

I think mini-gmp is wrong here, I'll have to investigate. Thanks for
reporting.

> The docs state: "abs(s) < abs(b) / (2 g) and abs(t) < abs(a) / (2 g)"
> I think that should be: "abs(s) <= abs(b) / (2 g) and abs(t) <= abs(a) / (2
> g)"
> (< should be <=)

These strict equalities don't hold for the documented "exceptional"
cases. IT seems both your examples have abs(B) = 2G, so that's the case
that applies.

It is a bit subtle. And even if one relaxes the inequalities as you
suggest, 1 = gcd(1,1) would still be an exception (since in that case,
we'd get |s|, |t| < 1/2, i.e., s = t = 0, which isn't reasonable).

Can you suggest any doc change to make it clearer?

Regards,
/Niels

-- 
Niels Möller. PGP key CB4962D070D77D7FCB8BA36271D8F1FF368C6677.
Internet email is subject to wholesale government surveillance.
___
gmp-bugs mailing list
gmp-bugs@gmplib.org
https://gmplib.org/mailman/listinfo/gmp-bugs


mini-gmp mpz_gcdext Bézout coefficients do not match documentation

2024-02-17 Thread Guido Vranken
Computing extended gcd using mpz_gcdext where a = 1, b = 2:

libgmp: g: 1, s: 1, t: 0
mini-gmp: g: 1, s: -1, t: 1
This violates the docs: "s = sgn(a) if abs(b) = 2", e.g. s must be 1

Computing extended gcd using mpz_gcdext where a = 6, b = 4:
libgmp: g: 2, s: 1, t: -1
mini-libgmp: g: 2, s: -1, t: 2

The docs state: "abs(s) < abs(b) / (2 g) and abs(t) < abs(a) / (2 g)"
I think that should be: "abs(s) <= abs(b) / (2 g) and abs(t) <= abs(a) / (2
g)"
(< should be <=)
Then, for libgmp this is true: s: 1 <= 4 / 4, t: 1 <= 6 / 4
For mini-gmp this is false for t: s: 1 <= 4 / 4, t: 2 <= 6 / 4

Reproducer:

//#include 
#include "mini-gmp.c"

#include 

static void test(
const unsigned long int a_ui,
const unsigned long int b_ui)
{
mpz_t a, b, g, s, t;
mpz_init(a);
mpz_init(b);
mpz_init(g);
mpz_init(s);
mpz_init(t);

mpz_set_ui(a, a_ui);
mpz_set_ui(b, b_ui);

mpz_gcdext(g, s, t, a, b);

printf("g: %s, s: %s, t: %s\n",
mpz_get_str(NULL, 10, g),
mpz_get_str(NULL, 10, s),
mpz_get_str(NULL, 10, t));
}

int main()
{
test(1, 2);
test(6, 4);
}

Tested on Linux x64 using the latest repository code.
___
gmp-bugs mailing list
gmp-bugs@gmplib.org
https://gmplib.org/mailman/listinfo/gmp-bugs