Small operands gcd improvements

[Niels Möller]

tg at gmplib.org (Torbjörn Granlund) writes:

> Now pushed gcd_11.asm for all gcd_1.asm except the one for ia64 and the
> one for AMD k6.  

Excellent!

> I might convert the former, but I am tempted to simply delete the k6
> gcd_1.asm.

The k6 code uses a branch on the u - v sign? Might be slower than the
current C code.  What hardware is it used for?

I also see that there's no x86/gcd_1.asm.

Below a gcd_11 unit test. Not terribly interesting, but we'll need
something similar for testing gcd_22. And we get a place to add tests
for any problematic corner cases.

Regards,
/Niels

diff -Nprc2 gmp-gcd_11.fc77a3dc184b/tests/mpn/Makefile.am gmp-gcd_11/tests/mpn/Makefile.am
*** gmp-gcd_11.fc77a3dc184b/tests/mpn/Makefile.am	2019-08-07 17:21:54.000000000 +0200
--- gmp-gcd_11/tests/mpn/Makefile.am	2019-08-07 23:23:35.549830000 +0200
*************** check_PROGRAMS = t-asmtype t-aors_1 t-di
*** 30,34 ****
    t-div t-mul t-mullo t-sqrlo t-mulmod_bnm1 t-sqrmod_bnm1 t-mulmid	\
    t-hgcd t-hgcd_appr t-matrix22 t-invert t-bdiv t-fib2m			\
!   t-broot t-brootinv t-minvert t-sizeinbase
  
  EXTRA_DIST = toom-shared.h toom-sqr-shared.h
--- 30,34 ----
    t-div t-mul t-mullo t-sqrlo t-mulmod_bnm1 t-sqrmod_bnm1 t-mulmid	\
    t-hgcd t-hgcd_appr t-matrix22 t-invert t-bdiv t-fib2m			\
!   t-broot t-brootinv t-minvert t-sizeinbase t-gcd_11
  
  EXTRA_DIST = toom-shared.h toom-sqr-shared.h
