[Gmp-commit] /var/hg/gmp: 2 new changesets

mercurial at gmplib.org mercurial at gmplib.org
Mon Jun 15 19:57:47 UTC 2015


details:   /var/hg/gmp/rev/87c0b11998c5
changeset: 16715:87c0b11998c5
user:      Marco Bodrato <bodrato at mail.dm.unipi.it>
date:      Mon Jun 15 21:53:23 2015 +0200
description:
tests/devel/try.c: Support mpn_sqrt (sqrtrem with remainder = NULL)

details:   /var/hg/gmp/rev/92af7d28c9af
changeset: 16716:92af7d28c9af
user:      Marco Bodrato <bodrato at mail.dm.unipi.it>
date:      Mon Jun 15 21:57:38 2015 +0200
description:
ChangeLog

diffstat:

 ChangeLog             |   4 +++
 mpn/generic/sqrtrem.c |  25 ++++++++++---------
 tests/devel/try.c     |  62 +++++++++++++++++++++++++++++++++++++++++++++++++-
 tests/mpz/t-sqrtrem.c |   7 +++++
 4 files changed, 84 insertions(+), 14 deletions(-)

diffs (186 lines):

diff -r 5423329f78dd -r 92af7d28c9af ChangeLog
--- a/ChangeLog	Mon Jun 15 18:05:00 2015 +0200
+++ b/ChangeLog	Mon Jun 15 21:57:38 2015 +0200
@@ -1,3 +1,7 @@
+2015-06-15 Marco Bodrato <bodrato at mail.dm.unipi.it>
+
+	* tests/devel/try.c: Support mpn_sqrt (sqrtrem with remainder = NULL).
+
 2015-06-15  Torbjörn Granlund  <torbjorng at google.com>
 
 	* config.guess: Rewrite code for AVX handling to deal with broken cpuid
diff -r 5423329f78dd -r 92af7d28c9af mpn/generic/sqrtrem.c
--- a/mpn/generic/sqrtrem.c	Mon Jun 15 18:05:00 2015 +0200
+++ b/mpn/generic/sqrtrem.c	Mon Jun 15 21:57:38 2015 +0200
@@ -292,30 +292,31 @@
   ASSERT (! MPN_OVERLAP_P (sp, (nn + 1) / 2, np, nn));
 
   high = np[nn - 1];
-  if (nn == 1)
-    if (high & (GMP_NUMB_HIGHBIT | (GMP_NUMB_HIGHBIT >> 1)))
+  if (high & (GMP_NUMB_HIGHBIT | (GMP_NUMB_HIGHBIT / 2)))
+    c = 0;
+  else
+    {
+      count_leading_zeros (c, high);
+      c -= GMP_NAIL_BITS;
+
+      c = c / 2; /* we have to shift left by 2c bits to normalize {np, nn} */
+    }
+  if (nn == 1) {
+    if (c == 0)
       {
 	sp[0] = mpn_sqrtrem1 (&rl, high);
 	if (rp != NULL)
 	  rp[0] = rl;
-	return rl != 0;
       }
     else
       {
-	count_leading_zeros (c, high);
-	c -= GMP_NAIL_BITS;
-
-	c = c / 2; /* we have to shift left by 2c bits to normalize {np, nn} */
 	cc = mpn_sqrtrem1 (&rl, high << (2*c)) >> c;
 	sp[0] = cc;
 	if (rp != NULL)
 	  rp[0] = rl = high - cc*cc;
-	return rl != 0;
       }
-  count_leading_zeros (c, high);
-  c -= GMP_NAIL_BITS;
-
-  c = c / 2; /* we have to shift left by 2c bits to normalize {np, nn} */
+    return rl != 0;
+  }
   tn = (nn + 1) / 2; /* 2*tn is the smallest even integer >= nn */
 
   TMP_MARK;
diff -r 5423329f78dd -r 92af7d28c9af tests/devel/try.c
--- a/tests/devel/try.c	Mon Jun 15 18:05:00 2015 +0200
+++ b/tests/devel/try.c	Mon Jun 15 21:57:38 2015 +0200
@@ -618,6 +618,46 @@
     validate_fail ();
 }
 
