Add libgmp 6.2.0 to third_party
Don't build it yet. That will come in the next review.
Change-Id: Idf3266558165e5ab45f4a41c98cc8c838c8244d5
diff --git a/third_party/gmp/tests/refmpz.c b/third_party/gmp/tests/refmpz.c
new file mode 100644
index 0000000..167799f
--- /dev/null
+++ b/third_party/gmp/tests/refmpz.c
@@ -0,0 +1,343 @@
+/* Reference mpz functions.
+
+Copyright 1997, 1999-2002 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/. */
+
+/* always do assertion checking */
+#define WANT_ASSERT 1
+
+#include <stdio.h>
+#include <stdlib.h> /* for free */
+#include "gmp-impl.h"
+#include "longlong.h"
+#include "tests.h"
+
+
+/* Change this to "#define TRACE(x) x" for some traces. */
+#define TRACE(x)
+
+
+/* FIXME: Shouldn't use plain mpz functions in a reference routine. */
+void
+refmpz_combit (mpz_ptr r, unsigned long bit)
+{
+ if (mpz_tstbit (r, bit))
+ mpz_clrbit (r, bit);
+ else
+ mpz_setbit (r, bit);
+}
+
+
+unsigned long
+refmpz_hamdist (mpz_srcptr x, mpz_srcptr y)
+{
+ mp_size_t xsize, ysize, tsize;
+ mp_ptr xp, yp;
+ unsigned long ret;
+
+ if ((SIZ(x) < 0 && SIZ(y) >= 0)
+ || (SIZ(y) < 0 && SIZ(x) >= 0))
+ return ULONG_MAX;
+
+ xsize = ABSIZ(x);
+ ysize = ABSIZ(y);
+ tsize = MAX (xsize, ysize);
+
+ xp = refmpn_malloc_limbs (tsize);
+ refmpn_zero (xp, tsize);
+ refmpn_copy (xp, PTR(x), xsize);
+
+ yp = refmpn_malloc_limbs (tsize);
+ refmpn_zero (yp, tsize);
+ refmpn_copy (yp, PTR(y), ysize);
+
+ if (SIZ(x) < 0)
+ refmpn_neg (xp, xp, tsize);
+
+ if (SIZ(x) < 0)
+ refmpn_neg (yp, yp, tsize);
+
+ ret = refmpn_hamdist (xp, yp, tsize);
+
+ free (xp);
+ free (yp);
+ return ret;
+}
+
+void
+refmpz_gcd (mpz_ptr g, mpz_srcptr a_orig, mpz_srcptr b_orig)
+{
+ mp_bitcnt_t a_twos, b_twos, common_twos;
+ mpz_t a;
+ mpz_t b;
+ mpz_init (a);
+ mpz_init (b);
+ mpz_abs (a, a_orig);
+ mpz_abs (b, b_orig);
+
+ if (mpz_sgn (a) == 0)
+ {
+ mpz_set (g, b);
+ return;
+ }
+ if (mpz_sgn (b) == 0)
+ {
+ mpz_set (g, a);
+ return;
+ }
+ a_twos = mpz_scan1 (a, 0);
+ mpz_tdiv_q_2exp (a, a, a_twos);
+
+ b_twos = mpz_scan1 (b, 0);
+ mpz_tdiv_q_2exp (b, b, b_twos);
+
+ common_twos = MIN(a_twos, b_twos);
+ for (;;)
+ {
+ int c;
+ mp_bitcnt_t twos;
+ c = mpz_cmp (a, b);
+ if (c == 0)
+ break;
+ if (c < 0)
+ mpz_swap (a, b);
+ mpz_sub (a, a, b);
+ twos = mpz_scan1 (a, 0);
+ mpz_tdiv_q_2exp (a, a, twos);
+ }
+ mpz_mul_2exp (g, a, common_twos);
+
+ mpz_clear (a);
+ mpz_clear (b);
+}
+
+/* (0/b), with mpz b; is 1 if b=+/-1, 0 otherwise */
+#define JACOBI_0Z(b) JACOBI_0LS (PTR(b)[0], SIZ(b))
+
+/* (a/b) effect due to sign of b: mpz/mpz */
+#define JACOBI_BSGN_ZZ_BIT1(a, b) JACOBI_BSGN_SS_BIT1 (SIZ(a), SIZ(b))
+
+/* (a/b) effect due to sign of a: mpz/unsigned-mpz, b odd;
+ is (-1/b) if a<0, or +1 if a>=0 */
+#define JACOBI_ASGN_ZZU_BIT1(a, b) JACOBI_ASGN_SU_BIT1 (SIZ(a), PTR(b)[0])
+
+int
+refmpz_kronecker (mpz_srcptr a_orig, mpz_srcptr b_orig)
+{
+ unsigned long twos;
+ mpz_t a, b;
+ int result_bit1 = 0;
+
+ if (mpz_sgn (b_orig) == 0)
+ return JACOBI_Z0 (a_orig); /* (a/0) */
+
+ if (mpz_sgn (a_orig) == 0)
+ return JACOBI_0Z (b_orig); /* (0/b) */
+
+ if (mpz_even_p (a_orig) && mpz_even_p (b_orig))
+ return 0;
+
+ if (mpz_cmp_ui (b_orig, 1) == 0)
+ return 1;
+
+ mpz_init_set (a, a_orig);
+ mpz_init_set (b, b_orig);
+
+ if (mpz_sgn (b) < 0)
+ {
+ result_bit1 ^= JACOBI_BSGN_ZZ_BIT1 (a, b);
+ mpz_neg (b, b);
+ }
+ if (mpz_even_p (b))
+ {
+ twos = mpz_scan1 (b, 0L);
+ mpz_tdiv_q_2exp (b, b, twos);
+ result_bit1 ^= JACOBI_TWOS_U_BIT1 (twos, PTR(a)[0]);
+ }
+
+ if (mpz_sgn (a) < 0)
+ {
+ result_bit1 ^= JACOBI_N1B_BIT1 (PTR(b)[0]);
+ mpz_neg (a, a);
+ }
+ if (mpz_even_p (a))
+ {
+ twos = mpz_scan1 (a, 0L);
+ mpz_tdiv_q_2exp (a, a, twos);
+ result_bit1 ^= JACOBI_TWOS_U_BIT1 (twos, PTR(b)[0]);
+ }
+
+ for (;;)
+ {
+ ASSERT (mpz_odd_p (a));
+ ASSERT (mpz_odd_p (b));
+ ASSERT (mpz_sgn (a) > 0);
+ ASSERT (mpz_sgn (b) > 0);
+
+ TRACE (printf ("top\n");
+ mpz_trace (" a", a);
+ mpz_trace (" b", b));
+
+ if (mpz_cmp (a, b) < 0)
+ {
+ TRACE (printf ("swap\n"));
+ mpz_swap (a, b);
+ result_bit1 ^= JACOBI_RECIP_UU_BIT1 (PTR(a)[0], PTR(b)[0]);
+ }
+
+ if (mpz_cmp_ui (b, 1) == 0)
+ break;
+
+ mpz_sub (a, a, b);
+ TRACE (printf ("sub\n");
+ mpz_trace (" a", a));
+ if (mpz_sgn (a) == 0)
+ goto zero;
+
+ twos = mpz_scan1 (a, 0L);
+ mpz_fdiv_q_2exp (a, a, twos);
+ TRACE (printf ("twos %lu\n", twos);
+ mpz_trace (" a", a));
+ result_bit1 ^= JACOBI_TWOS_U_BIT1 (twos, PTR(b)[0]);
+ }
+
+ mpz_clear (a);
+ mpz_clear (b);
+ return JACOBI_BIT1_TO_PN (result_bit1);
+
+ zero:
+ mpz_clear (a);
+ mpz_clear (b);
+ return 0;
+}
+
+/* Same as mpz_kronecker, but ignoring factors of 2 on b */
+int
+refmpz_jacobi (mpz_srcptr a, mpz_srcptr b)
+{
+ ASSERT_ALWAYS (mpz_sgn (b) > 0);
+ ASSERT_ALWAYS (mpz_odd_p (b));
+
+ return refmpz_kronecker (a, b);
+}
+
+/* Legendre symbol via powm. p must be an odd prime. */
+int
+refmpz_legendre (mpz_srcptr a, mpz_srcptr p)
+{
+ int res;
+
+ mpz_t r;
+ mpz_t e;
+
+ ASSERT_ALWAYS (mpz_sgn (p) > 0);
+ ASSERT_ALWAYS (mpz_odd_p (p));
+
+ mpz_init (r);
+ mpz_init (e);
+
+ mpz_fdiv_r (r, a, p);
+
+ mpz_set (e, p);
+ mpz_sub_ui (e, e, 1);
+ mpz_fdiv_q_2exp (e, e, 1);
+ mpz_powm (r, r, e, p);
+
+ /* Normalize to a more or less symmetric range around zero */
+ if (mpz_cmp (r, e) > 0)
+ mpz_sub (r, r, p);
+
+ ASSERT_ALWAYS (mpz_cmpabs_ui (r, 1) <= 0);
+
+ res = mpz_sgn (r);
+
+ mpz_clear (r);
+ mpz_clear (e);
+
+ return res;
+}
+
+
+int
+refmpz_kronecker_ui (mpz_srcptr a, unsigned long b)
+{
+ mpz_t bz;
+ int ret;
+ mpz_init_set_ui (bz, b);
+ ret = refmpz_kronecker (a, bz);
+ mpz_clear (bz);
+ return ret;
+}
+
+int
+refmpz_kronecker_si (mpz_srcptr a, long b)
+{
+ mpz_t bz;
+ int ret;
+ mpz_init_set_si (bz, b);
+ ret = refmpz_kronecker (a, bz);
+ mpz_clear (bz);
+ return ret;
+}
+
+int
+refmpz_ui_kronecker (unsigned long a, mpz_srcptr b)
+{
+ mpz_t az;
+ int ret;
+ mpz_init_set_ui (az, a);
+ ret = refmpz_kronecker (az, b);
+ mpz_clear (az);
+ return ret;
+}
+
+int
+refmpz_si_kronecker (long a, mpz_srcptr b)
+{
+ mpz_t az;
+ int ret;
+ mpz_init_set_si (az, a);
+ ret = refmpz_kronecker (az, b);
+ mpz_clear (az);
+ return ret;
+}
+
+
+void
+refmpz_pow_ui (mpz_ptr w, mpz_srcptr b, unsigned long e)
+{
+ mpz_t s, t;
+ unsigned long i;
+
+ mpz_init_set_ui (t, 1L);
+ mpz_init_set (s, b);
+
+ if ((e & 1) != 0)
+ mpz_mul (t, t, s);
+
+ for (i = 2; i <= e; i <<= 1)
+ {
+ mpz_mul (s, s, s);
+ if ((i & e) != 0)
+ mpz_mul (t, t, s);
+ }
+
+ mpz_set (w, t);
+
+ mpz_clear (s);
+ mpz_clear (t);
+}