Re: [bug-gawk] GAWK 4.1 SQRT() bug

[Aharon Robbins]

Greetings. Re this:

> Date: Sat, 02 Aug 2014 12:27:00 -0400
> To: address@hidden
> From: Katherine Wasserman <address@hidden>
> Subject: [bug-gawk] GAWK 4.1 SQRT() bug
>
> In version 4.1.60 of GAWK the sqrt() function does not work correctly on 
> bignums.

Turns out to be true on the stable version also.

> Here's a demo of the problem along with, a function that does work correctly.
>
> Running this program (sqrt-bug.awk):
> --------------------------------------------------------------------
> BEGIN {
> a=11111111111111111111111111111111111111111111111111111111111
> print sqrt(a^2)
> print sq_root(a^2)
> }
>
>
> function sq_root(x, temp,r,z)
> {  temp=substr(x,1,length(x)/2) + 0 # a good first guess
>    z=0
>    while (abs(z-temp)>1)
>     { z=temp
>       div(x,temp,r)
>       temp=r["quotient"] + temp
>       div(temp,2,r)
>       temp=r["quotient"]
>     }
>    return temp
> }
>
> function abs(x)
> { return (x<0 ? -x : x)
> }
>
> --------------------------------------------------------------------
> gawk -M -f sqrt-bug.awk
>
> results in:
> 11111111111111111261130863809439559987542611609749437808640
> 11111111111111111111111111111111111111111111111111111111111
>
> Thanks,
> Katie

This is indeed a bug. At least for sqrt(), the precision of the
result needs to match that of the original argument. Doing this
for all functions doesn't seem to hurt.

Here is the patch that I will be pushing shortly. I will be adding
your program to the test suite, also.

Thanks!

Arnold
------------------------------------------------------------
diff --git a/mpfr.c b/mpfr.c
index cdd64b5..66e2b04 100644
--- a/mpfr.c
+++ b/mpfr.c
@@ -697,22 +697,35 @@ do_mpfr_atan2(int nargs)
        return res;
 }
 
+/* do_mpfr_func --- run an MPFR function - not inline, for debugging */
 
-#define SPEC_MATH(X)                                           \
-NODE *t1, *res;                                                        \
-mpfr_ptr p1;                                                   \
-int tval;                                                      \
-t1 = POP_SCALAR();                                             \
-if (do_lint && (t1->flags & (NUMCUR|NUMBER)) == 0)             \
-       lintwarn(_("%s: received non-numeric argument"), #X);   \
-force_number(t1);                                              \
-p1 = MP_FLOAT(t1);                                             \
-res = mpg_float();                                             \
-tval = mpfr_##X(res->mpg_numbr, p1, ROUND_MODE);                       \
-IEEE_FMT(res->mpg_numbr, tval);                                        \
-DEREF(t1);                                                     \
-return res
+static NODE *
+do_mpfr_func(const char *name,
+               int (*mpfr_func)(),     /* putting argument types just gets the 
compiler confused */
+               int nargs)
+{
+       NODE *t1, *res;
+       mpfr_ptr p1;
+       int tval;
 
+       t1 = POP_SCALAR();
+       if (do_lint && (t1->flags & (NUMCUR|NUMBER)) == 0)
+               lintwarn(_("%s: received non-numeric argument"), name);
+
+       force_number(t1);
+       p1 = MP_FLOAT(t1);
+       res = mpg_float();
+       mpfr_set_prec(res->mpg_numbr, mpfr_get_prec(p1));       /* needed at 
least for sqrt() */
+       tval = mpfr_func(res->mpg_numbr, p1, ROUND_MODE);
+       IEEE_FMT(res->mpg_numbr, tval);
+       DEREF(t1);
+       return res;
+}
+
+#define SPEC_MATH(X)                           \
+NODE *result;                                  \
+result = do_mpfr_func(#X, mpfr_##X, nargs);    \
+return result
 
 /* do_mpfr_sin --- do the sin function */