diff -Nprc2 gmp-gcd_11.fc77a3dc184b/tests/mpn/t-gcd_11.c gmp-gcd_11/tests/mpn/t-gcd_11.c
*** gmp-gcd_11.fc77a3dc184b/tests/mpn/t-gcd_11.c	1970-01-01 01:00:00.000000000 +0100
--- gmp-gcd_11/tests/mpn/t-gcd_11.c	2019-08-07 23:23:35.549830000 +0200
***************
*** 0 ****
--- 1,82 ----
+ /* Test mpn_gcd_11.
+ 
+ Copyright 2019 Free Software Foundation, Inc.
+ 
+ This file is part of the GNU MP Library test suite.
+ 
+ The GNU MP Library test suite is free software; you can redistribute it
+ and/or modify it under the terms of the GNU General Public License as
+ published by the Free Software Foundation; either version 3 of the License,
+ or (at your option) any later version.
+ 
+ The GNU MP Library test suite is distributed in the hope that it will be
+ useful, but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General
+ Public License for more details.
+ 
+ You should have received a copy of the GNU General Public License along with
+ the GNU MP Library test suite.  If not, see https://www.gnu.org/licenses/.  */
+ 
+ #include <stdio.h>
+ #include <stdlib.h>
+ 
+ #include "gmp-impl.h"
+ #include "tests.h"
+ 
+ #ifndef COUNT
+ #define COUNT 10000
+ #endif
+ 
+ static void
+ one_test (mp_limb_t a, mp_limb_t b, mp_limb_t ref)
+ {
+   mp_limb_t r = mpn_gcd_11 (a, b);
+   if (r != ref)
+     {
+       gmp_fprintf (stderr,
+ 		   "gcd_11 (0x%Mx, 0x%Mx) failed, got: 0x%Mx, ref: 0x%Mx\n",
+ 		   a, b, r, ref);
+       abort();
+     }
+ }
+ 
+ int
+ main (int argc, char **argv)
+ {
+   mpz_t a, b;
+   int count = COUNT;
+   int test;
+   gmp_randstate_ptr rands;
+ 
+   TESTS_REPS (count, argv, argc);
+ 
+   tests_start ();
+   rands = RANDS;
+ 
+   mpz_init (a);
+   mpz_init (b);
+   for (test = 0; test < count; test++)
+     {
+       mp_limb_t al, bl;
+       mp_bitcnt_t size = 1 + gmp_urandomm_ui(rands, GMP_NUMB_BITS);
+       if (test & 1)
+ 	{
+ 	  mpz_urandomb (a, rands, size);
+ 	  mpz_urandomb (b, rands, size);
+ 	}
+       else
+ 	{
+ 	  mpz_rrandomb (a, rands, size);
+ 	  mpz_rrandomb (b, rands, size);
+ 	}
+ 
+       mpz_setbit (a, 0);
+       mpz_setbit (b, 0);
+       al = mpz_getlimbn (a, 0);
+       bl = mpz_getlimbn (b, 0);
+       one_test (al, bl, refmpn_gcd_11 (al, bl));
+     }
+ 
+   mpz_clear (a);
+   mpz_clear (b);
+ }
diff -Nprc2 gmp-gcd_11.fc77a3dc184b/tests/refmpn.c gmp-gcd_11/tests/refmpn.c
*** gmp-gcd_11.fc77a3dc184b/tests/refmpn.c	2019-08-07 17:21:54.000000000 +0200
--- gmp-gcd_11/tests/refmpn.c	2019-08-07 23:23:35.549830000 +0200
*************** refmpn_mul_any (mp_ptr prodp,
*** 1981,1984 ****
--- 1981,2002 ----
  
  mp_limb_t
+ refmpn_gcd_11 (mp_limb_t x, mp_limb_t y)
+ {
+   do
+     {
+       while ((x & 1) == 0)  x >>= 1;
+       while ((y & 1) == 0)  y >>= 1;
+ 
+       if (x < y)
+ 	MP_LIMB_T_SWAP (x, y);
+ 
+       x -= y;
+     }
+   while (x != 0);
+ 
+   return y;
+ }
+ 
+ mp_limb_t
  refmpn_gcd_1 (mp_srcptr xp, mp_size_t xsize, mp_limb_t y)
  {
*************** refmpn_gcd_1 (mp_srcptr xp, mp_size_t xs
*** 2003,2020 ****
      }
  
!   for (;;)
!     {
!       while ((x & 1) == 0)  x >>= 1;
!       while ((y & 1) == 0)  y >>= 1;
! 
!       if (x < y)
! 	MP_LIMB_T_SWAP (x, y);
! 
!       x -= y;
!       if (x == 0)
! 	break;
!     }
! 
!   return y << twos;
  }
  
--- 2021,2025 ----
      }
  
!   return refmpn_gcd_11 (x, y) << twos;
  }
  
diff -Nprc2 gmp-gcd_11.fc77a3dc184b/tests/tests.h gmp-gcd_11/tests/tests.h
*** gmp-gcd_11.fc77a3dc184b/tests/tests.h	2019-08-07 17:21:54.000000000 +0200
--- gmp-gcd_11/tests/tests.h	2019-08-07 23:23:35.549830000 +0200
*************** int refmpn_equal_anynail (mp_srcptr, mp_
*** 231,234 ****
--- 231,235 ----
  void refmpn_fill (mp_ptr, mp_size_t, mp_limb_t);
  
+ mp_limb_t refmpn_gcd_11 (mp_limb_t, mp_limb_t);
  mp_limb_t refmpn_gcd_1 (mp_srcptr, mp_size_t, mp_limb_t);
  mp_limb_t refmpn_gcd (mp_ptr, mp_ptr, mp_size_t, mp_ptr, mp_size_t);

-- 
Niels Möller. PGP-encrypted email is preferred. Keyid 368C6677.
Internet email is subject to wholesale government surveillance.