Re: [GMP kernel] using mpn_sqrtrem for sqrtr

[Bill Allombert <allomber@math.u-bordeaux.fr>]

On Mon, Sep 20, 2004 at 10:52:02PM +0200, Bill Allombert wrote:
> Hello PARI-dev,
> 
> Here a patch that implement sqrtr using the GMP function
> mpn_sqrtrem() (when the GMP kernel is used).

Here a patch for the native kernel that implement the same
algorithm as sqrtr_with_gmp (in the function sqrtr_without_gmp)
With that patch sqrtr_with_gmp and sqrtr_without_gmp should
produce the exact same result.

The problem is that sqrtr_without_gmp is 2.5 times slower than
sqrtr_abs(), which probably means that isqrti is 2.5 times
slower than sqrtr_abs() with the native kernel, which is a bug.

Cheers,
Bill.

Index: src/kernel/none/mp.c
===================================================================
RCS file: /home/cvs/pari/src/kernel/none/mp.c,v
retrieving revision 1.138
diff -u -r1.138 mp.c
--- src/kernel/none/mp.c	14 Sep 2004 20:41:22 -0000	1.138
+++ src/kernel/none/mp.c	21 Sep 2004 10:42:37 -0000
@@ -1505,6 +1505,28 @@
 
 /* Return trunc(sqrt(a))). a must be an non-negative integer*/
 GEN isqrti(GEN a) {return racine_r(a,lgefint(a));}
+
+GEN sqrtr_without_gmp(GEN a)
+{
+  GEN res, b, c;
+  long pr=lg(a)-2;
+  long e=expo(a),er=e>>1;
+  int i;
+  if (!signe(a)) return realzero_bit(er);
+  res = cgetr(2 + pr);
+  b = new_chunk(pr + 2);
+  for (i=0; i<pr+2; i++) b[i]=0;
+  c = ishiftr_spec(a,pr+2,(e&1)?0:-1);
+  setsigne(c,1);
+  c = racine_r(c,(pr<<1)+4);
+  i=pr+2;
+  if (c[i]<0)
+    while( !++c[--i]);
+  for(i=2;i<pr+2;i++) res[i]=c[i];
+  res[1] = evalsigne(1) | evalexpo(er);
+  avma=(pari_sp)res;
+  return res;
+}
 
 /* Normalize a non-negative integer.  */
 GEN