This is the mail archive of the libc-hacker@sourceware.org mailing list for the glibc project.
Note that libc-hacker is a closed list. You may look at the archives of this list, but subscription and posting are not open.
| Index Nav: | [Date Index] [Subject Index] [Author Index] [Thread Index] | |
|---|---|---|
| Message Nav: | [Date Prev] [Date Next] | [Thread Prev] [Thread Next] |
| Other format: | [Raw text] | |
Hi!
http://www.opengroup.org/onlinepubs/009695399/functions/pow.html
says that pow{,f,l} should return NaN when either x or y is NaN,
except for the few exceptions (pow (1, NAN) and pow (NAN, +-0) should
return 1.0). But dbl-64 generic pow was checking for NAN too late
and missed a couple of cases. While adding tests I noticed the i386/x86_64
assembly pow* raise invalid exception when they should not, so that's
fixed in the patch as well.
2007-03-01 Jakub Jelinek <jakub@redhat.com>
[BZ #4096]
* sysdeps/ieee754/dbl-64/e_pow.c (__ieee754_pow): Check for NaN earlier.
* math/libm-test.inc (pow_test): Add more tests involving NaNs.
* sysdeps/i386/fpu/e_powf.S (__ieee754_powf): Avoid invalid exception
for x qNaN and y either +-inf or non-integer value.
* sysdeps/i386/fpu/e_pow.S (__ieee754_pow): Likewise.
* sysdeps/i386/fpu/e_powl.S (__ieee754_powl): Likewise.
* sysdeps/x86_64/fpu/e_powl.S (__ieee754_powl): Likewise.
--- libc/math/libm-test.inc.jj 2006-09-20 21:36:26.000000000 +0200
+++ libc/math/libm-test.inc 2007-03-01 15:06:05.000000000 +0100
@@ -4614,6 +4614,17 @@ pow_test (void)
/* pow (x, NaN) == NaN. */
TEST_ff_f (pow, 3.0, nan_value, nan_value);
+ TEST_ff_f (pow, minus_zero, nan_value, nan_value);
+ TEST_ff_f (pow, plus_infty, nan_value, nan_value);
+ TEST_ff_f (pow, -3.0, nan_value, nan_value);
+ TEST_ff_f (pow, minus_infty, nan_value, nan_value);
+
+ TEST_ff_f (pow, nan_value, 3.0, nan_value);
+ TEST_ff_f (pow, nan_value, -3.0, nan_value);
+ TEST_ff_f (pow, nan_value, plus_infty, nan_value);
+ TEST_ff_f (pow, nan_value, minus_infty, nan_value);
+ TEST_ff_f (pow, nan_value, 2.5, nan_value);
+ TEST_ff_f (pow, nan_value, -2.5, nan_value);
TEST_ff_f (pow, 1, plus_infty, 1);
TEST_ff_f (pow, -1, plus_infty, 1);
--- libc/sysdeps/i386/fpu/e_pow.S.jj 2005-05-04 19:45:15.000000000 +0200
+++ libc/sysdeps/i386/fpu/e_pow.S 2007-03-01 16:10:28.000000000 +0100
@@ -1,5 +1,5 @@
/* ix87 specific implementation of pow function.
- Copyright (C) 1996, 1997, 1998, 1999, 2001, 2004, 2005
+ Copyright (C) 1996, 1997, 1998, 1999, 2001, 2004, 2005, 2007
Free Software Foundation, Inc.
This file is part of the GNU C Library.
Contributed by Ulrich Drepper <drepper@cygnus.com>, 1996.
@@ -161,10 +161,11 @@ ENTRY(__ieee754_pow)
2: /* y is a real number. */
fxch // x : y
fldl MO(one) // 1.0 : x : y
- fld %st(1) // x : 1.0 : x : y
- fsub %st(1) // x-1 : 1.0 : x : y
- fabs // |x-1| : 1.0 : x : y
- fcompl MO(limit) // 1.0 : x : y
+ fldl MO(limit) // 0.29 : 1.0 : x : y
+ fld %st(2) // x : 0.29 : 1.0 : x : y
+ fsub %st(2) // x-1 : 0.29 : 1.0 : x : y
+ fabs // |x-1| : 0.29 : 1.0 : x : y
+ fucompp // 1.0 : x : y
fnstsw
fxch // x : 1.0 : y
sahf
@@ -197,9 +198,10 @@ ENTRY(__ieee754_pow)
// y == ±inf
.align ALIGNARG(4)
12: fstp %st(0) // pop y
- fldl 4(%esp) // x
- fabs
- fcompl MO(one) // < 1, == 1, or > 1
+ fldl MO(one) // 1
+ fldl 4(%esp) // x : 1
+ fabs // abs(x) : 1
+ fucompp // < 1, == 1, or > 1
fnstsw
andb $0x45, %ah
cmpb $0x45, %ah
--- libc/sysdeps/i386/fpu/e_powl.S.jj 2005-05-04 19:45:15.000000000 +0200
+++ libc/sysdeps/i386/fpu/e_powl.S 2007-03-01 16:02:21.000000000 +0100
@@ -1,5 +1,5 @@
/* ix87 specific implementation of pow function.
- Copyright (C) 1996, 1997, 1998, 1999, 2001, 2004, 2005
+ Copyright (C) 1996, 1997, 1998, 1999, 2001, 2004, 2005, 2007
Free Software Foundation, Inc.
This file is part of the GNU C Library.
Contributed by Ulrich Drepper <drepper@cygnus.com>, 1996.
@@ -161,10 +161,11 @@ ENTRY(__ieee754_powl)
2: /* y is a real number. */
fxch // x : y
fldl MO(one) // 1.0 : x : y
- fld %st(1) // x : 1.0 : x : y
- fsub %st(1) // x-1 : 1.0 : x : y
- fabs // |x-1| : 1.0 : x : y
- fcompl MO(limit) // 1.0 : x : y
+ fldl MO(limit) // 0.29 : 1.0 : x : y
+ fld %st(2) // x : 0.29 : 1.0 : x : y
+ fsub %st(2) // x-1 : 0.29 : 1.0 : x : y
+ fabs // |x-1| : 0.29 : 1.0 : x : y
+ fucompp // 1.0 : x : y
fnstsw
fxch // x : 1.0 : y
sahf
@@ -210,9 +211,10 @@ ENTRY(__ieee754_powl)
// y == ±inf
.align ALIGNARG(4)
12: fstp %st(0) // pop y
- fldt 4(%esp) // x
- fabs
- fcompl MO(one) // < 1, == 1, or > 1
+ fldl MO(one) // 1
+ fldt 4(%esp) // x : 1
+ fabs // abs(x) : 1
+ fucompp // < 1, == 1, or > 1
fnstsw
andb $0x45, %ah
cmpb $0x45, %ah
--- libc/sysdeps/i386/fpu/e_powf.S.jj 2005-05-04 19:45:15.000000000 +0200
+++ libc/sysdeps/i386/fpu/e_powf.S 2007-03-01 16:11:59.000000000 +0100
@@ -1,5 +1,5 @@
/* ix87 specific implementation of pow function.
- Copyright (C) 1996, 1997, 1999, 2001, 2004, 2005
+ Copyright (C) 1996, 1997, 1999, 2001, 2004, 2005, 2007
Free Software Foundation, Inc.
This file is part of the GNU C Library.
Contributed by Ulrich Drepper <drepper@cygnus.com>, 1996.
@@ -155,10 +155,11 @@ ENTRY(__ieee754_powf)
2: /* y is a real number. */
fxch // x : y
fldl MO(one) // 1.0 : x : y
- fld %st(1) // x : 1.0 : x : y
- fsub %st(1) // x-1 : 1.0 : x : y
- fabs // |x-1| : 1.0 : x : y
- fcompl MO(limit) // 1.0 : x : y
+ fldl MO(limit) // 0.29 : 1.0 : x : y
+ fld %st(2) // x : 0.29 : 1.0 : x : y
+ fsub %st(2) // x-1 : 0.29 : 1.0 : x : y
+ fabs // |x-1| : 0.29 : 1.0 : x : y
+ fucompp // 1.0 : x : y
fnstsw
fxch // x : 1.0 : y
sahf
@@ -191,9 +192,10 @@ ENTRY(__ieee754_powf)
// y == ±inf
.align ALIGNARG(4)
12: fstp %st(0) // pop y
- flds 4(%esp) // x
- fabs
- fcompl MO(one) // < 1, == 1, or > 1
+ fldl MO(one) // 1
+ flds 4(%esp) // x : 1
+ fabs // abs(x) : 1
+ fucompp // < 1, == 1, or > 1
fnstsw
andb $0x45, %ah
cmpb $0x45, %ah
--- libc/sysdeps/ieee754/dbl-64/e_pow.c.jj 2004-01-23 14:16:24.000000000 +0100
+++ libc/sysdeps/ieee754/dbl-64/e_pow.c 2007-03-01 14:23:31.000000000 +0100
@@ -106,20 +106,28 @@ double __ieee754_pow(double x, double y)
else
return y < 0 ? 1.0/ABS(x) : 0.0; /* return 0 */
}
+
+ qx = u.i[HIGH_HALF]&0x7fffffff; /* no sign */
+ qy = v.i[HIGH_HALF]&0x7fffffff; /* no sign */
+
+ if (qx >= 0x7ff00000 && (qx > 0x7ff00000 || u.i[LOW_HALF] != 0)) return NaNQ.x;
+ if (qy >= 0x7ff00000 && (qy > 0x7ff00000 || v.i[LOW_HALF] != 0))
+ return x == 1.0 ? 1.0 : NaNQ.x;
+
/* if x<0 */
if (u.i[HIGH_HALF] < 0) {
k = checkint(y);
if (k==0) {
- if ((v.i[HIGH_HALF] & 0x7fffffff) == 0x7ff00000 && v.i[LOW_HALF] == 0) {
+ if (qy == 0x7ff00000) {
if (x == -1.0) return 1.0;
else if (x > -1.0) return v.i[HIGH_HALF] < 0 ? INF.x : 0.0;
else return v.i[HIGH_HALF] < 0 ? 0.0 : INF.x;
}
- else if (u.i[HIGH_HALF] == 0xfff00000 && u.i[LOW_HALF] == 0)
+ else if (qx == 0x7ff00000)
return y < 0 ? 0.0 : INF.x;
return NaNQ.x; /* y not integer and x<0 */
}
- else if (u.i[HIGH_HALF] == 0xfff00000 && u.i[LOW_HALF] == 0)
+ else if (qx == 0x7ff00000)
{
if (k < 0)
return y < 0 ? nZERO.x : nINF.x;
@@ -129,14 +137,6 @@ double __ieee754_pow(double x, double y)
return (k==1)?__ieee754_pow(-x,y):-__ieee754_pow(-x,y); /* if y even or odd */
}
/* x>0 */
- qx = u.i[HIGH_HALF]&0x7fffffff; /* no sign */
- qy = v.i[HIGH_HALF]&0x7fffffff; /* no sign */
-
- if (qx > 0x7ff00000 || (qx == 0x7ff00000 && u.i[LOW_HALF] != 0)) return NaNQ.x;
- /* if 0<x<2^-0x7fe */
- if (qy > 0x7ff00000 || (qy == 0x7ff00000 && v.i[LOW_HALF] != 0))
- return x == 1.0 ? 1.0 : NaNQ.x;
- /* if y<2^-0x7fe */
if (qx == 0x7ff00000) /* x= 2^-0x3ff */
{if (y == 0) return NaNQ.x;
--- libc/sysdeps/x86_64/fpu/e_powl.S.jj 2004-07-20 09:05:58.000000000 +0200
+++ libc/sysdeps/x86_64/fpu/e_powl.S 2007-03-01 15:30:42.000000000 +0100
@@ -1,5 +1,6 @@
/* ix87 specific implementation of pow function.
- Copyright (C) 1996, 1997, 1998, 1999, 2001, 2004 Free Software Foundation, Inc.
+ Copyright (C) 1996, 1997, 1998, 1999, 2001, 2004, 2007
+ Free Software Foundation, Inc.
This file is part of the GNU C Library.
Contributed by Ulrich Drepper <drepper@cygnus.com>, 1996.
@@ -146,10 +147,11 @@ ENTRY(__ieee754_powl)
2: /* y is a real number. */
fxch // x : y
fldl MO(one) // 1.0 : x : y
- fld %st(1) // x : 1.0 : x : y
- fsub %st(1) // x-1 : 1.0 : x : y
- fabs // |x-1| : 1.0 : x : y
- fcompl MO(limit) // 1.0 : x : y
+ fldl MO(limit) // 0.29 : 1.0 : x : y
+ fld %st(2) // x : 0.29 : 1.0 : x : y
+ fsub %st(2) // x-1 : 0.29 : 1.0 : x : y
+ fabs // |x-1| : 0.29 : 1.0 : x : y
+ fucompp // 1.0 : x : y
fnstsw
fxch // x : 1.0 : y
test $4500,%eax
@@ -190,9 +192,10 @@ ENTRY(__ieee754_powl)
// y == ±inf
.align ALIGNARG(4)
12: fstp %st(0) // pop y
- fldt 8(%rsp) // x
- fabs
- fcompl MO(one) // < 1, == 1, or > 1
+ fldl MO(one) // 1
+ fldt 8(%rsp) // x : 1
+ fabs // abs(x) : 1
+ fucompp // < 1, == 1, or > 1
fnstsw
andb $0x45, %ah
cmpb $0x45, %ah
Jakub
| Index Nav: | [Date Index] [Subject Index] [Author Index] [Thread Index] | |
|---|---|---|
| Message Nav: | [Date Prev] [Date Next] | [Thread Prev] [Thread Next] |