math: fix remquo.c when x==-y and a subnormal remainder bug as well
authornsz <nsz@port70.net>
Mon, 7 May 2012 22:22:56 +0000 (00:22 +0200)
committernsz <nsz@port70.net>
Mon, 7 May 2012 22:22:56 +0000 (00:22 +0200)
backported fix from freebsd:
http://svnweb.FreeBSD.org/base?view=revision&revision=233973

src/math/remquo.c
src/math/remquof.c
src/math/remquol.c

index e92984edf2ec018bdc801738ea9e578047d4ff7c..085466e6a2f84ba16bdabfe8eb6d94bfa7e8a7b1 100644 (file)
@@ -44,7 +44,7 @@ double remquo(double x, double y, int *quo)
                        goto fixup;
                }
                if (lx == ly) {            /* |x| = |y| return x*0 */
-                       *quo = 1;
+                       *quo = sxy ? -1 : 1;
                        return Zero[(uint32_t)sx>>31];
                }
        }
@@ -127,6 +127,7 @@ double remquo(double x, double y, int *quo)
 
        /* convert back to floating value and restore the sign */
        if ((hx|lx) == 0) {  /* return sign(x)*0 */
+               q &= 0x7fffffff;
                *quo = sxy ? -q : q;
                return Zero[(uint32_t)sx>>31];
        }
@@ -144,10 +145,10 @@ double remquo(double x, double y, int *quo)
                        hx >>= n;
                } else if (n <= 31) {
                        lx = (hx<<(32-n))|(lx>>n);
-                       hx = sx;
+                       hx = 0;
                } else {
                        lx = hx>>(n-32);
-                       hx = sx;
+                       hx = 0;
                }
        }
 fixup:
index 11569ce8a386a892913c1121bc66585f1fd9434a..536a050aa596558158eb825a8dd58f799412beb2 100644 (file)
@@ -41,7 +41,7 @@ float remquof(float x, float y, int *quo)
                q = 0;
                goto fixup;
        } else if(hx==hy) {  /* |x| = |y| return x*0*/
-               *quo = 1;
+               *quo = sxy ? -1 : 1;
                return Zero[(uint32_t)sx>>31];
        }
 
@@ -92,6 +92,7 @@ float remquof(float x, float y, int *quo)
 
        /* convert back to floating value and restore the sign */
        if (hx == 0) {                             /* return sign(x)*0 */
+               q &= 0x7fffffff;
                *quo = sxy ? -q : q;
                return Zero[(uint32_t)sx>>31];
        }
index 721231b4c02184e98dcb9ed21f3d4eb781995651..a2e11728169d3c74dcd3676f3919b62c10cd939a 100644 (file)
@@ -94,7 +94,7 @@ long double remquol(long double x, long double y, int *quo)
                        goto fixup;       /* |x|<|y| return x or x-y */
                }
                if (ux.bits.manh == uy.bits.manh && ux.bits.manl == uy.bits.manl) {
-                       *quo = 1;
+                       *quo = sxy ? -1 : 1;
                        return Zero[sx];  /* |x|=|y| return x*0*/
                }
        }
@@ -152,6 +152,7 @@ long double remquol(long double x, long double y, int *quo)
 
        /* convert back to floating value and restore the sign */
        if ((hx|lx) == 0) {  /* return sign(x)*0 */
+               q &= 0x7fffffff;
                *quo = sxy ? -q : q;
                return Zero[sx];
        }