diff options
| -rw-r--r-- | ChangeLog | 15 | ||||
| -rw-r--r-- | NEWS | 20 | ||||
| -rw-r--r-- | math/libm-test.inc | 253 | ||||
| -rw-r--r-- | sysdeps/i386/fpu/e_pow.S | 9 | ||||
| -rw-r--r-- | sysdeps/i386/fpu/e_powf.S | 9 | ||||
| -rw-r--r-- | sysdeps/i386/fpu/e_powl.S | 56 | ||||
| -rw-r--r-- | sysdeps/ieee754/k_standard.c | 54 | ||||
| -rw-r--r-- | sysdeps/x86_64/fpu/e_powl.S | 47 |
8 files changed, 436 insertions, 27 deletions
@@ -1,5 +1,20 @@ 2012-03-28 Joseph Myers <joseph@codesourcery.com> + [BZ #369] + [BZ #2678] + [BZ #3866] + * sysdeps/i386/fpu/e_pow.S (__ieee754_pow): Take absolute value of + x for large integer exponent. + * sysdeps/i386/fpu/e_powf.S (__ieee754_powf): Likewise. + * sysdeps/i386/fpu/e_powl.S (__ieee754_powl): Likewise. Adjust + sign of result as needed afterwards. + * sysdeps/x86_64/fpu/e_powl.S (__ieee754_powl): Likewise. + * sysdeps/ieee754/k_standard.c (__kernel_standard): Handle sign of + result for underflowing pow the same as for overflow. + (__kernel_standard_l): Handle powl overflow and underflow here + rather than calling __kernel_standard. + * math/libm-test.inc (pow_test): Add more tests. + [BZ #3868] [BZ #13879] [BZ #13910] @@ -9,16 +9,16 @@ Version 2.16 * The following bugs are resolved with this release: - 174, 350, 411, 2541, 2547, 2548, 2551, 2552, 2553, 2554, 2562, 2563, 2565, - 2566, 2576, 3335, 3868, 3976, 3992, 4026, 4108, 4596, 4822, 5077, 5461, - 5805, 5993, 6471, 6730, 6884, 6907, 6911, 9739, 9902, 10110, 10135, 10140, - 10210, 10545, 10716, 11174, 11322, 11365, 11451, 11494, 12047, 13058, - 13525, 13526, 13527, 13528, 13529, 13530, 13531, 13532, 13533, 13547, - 13551, 13552, 13553, 13555, 13559, 13566, 13583, 13618, 13637, 13656, - 13658, 13673, 13695, 13704, 13706, 13726, 13738, 13786, 13792, 13806, - 13824, 13840, 13841, 13844, 13846, 13851, 13852, 13854, 13871, 13879, - 13883, 13892, 13910, 13911, 13912, 13913, 13915, 13916, 13917, 13918, - 13919, 13920, 13921 + 174, 350, 369, 411, 2541, 2547, 2548, 2551, 2552, 2553, 2554, 2562, 2563, + 2565, 2566, 2576, 2678, 3335, 3866, 3868, 3976, 3992, 4026, 4108, 4596, + 4822, 5077, 5461, 5805, 5993, 6471, 6730, 6884, 6907, 6911, 9739, 9902, + 10110, 10135, 10140, 10210, 10545, 10716, 11174, 11322, 11365, 11451, + 11494, 12047, 13058, 13525, 13526, 13527, 13528, 13529, 13530, 13531, + 13532, 13533, 13547, 13551, 13552, 13553, 13555, 13559, 13566, 13583, + 13618, 13637, 13656, 13658, 13673, 13695, 13704, 13706, 13726, 13738, + 13786, 13792, 13806, 13824, 13840, 13841, 13844, 13846, 13851, 13852, + 13854, 13871, 13879, 13883, 13892, 13910, 13911, 13912, 13913, 13915, + 13916, 13917, 13918, 13919, 13920, 13921 * ISO C11 support: diff --git a/math/libm-test.inc b/math/libm-test.inc index 68f6ef2..32bce45 100644 --- a/math/libm-test.inc +++ b/math/libm-test.inc @@ -5818,6 +5818,259 @@ pow_test (void) TEST_ff_f (pow, -7.49321e+133, -9.80818e+16, 0); #endif + TEST_ff_f (pow, -1.0, -0xffffff, -1.0); + TEST_ff_f (pow, -1.0, -0x1fffffe, 1.0); +#ifndef TEST_FLOAT + TEST_ff_f (pow, -1.0, -0x1.fffffffffffffp+52L, -1.0); + TEST_ff_f (pow, -1.0, -0x1.fffffffffffffp+53L, 1.0); +#endif +#ifdef TEST_LDOUBLE +# if LDBL_MANT_DIG >= 64 + TEST_ff_f (pow, -1.0, -0x1.fffffffffffffffep+63L, -1.0); + TEST_ff_f (pow, -1.0, -0x1.fffffffffffffffep+64L, 1.0); +# endif +# if LDBL_MANT_DIG >= 106 + TEST_ff_f (pow, -1.0, -0x1.ffffffffffffffffffffffffff8p+105L, -1.0); + TEST_ff_f (pow, -1.0, -0x1.ffffffffffffffffffffffffff8p+106L, 1.0); +# endif +# if LDBL_MANT_DIG >= 113 + TEST_ff_f (pow, -1.0, -0x1.ffffffffffffffffffffffffffffp+112L, -1.0); + TEST_ff_f (pow, -1.0, -0x1.ffffffffffffffffffffffffffffp+113L, 1.0); +# endif +#endif + TEST_ff_f (pow, -1.0, -max_value, 1.0); + + TEST_ff_f (pow, -1.0, 0xffffff, -1.0); + TEST_ff_f (pow, -1.0, 0x1fffffe, 1.0); +#ifndef TEST_FLOAT + TEST_ff_f (pow, -1.0, 0x1.fffffffffffffp+52L, -1.0); + TEST_ff_f (pow, -1.0, 0x1.fffffffffffffp+53L, 1.0); +#endif +#ifdef TEST_LDOUBLE +# if LDBL_MANT_DIG >= 64 + TEST_ff_f (pow, -1.0, 0x1.fffffffffffffffep+63L, -1.0); + TEST_ff_f (pow, -1.0, 0x1.fffffffffffffffep+64L, 1.0); +# endif +# if LDBL_MANT_DIG >= 106 + TEST_ff_f (pow, -1.0, 0x1.ffffffffffffffffffffffffff8p+105L, -1.0); + TEST_ff_f (pow, -1.0, 0x1.ffffffffffffffffffffffffff8p+106L, 1.0); +# endif +# if LDBL_MANT_DIG >= 113 + TEST_ff_f (pow, -1.0, 0x1.ffffffffffffffffffffffffffffp+112L, -1.0); + TEST_ff_f (pow, -1.0, 0x1.ffffffffffffffffffffffffffffp+113L, 1.0); +# endif +#endif + TEST_ff_f (pow, -1.0, max_value, 1.0); + + TEST_ff_f (pow, -2.0, 126, 0x1p126); + TEST_ff_f (pow, -2.0, 127, -0x1p127); + TEST_ff_f (pow, -2.0, -126, 0x1p-126); + TEST_ff_f (pow, -2.0, -127, -0x1p-127); + + TEST_ff_f (pow, -2.0, -0xffffff, minus_zero); + TEST_ff_f (pow, -2.0, -0x1fffffe, plus_zero); +#ifndef TEST_FLOAT + TEST_ff_f (pow, -2.0, -0x1.fffffffffffffp+52L, minus_zero); + TEST_ff_f (pow, -2.0, -0x1.fffffffffffffp+53L, plus_zero); +#endif +#ifdef TEST_LDOUBLE +# if LDBL_MANT_DIG >= 64 + TEST_ff_f (pow, -2.0, -0x1.fffffffffffffffep+63L, minus_zero); + TEST_ff_f (pow, -2.0, -0x1.fffffffffffffffep+64L, plus_zero); +# endif +# if LDBL_MANT_DIG >= 106 + TEST_ff_f (pow, -2.0, -0x1.ffffffffffffffffffffffffff8p+105L, minus_zero); + TEST_ff_f (pow, -2.0, -0x1.ffffffffffffffffffffffffff8p+106L, plus_zero); +# endif +# if LDBL_MANT_DIG >= 113 + TEST_ff_f (pow, -2.0, -0x1.ffffffffffffffffffffffffffffp+112L, minus_zero); + TEST_ff_f (pow, -2.0, -0x1.ffffffffffffffffffffffffffffp+113L, plus_zero); +# endif +#endif + TEST_ff_f (pow, -2.0, -max_value, plus_zero); + + TEST_ff_f (pow, -2.0, 0xffffff, minus_infty, OVERFLOW_EXCEPTION); + TEST_ff_f (pow, -2.0, 0x1fffffe, plus_infty, OVERFLOW_EXCEPTION); +#ifndef TEST_FLOAT + TEST_ff_f (pow, -2.0, 0x1.fffffffffffffp+52L, minus_infty, OVERFLOW_EXCEPTION); + TEST_ff_f (pow, -2.0, 0x1.fffffffffffffp+53L, plus_infty, OVERFLOW_EXCEPTION); +#endif +#ifdef TEST_LDOUBLE +# if LDBL_MANT_DIG >= 64 + TEST_ff_f (pow, -2.0, 0x1.fffffffffffffffep+63L, minus_infty, OVERFLOW_EXCEPTION); + TEST_ff_f (pow, -2.0, 0x1.fffffffffffffffep+64L, plus_infty, OVERFLOW_EXCEPTION); +# endif +# if LDBL_MANT_DIG >= 106 + TEST_ff_f (pow, -2.0, 0x1.ffffffffffffffffffffffffff8p+105L, minus_infty, OVERFLOW_EXCEPTION); + TEST_ff_f (pow, -2.0, 0x1.ffffffffffffffffffffffffff8p+106L, plus_infty, OVERFLOW_EXCEPTION); +# endif +# if LDBL_MANT_DIG >= 113 + TEST_ff_f (pow, -2.0, 0x1.ffffffffffffffffffffffffffffp+112L, minus_infty, OVERFLOW_EXCEPTION); + TEST_ff_f (pow, -2.0, 0x1.ffffffffffffffffffffffffffffp+113L, plus_infty, OVERFLOW_EXCEPTION); +# endif +#endif + /* Bug 13873: OVERFLOW exception may be missing. */ + TEST_ff_f (pow, -2.0, max_value, plus_infty, OVERFLOW_EXCEPTION_OK); + + TEST_ff_f (pow, -max_value, 0.5, nan_value, INVALID_EXCEPTION); + TEST_ff_f (pow, -max_value, 1.5, nan_value, INVALID_EXCEPTION); + TEST_ff_f (pow, -max_value, 1000.5, nan_value, INVALID_EXCEPTION); + TEST_ff_f (pow, -max_value, -2, plus_zero); + TEST_ff_f (pow, -max_value, -3, minus_zero); + TEST_ff_f (pow, -max_value, 2, plus_infty, OVERFLOW_EXCEPTION); + TEST_ff_f (pow, -max_value, 3, minus_infty, OVERFLOW_EXCEPTION); + + TEST_ff_f (pow, -max_value, -0xffffff, minus_zero); + TEST_ff_f (pow, -max_value, -0x1fffffe, plus_zero); +#ifndef TEST_FLOAT + TEST_ff_f (pow, -max_value, -0x1.fffffffffffffp+52L, minus_zero); + TEST_ff_f (pow, -max_value, -0x1.fffffffffffffp+53L, plus_zero); +#endif +#ifdef TEST_LDOUBLE +# if LDBL_MANT_DIG >= 64 + TEST_ff_f (pow, -max_value, -0x1.fffffffffffffffep+63L, minus_zero); + TEST_ff_f (pow, -max_value, -0x1.fffffffffffffffep+64L, plus_zero); +# endif +# if LDBL_MANT_DIG >= 106 + TEST_ff_f (pow, -max_value, -0x1.ffffffffffffffffffffffffff8p+105L, minus_zero); + TEST_ff_f (pow, -max_value, -0x1.ffffffffffffffffffffffffff8p+106L, plus_zero); +# endif +# if LDBL_MANT_DIG >= 113 + TEST_ff_f (pow, -max_value, -0x1.ffffffffffffffffffffffffffffp+112L, minus_zero); + TEST_ff_f (pow, -max_value, -0x1.ffffffffffffffffffffffffffffp+113L, plus_zero); +# endif +#endif + /* Bug 13872: spurious OVERFLOW exception may be present. */ + TEST_ff_f (pow, -max_value, -max_value, plus_zero, OVERFLOW_EXCEPTION_OK); + + TEST_ff_f (pow, -max_value, 0xffffff, minus_infty, OVERFLOW_EXCEPTION); + TEST_ff_f (pow, -max_value, 0x1fffffe, plus_infty, OVERFLOW_EXCEPTION); +#ifndef TEST_FLOAT + TEST_ff_f (pow, -max_value, 0x1.fffffffffffffp+52L, minus_infty, OVERFLOW_EXCEPTION); + TEST_ff_f (pow, -max_value, 0x1.fffffffffffffp+53L, plus_infty, OVERFLOW_EXCEPTION); +#endif +#ifdef TEST_LDOUBLE +# if LDBL_MANT_DIG >= 64 + TEST_ff_f (pow, -max_value, 0x1.fffffffffffffffep+63L, minus_infty, OVERFLOW_EXCEPTION); + TEST_ff_f (pow, -max_value, 0x1.fffffffffffffffep+64L, plus_infty, OVERFLOW_EXCEPTION); +# endif +# if LDBL_MANT_DIG >= 106 + TEST_ff_f (pow, -max_value, 0x1.ffffffffffffffffffffffffff8p+105L, minus_infty, OVERFLOW_EXCEPTION); + TEST_ff_f (pow, -max_value, 0x1.ffffffffffffffffffffffffff8p+106L, plus_infty, OVERFLOW_EXCEPTION); +# endif +# if LDBL_MANT_DIG >= 113 + TEST_ff_f (pow, -max_value, 0x1.ffffffffffffffffffffffffffffp+112L, minus_infty, OVERFLOW_EXCEPTION); + TEST_ff_f (pow, -max_value, 0x1.ffffffffffffffffffffffffffffp+113L, plus_infty, OVERFLOW_EXCEPTION); +# endif +#endif + /* Bug 13873: OVERFLOW exception may be missing. */ + TEST_ff_f (pow, -max_value, max_value, plus_infty, OVERFLOW_EXCEPTION_OK); + + TEST_ff_f (pow, -0.5, 126, 0x1p-126); + TEST_ff_f (pow, -0.5, 127, -0x1p-127); + TEST_ff_f (pow, -0.5, -126, 0x1p126); + TEST_ff_f (pow, -0.5, -127, -0x1p127); + + TEST_ff_f (pow, -0.5, -0xffffff, minus_infty, OVERFLOW_EXCEPTION); + TEST_ff_f (pow, -0.5, -0x1fffffe, plus_infty, OVERFLOW_EXCEPTION); +#ifndef TEST_FLOAT + TEST_ff_f (pow, -0.5, -0x1.fffffffffffffp+52L, minus_infty, OVERFLOW_EXCEPTION); + TEST_ff_f (pow, -0.5, -0x1.fffffffffffffp+53L, plus_infty, OVERFLOW_EXCEPTION); +#endif +#ifdef TEST_LDOUBLE +# if LDBL_MANT_DIG >= 64 + TEST_ff_f (pow, -0.5, -0x1.fffffffffffffffep+63L, minus_infty, OVERFLOW_EXCEPTION); + TEST_ff_f (pow, -0.5, -0x1.fffffffffffffffep+64L, plus_infty, OVERFLOW_EXCEPTION); +# endif +# if LDBL_MANT_DIG >= 106 + TEST_ff_f (pow, -0.5, -0x1.ffffffffffffffffffffffffff8p+105L, minus_infty, OVERFLOW_EXCEPTION); + TEST_ff_f (pow, -0.5, -0x1.ffffffffffffffffffffffffff8p+106L, plus_infty, OVERFLOW_EXCEPTION); +# endif +# if LDBL_MANT_DIG >= 113 + TEST_ff_f (pow, -0.5, -0x1.ffffffffffffffffffffffffffffp+112L, minus_infty, OVERFLOW_EXCEPTION); + TEST_ff_f (pow, -0.5, -0x1.ffffffffffffffffffffffffffffp+113L, plus_infty, OVERFLOW_EXCEPTION); +# endif +#endif + /* Bug 13873: OVERFLOW exception may be missing. */ + TEST_ff_f (pow, -0.5, -max_value, plus_infty, OVERFLOW_EXCEPTION_OK); + + TEST_ff_f (pow, -0.5, 0xffffff, minus_zero); + TEST_ff_f (pow, -0.5, 0x1fffffe, plus_zero); +#ifndef TEST_FLOAT + TEST_ff_f (pow, -0.5, 0x1.fffffffffffffp+52L, minus_zero); + TEST_ff_f (pow, -0.5, 0x1.fffffffffffffp+53L, plus_zero); +#endif +#ifdef TEST_LDOUBLE +# if LDBL_MANT_DIG >= 64 + TEST_ff_f (pow, -0.5, 0x1.fffffffffffffffep+63L, minus_zero); + TEST_ff_f (pow, -0.5, 0x1.fffffffffffffffep+64L, plus_zero); +# endif +# if LDBL_MANT_DIG >= 106 + TEST_ff_f (pow, -0.5, 0x1.ffffffffffffffffffffffffff8p+105L, minus_zero); + TEST_ff_f (pow, -0.5, 0x1.ffffffffffffffffffffffffff8p+106L, plus_zero); +# endif +# if LDBL_MANT_DIG >= 113 + TEST_ff_f (pow, -0.5, 0x1.ffffffffffffffffffffffffffffp+112L, minus_zero); + TEST_ff_f (pow, -0.5, 0x1.ffffffffffffffffffffffffffffp+113L, plus_zero); +# endif +#endif + TEST_ff_f (pow, -0.5, max_value, plus_zero); + + TEST_ff_f (pow, -min_value, 0.5, nan_value, INVALID_EXCEPTION); + TEST_ff_f (pow, -min_value, 1.5, nan_value, INVALID_EXCEPTION); + TEST_ff_f (pow, -min_value, 1000.5, nan_value, INVALID_EXCEPTION); + TEST_ff_f (pow, -min_value, -2, plus_infty, OVERFLOW_EXCEPTION); + TEST_ff_f (pow, -min_value, -3, minus_infty, OVERFLOW_EXCEPTION); + TEST_ff_f (pow, -min_value, 1, -min_value); + TEST_ff_f (pow, -min_value, 2, plus_zero); + TEST_ff_f (pow, -min_value, 3, minus_zero); + + TEST_ff_f (pow, -min_value, -0xffffff, minus_infty, OVERFLOW_EXCEPTION); + TEST_ff_f (pow, -min_value, -0x1fffffe, plus_infty, OVERFLOW_EXCEPTION); +#ifndef TEST_FLOAT + TEST_ff_f (pow, -min_value, -0x1.fffffffffffffp+52L, minus_infty, OVERFLOW_EXCEPTION); + TEST_ff_f (pow, -min_value, -0x1.fffffffffffffp+53L, plus_infty, OVERFLOW_EXCEPTION); +#endif +#ifdef TEST_LDOUBLE +# if LDBL_MANT_DIG >= 64 + TEST_ff_f (pow, -min_value, -0x1.fffffffffffffffep+63L, minus_infty, OVERFLOW_EXCEPTION); + TEST_ff_f (pow, -min_value, -0x1.fffffffffffffffep+64L, plus_infty, OVERFLOW_EXCEPTION); +# endif +# if LDBL_MANT_DIG >= 106 + TEST_ff_f (pow, -min_value, -0x1.ffffffffffffffffffffffffff8p+105L, minus_infty, OVERFLOW_EXCEPTION); + TEST_ff_f (pow, -min_value, -0x1.ffffffffffffffffffffffffff8p+106L, plus_infty, OVERFLOW_EXCEPTION); +# endif +# if LDBL_MANT_DIG >= 113 + TEST_ff_f (pow, -min_value, -0x1.ffffffffffffffffffffffffffffp+112L, minus_infty, OVERFLOW_EXCEPTION); + TEST_ff_f (pow, -min_value, -0x1.ffffffffffffffffffffffffffffp+113L, plus_infty, OVERFLOW_EXCEPTION); +# endif +#endif + /* Bug 13873: OVERFLOW exception may be missing. */ + TEST_ff_f (pow, -min_value, -max_value, plus_infty, OVERFLOW_EXCEPTION_OK); + + TEST_ff_f (pow, -min_value, 0xffffff, minus_zero); + TEST_ff_f (pow, -min_value, 0x1fffffe, plus_zero); +#ifndef TEST_FLOAT + TEST_ff_f (pow, -min_value, 0x1.fffffffffffffp+52L, minus_zero); + TEST_ff_f (pow, -min_value, 0x1.fffffffffffffp+53L, plus_zero); +#endif +#ifdef TEST_LDOUBLE +# if LDBL_MANT_DIG >= 64 + TEST_ff_f (pow, -min_value, 0x1.fffffffffffffffep+63L, minus_zero); + TEST_ff_f (pow, -min_value, 0x1.fffffffffffffffep+64L, plus_zero); +# endif +# if LDBL_MANT_DIG >= 106 + TEST_ff_f (pow, -min_value, 0x1.ffffffffffffffffffffffffff8p+105L, minus_zero); + TEST_ff_f (pow, -min_value, 0x1.ffffffffffffffffffffffffff8p+106L, plus_zero); +# endif +# if LDBL_MANT_DIG >= 113 + TEST_ff_f (pow, -min_value, 0x1.ffffffffffffffffffffffffffffp+112L, minus_zero); + TEST_ff_f (pow, -min_value, 0x1.ffffffffffffffffffffffffffffp+113L, plus_zero); +# endif +#endif + /* Bug 13872: spurious OVERFLOW exception may be present. */ + TEST_ff_f (pow, -min_value, max_value, plus_zero, OVERFLOW_EXCEPTION_OK); + END (pow); } diff --git a/sysdeps/i386/fpu/e_pow.S b/sysdeps/i386/fpu/e_pow.S index 1abedf6..b61a946 100644 --- a/sysdeps/i386/fpu/e_pow.S +++ b/sysdeps/i386/fpu/e_pow.S @@ -114,7 +114,7 @@ ENTRY(__ieee754_pow) fucomp %st(1) // y : x fnstsw sahf - jne 2f + jne 3f /* OK, we have an integer value for y. */ popl %eax @@ -157,7 +157,12 @@ ENTRY(__ieee754_pow) cfi_adjust_cfa_offset (8) .align ALIGNARG(4) -2: /* y is a real number. */ +2: /* y is a large integer (so even). */ + fxch // x : y + fabs // |x| : y + fxch // y : x + .align ALIGNARG(4) +3: /* y is a real number. */ fxch // x : y fldl MO(one) // 1.0 : x : y fldl MO(limit) // 0.29 : 1.0 : x : y diff --git a/sysdeps/i386/fpu/e_powf.S b/sysdeps/i386/fpu/e_powf.S index aa58ed2..529a96f 100644 --- a/sysdeps/i386/fpu/e_powf.S +++ b/sysdeps/i386/fpu/e_powf.S @@ -114,7 +114,7 @@ ENTRY(__ieee754_powf) fucomp %st(1) // y : x fnstsw sahf - jne 2f + jne 3f /* OK, we have an integer value for y. */ popl %edx @@ -151,7 +151,12 @@ ENTRY(__ieee754_powf) cfi_adjust_cfa_offset (4) .align ALIGNARG(4) -2: /* y is a real number. */ +2: /* y is a large integer (so even). */ + fxch // x : y + fabs // |x| : y + fxch // y : x + .align ALIGNARG(4) +3: /* y is a real number. */ fxch // x : y fldl MO(one) // 1.0 : x : y fldl MO(limit) // 0.29 : 1.0 : x : y diff --git a/sysdeps/i386/fpu/e_powl.S b/sysdeps/i386/fpu/e_powl.S index c0aa194c..0e7c05b 100644 --- a/sysdeps/i386/fpu/e_powl.S +++ b/sysdeps/i386/fpu/e_powl.S @@ -117,7 +117,7 @@ ENTRY(__ieee754_powl) fucomp %st(1) // y : x fnstsw sahf - jne 2f + jne 3f /* OK, we have an integer value for y. */ popl %eax @@ -160,7 +160,14 @@ ENTRY(__ieee754_powl) cfi_adjust_cfa_offset (8) .align ALIGNARG(4) -2: /* y is a real number. */ +2: // y is a large integer (absolute value at least 1L<<63), but + // may be odd unless at least 1L<<64. So it may be necessary + // to adjust the sign of a negative result afterwards. + fxch // x : y + fabs // |x| : y + fxch // y : |x| + .align ALIGNARG(4) +3: /* y is a real number. */ fxch // x : y fldl MO(one) // 1.0 : x : y fldl MO(limit) // 0.29 : 1.0 : x : y @@ -190,18 +197,51 @@ ENTRY(__ieee754_powl) f2xm1 // 2^fract(y*log2(x))-1 : int(y*log2(x)) faddl MO(one) // 2^fract(y*log2(x)) : int(y*log2(x)) fscale // 2^fract(y*log2(x))*2^int(y*log2(x)) : int(y*log2(x)) - addl $8, %esp - cfi_adjust_cfa_offset (-8) fstp %st(1) // 2^fract(y*log2(x))*2^int(y*log2(x)) - ret + jmp 29f - cfi_adjust_cfa_offset (8) 28: fstp %st(1) // y*log2(x) fldl MO(one) // 1 : y*log2(x) fscale // 2^(y*log2(x)) : y*log2(x) - addl $8, %esp - cfi_adjust_cfa_offset (-8) fstp %st(1) // 2^(y*log2(x)) +29: testb $2, %dh + jz 292f + // x is negative. If y is an odd integer, negate the result. + fldt 24(%esp) // y : abs(result) + fld %st // y : y : abs(result) + fabs // |y| : y : abs(result) + fcompl MO(p64) // y : abs(result) + fnstsw + sahf + jnc 291f + fldl MO(p63) // p63 : y : abs(result) + fxch // y : p63 : abs(result) + fprem // y%p63 : p63 : abs(result) + fstp %st(1) // y%p63 : abs(result) + + // We must find out whether y is an odd integer. + fld %st // y : y : abs(result) + fistpll (%esp) // y : abs(result) + fildll (%esp) // int(y) : y : abs(result) + fucompp // abs(result) + fnstsw + sahf + jne 292f + + // OK, the value is an integer, but is it odd? + popl %eax + cfi_adjust_cfa_offset (-4) + popl %edx + cfi_adjust_cfa_offset (-4) + andb $1, %al + jz 290f // jump if not odd + // It's an odd integer. + fchs +290: ret + cfi_adjust_cfa_offset (8) +291: fstp %st(0) // abs(result) +292: addl $8, %esp + cfi_adjust_cfa_offset (-8) ret // pow(x,±0) = 1 diff --git a/sysdeps/ieee754/k_standard.c b/sysdeps/ieee754/k_standard.c index c3326d9..4e65bb1 100644 --- a/sysdeps/ieee754/k_standard.c +++ b/sysdeps/ieee754/k_standard.c @@ -508,6 +508,9 @@ __kernel_standard(double x, double y, int type) exc.type = UNDERFLOW; exc.name = type < 100 ? "pow" : (type < 200 ? "powf" : "powl"); exc.retval = zero; + y *= 0.5; + if (x < zero && __rint (y) != y) + exc.retval = -zero; if (_LIB_VERSION == _POSIX_) __set_errno (ERANGE); else if (!matherr(&exc)) { @@ -1004,6 +1007,8 @@ long double __kernel_standard_l (long double x, long double y, int type) { double dx, dy; + struct exception exc; + if (isfinite (x)) { long double ax = fabsl (x); @@ -1028,5 +1033,52 @@ __kernel_standard_l (long double x, long double y, int type) } else dy = y; - return __kernel_standard (dx, dy, type); + + switch (type) + { + case 221: + /* powl (x, y) overflow. */ + exc.arg1 = dx; + exc.arg2 = dy; + exc.type = OVERFLOW; + exc.name = "powl"; + if (_LIB_VERSION == _SVID_) + { + exc.retval = HUGE; + y *= 0.5; + if (x < zero && __rintl (y) != y) + exc.retval = -HUGE; + } + else + { + exc.retval = HUGE_VAL; + y *= 0.5; + if (x < zero && __rintl (y) != y) + exc.retval = -HUGE_VAL; + } + if (_LIB_VERSION == _POSIX_) + __set_errno (ERANGE); + else if (!matherr (&exc)) + __set_errno (ERANGE); + return exc.retval; + + case 222: + /* powl (x, y) underflow. */ + exc.arg1 = dx; + exc.arg2 = dy; + exc.type = UNDERFLOW; + exc.name = "powl"; + exc.retval = zero; + y *= 0.5; + if (x < zero && __rintl (y) != y) + exc.retval = -zero; + if (_LIB_VERSION == _POSIX_) + __set_errno (ERANGE); + else if (!matherr (&exc)) + __set_errno (ERANGE); + return exc.retval; + + default: + return __kernel_standard (dx, dy, type); + } } diff --git a/sysdeps/x86_64/fpu/e_powl.S b/sysdeps/x86_64/fpu/e_powl.S index 562791d..0626ce4 100644 --- a/sysdeps/x86_64/fpu/e_powl.S +++ b/sysdeps/x86_64/fpu/e_powl.S @@ -107,7 +107,7 @@ ENTRY(__ieee754_powl) fistpll -8(%rsp) // y : x fildll -8(%rsp) // int(y) : y : x fucomip %st(1),%st // y : x - jne 2f + jne 3f /* OK, we have an integer value for y. */ mov -8(%rsp),%eax @@ -145,7 +145,14 @@ ENTRY(__ieee754_powl) ret .align ALIGNARG(4) -2: /* y is a real number. */ +2: // y is a large integer (absolute value at least 1L<<63), but + // may be odd unless at least 1L<<64. So it may be necessary + // to adjust the sign of a negative result afterwards. + fxch // x : y + fabs // |x| : y + fxch // y : |x| + .align ALIGNARG(4) +3: /* y is a real number. */ fxch // x : y fldl MO(one) // 1.0 : x : y fldl MO(limit) // 0.29 : 1.0 : x : y @@ -176,13 +183,45 @@ ENTRY(__ieee754_powl) faddl MO(one) // 2^fract(y*log2(x)) : int(y*log2(x)) fscale // 2^fract(y*log2(x))*2^int(y*log2(x)) : int(y*log2(x)) fstp %st(1) // 2^fract(y*log2(x))*2^int(y*log2(x)) - ret + jmp 29f 28: fstp %st(1) // y*log2(x) fldl MO(one) // 1 : y*log2(x) fscale // 2^(y*log2(x)) : y*log2(x) fstp %st(1) // 2^(y*log2(x)) - ret +29: testb $2, %dh + jz 292f + // x is negative. If y is an odd integer, negate the result. + fldt 24(%rsp) // y : abs(result) + fldl MO(p64) // 1L<<64 : y : abs(result) + fld %st(1) // y : 1L<<64 : y : abs(result) + fabs // |y| : 1L<<64 : y : abs(result) + fcomip %st(1), %st // 1L<<64 : y : abs(result) + fstp %st(0) // y : abs(result) + jnc 291f + fldl MO(p63) // p63 : y : abs(result) + fxch // y : p63 : abs(result) + fprem // y%p63 : p63 : abs(result) + fstp %st(1) // y%p63 : abs(result) + + // We must find out whether y is an odd integer. + fld %st // y : y : abs(result) + fistpll -8(%rsp) // y : abs(result) + fildll -8(%rsp) // int(y) : y : abs(result) + fucomip %st(1),%st // y : abs(result) + ffreep %st // abs(result) + jne 292f + + // OK, the value is an integer, but is it odd? + mov -8(%rsp), %eax + mov -4(%rsp), %edx + andb $1, %al + jz 290f // jump if not odd + // It's an odd integer. + fchs +290: ret +291: fstp %st(0) // abs(result) +292: ret // pow(x,±0) = 1 .align ALIGNARG(4) |