+void
+validate_sqrt (void)
+{
+  mp_srcptr  orig_ptr = s[0].p;
+  mp_size_t  orig_size = size;
+  mp_size_t  root_size = (size+1)/2;
+  mp_srcptr  root_ptr = fun.d[0].p;
+  int        perf_pow = (fun.retval == 0);
+  mp_size_t  prod_size = 2*root_size;
+  mp_ptr     p;
+  mp_limb_t  cy;
+  int  error = 0;
+
+  p = refmpn_malloc_limbs (prod_size);
+
+  refmpn_sqr (p, root_ptr, root_size);
+  MPN_NORMALIZE (p, prod_size);
+  if (refmpn_cmp_twosizes (p,prod_size, orig_ptr,orig_size) != - !perf_pow)
+    {
+      printf ("root^2 bigger than original.\n");
+      mpn_trace ("prod", p, prod_size);
+      error = 1;
+    }
+
+  refmpn_sub (p, orig_ptr,orig_size, p,prod_size);
+  MPN_NORMALIZE (p, prod_size);
+  if (prod_size >= root_size &&
+      refmpn_sub (p, p,prod_size, root_ptr, root_size) == 0 &&
+      refmpn_cmp_twosizes (p, prod_size, root_ptr, root_size) > 0)
+    {
+      printf ("(root+1)^2 smaller than original.\n");
+      mpn_trace ("prod", p, prod_size);
+      error = 1;
+    }
+  free (p);
+
+  if (error)
+    validate_fail ();
+}
+
 
 /* These types are indexes into the param[] array and are arbitrary so long
    as they're all distinct and within the size of param[].  Renumber
@@ -680,7 +720,7 @@
 
   TYPE_SBPI1_DIV_QR, TYPE_TDIV_QR,
 
-  TYPE_SQRTREM, TYPE_ZERO, TYPE_GET_STR, TYPE_POPCOUNT, TYPE_HAMDIST,
+  TYPE_SQRTREM, TYPE_SQRT, TYPE_ZERO, TYPE_GET_STR, TYPE_POPCOUNT, TYPE_HAMDIST,
 
   TYPE_EXTRA
 };
@@ -1402,6 +1442,15 @@
   VALIDATE (validate_sqrtrem);
   REFERENCE (refmpn_sqrtrem);
 
+  p = &param[TYPE_SQRT];
+  p->retval = 1;
+  p->dst[0] = 1;
+  p->dst[1] = 0;
+  p->src[0] = 1;
+  p->dst_size[0] = SIZE_CEIL_HALF;
+  p->overlap = OVERLAP_NONE;
+  VALIDATE (validate_sqrt);
+
   p = &param[TYPE_ZERO];
   p->dst[0] = 1;
   p->size = SIZE_ALLOW_ZERO;
@@ -1669,6 +1718,9 @@
 MPN_ZERO_fun (mp_ptr ptr, mp_size_t size)
 { MPN_ZERO (ptr, size); }
 
+mp_size_t
+mpn_sqrt_fun (mp_ptr dst, mp_srcptr src, mp_size_t size)
+{ return mpn_sqrtrem (dst, NULL, src, size); }
 
 struct choice_t {
   const char  *name;
@@ -1970,7 +2022,8 @@
   { TRY(mpn_popcount),   TYPE_POPCOUNT },
   { TRY(mpn_hamdist),    TYPE_HAMDIST },
 
-  { TRY(mpn_sqrtrem),    TYPE_SQRTREM },
+  { TRY(mpn_sqrtrem),     TYPE_SQRTREM },
+  { TRY_FUNFUN(mpn_sqrt), TYPE_SQRT },
 
   { TRY_FUNFUN(MPN_ZERO), TYPE_ZERO },
 
@@ -2778,6 +2831,11 @@
       (e->d[0].p, e->d[1].p, e->s[0].p, size);
     break;
 
+  case TYPE_SQRT:
+    e->retval = (* (long (*)(ANYARGS)) CALLING_CONVENTIONS (function))
+      (e->d[0].p, e->s[0].p, size);
+    break;
+
   case TYPE_ZERO:
     CALLING_CONVENTIONS (function) (e->d[0].p, size);
     break;
diff -r 5423329f78dd -r 92af7d28c9af tests/mpz/t-sqrtrem.c
--- a/tests/mpz/t-sqrtrem.c	Mon Jun 15 18:05:00 2015 +0200
+++ b/tests/mpz/t-sqrtrem.c	Mon Jun 15 21:57:38 2015 +0200
@@ -64,10 +64,17 @@
 
       /* printf ("%ld\n", SIZ (x2)); */
 
+      mpz_sqrt (temp, x2);
+      MPZ_CHECK_FORMAT (temp);
+
       mpz_sqrtrem (x, rem, x2);
       MPZ_CHECK_FORMAT (x);
       MPZ_CHECK_FORMAT (rem);
 
+      /* Are results different?  */
+      if (mpz_cmp (temp, x) != 0)
+	dump_abort (x2, x, rem);
+
       mpz_mul (temp, x, x);
 
       /* Is square of result > argument?  */


More information about the gmp-commit mailing list