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

[mercurial at gmplib.org]

details:   /var/hg/gmp/rev/be4dbe9f533d
changeset: 16718:be4dbe9f533d
user:      Marco Bodrato <bodrato at mail.dm.unipi.it>
date:      Mon Jun 22 00:24:18 2015 +0200
description:
mpn/generic/sqrtrem.c (mpn_sqrtrem2): Remove branches, simplify

details:   /var/hg/gmp/rev/f953f9982d6e
changeset: 16719:f953f9982d6e
user:      Marco Bodrato <bodrato at mail.dm.unipi.it>
date:      Mon Jun 22 00:31:05 2015 +0200
description:
mpn/generic/sqrtrem.c (mpn_dc_sqrtrem): Move an add into a branch.

diffstat:

 mpn/generic/sqrtrem.c |  46 +++++++++++++++++-----------------------------
 1 files changed, 17 insertions(+), 29 deletions(-)

diffs (77 lines):

diff -r ac54869212ab -r f953f9982d6e mpn/generic/sqrtrem.c
--- a/mpn/generic/sqrtrem.c	Mon Jun 15 22:56:58 2015 +0200
+++ b/mpn/generic/sqrtrem.c	Mon Jun 22 00:31:05 2015 +0200
@@ -171,45 +171,33 @@
 static mp_limb_t
 mpn_sqrtrem2 (mp_ptr sp, mp_ptr rp, mp_srcptr np)
 {
-  mp_limb_t qhl, q, u, np0, sp0, rp0, q2;
+  mp_limb_t q, u, np0, sp0, rp0, q2;
   int cc;
 
   ASSERT (np[1] >= GMP_NUMB_HIGHBIT / 2);
 
   np0 = np[0];
   sp0 = mpn_sqrtrem1 (rp, np[1]);
-  qhl = 0;
   rp0 = rp[0];
-  while (rp0 >= sp0)
-    {
-      qhl++;
-      rp0 -= sp0;
-    }
-  /* now rp0 < sp0 < 2^Prec */
-  rp0 = (rp0 << Prec) + (np0 >> Prec);
-  u = 2 * sp0;
-  q = rp0 / u;
-  u = rp0 - q * u;
-  q += (qhl & 1) << (Prec - 1);
-  qhl >>= 1; /* if qhl=1, necessary q=0 as qhl*2^Prec + q <= 2^Prec */
-  /* now we have (initial rp0)<<Prec + np0>>Prec = (qhl<<Prec + q) * (2sp0) + u */
-  sp0 = ((sp0 + qhl) << Prec) + q;
-  cc = u >> Prec;
-  rp0 = ((u << Prec) & GMP_NUMB_MASK) + (np0 & ((CNST_LIMB (1) << Prec) - 1));
-  /* subtract q * q or qhl*2^(2*Prec) from rp */
+  /* rp0 <= 2*sp0 < 2^(Prec + 1) */
+  rp0 = (rp0 << (Prec - 1)) + (np0 >> (Prec + 1));
+  q = rp0 / sp0;
+  /* q <= 2^Prec, if q = 2^Prec, reduce the overestimate. */
+  q -= q >> Prec;
+  /* now we have q < 2^Prec */
+  u = rp0 - q * sp0;
+  /* now we have (rp[0]<<Prec + np0>>Prec)/2 = q * sp0 + u */
+  sp0 = (sp0 << Prec) | q;
+  cc = u >> (Prec - 1);
+  rp0 = ((u << (Prec + 1)) & GMP_NUMB_MASK) + (np0 & ((CNST_LIMB (1) << (Prec + 1)) - 1));
+  /* subtract q * q from rp */
   q2 = q * q;
-  cc -= (rp0 < q2) + qhl;
+  cc -= rp0 < q2;
   rp0 -= q2;
-  /* now subtract 2*q*2^Prec + 2^(2*Prec) if qhl is set */
   if (cc < 0)
     {
-      if (sp0 != 0)
-	{
-	  rp0 += sp0;
-	  cc += rp0 < sp0;
-	}
-      else
-	cc++;
+      rp0 += sp0;
+      cc += rp0 < sp0;
       --sp0;
       rp0 += sp0;
       cc += rp0 < sp0;
@@ -253,10 +241,10 @@
       mpn_sqr (np + n, sp, l);
       b = q + mpn_sub_n (np, np, np + n, 2 * l);
       c -= (l == h) ? b : mpn_sub_1 (np + 2 * l, np + 2 * l, 1, (mp_limb_t) b);
-      q = mpn_add_1 (sp + l, sp + l, h, q);
 
       if (c < 0)
 	{
+	  q = mpn_add_1 (sp + l, sp + l, h, q);
 #if HAVE_NATIVE_mpn_addlsh1_n
 	  c += mpn_addlsh1_n (np, np, sp, n) + 2 * q;
 #else