about summary refs log tree commit diff
path: root/sysdeps/ieee754/ldbl-128ibm
Commit message (Collapse)AuthorAgeFilesLines
* Fix uninitialized variable use in ldbl-128ibm nearbyintl.Joseph Myers2015-08-201-9/+5
| | | | | | | | | | | | | | | | | | | | | | | | | | | | Removing the use of -Wno-uninitialized for math/ shows errors for ldbl-128ibm: ../sysdeps/ieee754/ldbl-128ibm/s_nearbyintl.c: In function '__nearbyintl': ../sysdeps/ieee754/ldbl-128ibm/s_nearbyintl.c:119:34: error: 'low' may be used uninitialized in this function [-Werror=maybe-uninitialized] u.d[1].d = high - u.d[0].d + low; ^ ../sysdeps/ieee754/ldbl-128ibm/s_nearbyintl.c:119:23: error: 'high' may be used uninitialized in this function [-Werror=maybe-uninitialized] u.d[1].d = high - u.d[0].d + low; ^ These errors are correct: if the high part of the argument is a NaN, and the low part is nonzero but has absolute value less than 2^52, those variables can be used uninitialized. This patch rearranges the code so that the variables are always initialized with the natural values, and then possibly modified later, to avoid this uninitialized use. (Note that there are still other issues with this code and NaNs that are not fixed by this patch.) No bug filed in Bugzilla or testcase added for the uninitialized use since it wasn't user-visible with the compiler I tried (that is, I still got a NaN result). Tested for powerpc. * sysdeps/ieee754/ldbl-128ibm/s_nearbyintl.c: Always initialize variables for high and low parts before possibly modifying them.
* Fix tanh missing underflows (bug 16520).Joseph Myers2015-08-131-0/+8
| | | | | | | | | | | | | | | | | | | | | | | | | | | | | Similar to various other bugs in this area, some tanh implementations do not raise the underflow exception for subnormal arguments, when the result is tiny and inexact. This patch forces the exception in a similar way to previous fixes. Tested for x86_64, x86, mips64 and powerpc. [BZ #16520] * sysdeps/ieee754/dbl-64/s_tanh.c: Include <float.h>. (__tanh): Force underflow exception for arguments with small absolute value. * sysdeps/ieee754/flt-32/s_tanhf.c: Include <float.h>. (__tanhf): Force underflow exception for arguments with small absolute value. * sysdeps/ieee754/ldbl-128/s_tanhl.c: Include <float.h>. (__tanhl): Force underflow exception for arguments with small absolute value. * sysdeps/ieee754/ldbl-128ibm/s_tanhl.c: Include <float.h>. (__tanhl): Force underflow exception for arguments with small absolute value. * sysdeps/ieee754/ldbl-96/s_tanhl.c: Include <float.h>. (__tanhl): Force underflow exception for arguments with small absolute value. * math/auto-libm-test-in: Add more tests of tanh. * math/auto-libm-test-out: Regenerated. * sysdeps/i386/fpu/libm-test-ulps: Update.
* Fix ldbl-128ibm tanhl inaccuracy (bug 18790).Joseph Myers2015-08-101-5/+5
| | | | | | | | | | | | | | | | ldbl-128ibm tanhl uses a too-small threshold to decide when to return +/-1, resulting in large errors. This patch changes it to a more appropriate threshold (the requirement is for 2*exp(-2|x|) to be small in terms of ulps of 1). Tested for x86_64, x86 and powerpc. [BZ #18790] * sysdeps/ieee754/ldbl-128ibm/s_tanhl.c (__tanhl): Increase threshold for returning +/- 1. * math/auto-libm-test-in: Add more tests of tanh. * math/auto-libm-test-out: Regenerated. * sysdeps/i386/fpu/libm-test-ulps: Update.
* Fix ldbl-128ibm sinhl inaccuracy near 0 (bug 18789).Joseph Myers2015-08-101-1/+1
| | | | | | | | | | | | | | | ldbl-128ibm sinhl uses a too-big threshold to decide when to return the argument, resulting in large errors. This patch fixes it to use a more appropriate threshold. Tested for x86_64, x86 and powerpc. [BZ #18789] * sysdeps/ieee754/ldbl-128ibm/e_sinhl.c (__ieee754_sinhl): Use smaller threshold for returning the argument. * math/auto-libm-test-in: Add more tests of sinh. * math/auto-libm-test-out: Regenerated. * sysdeps/i386/fpu/libm-test-ulps: Update.
* Fix tan missing underflows (bug 16517).Joseph Myers2015-08-071-1/+11
| | | | | | | | | | | | | | | | | | | | | | | | | | | | Similar to various other bugs in this area, some tan implementations do not raise the underflow exception for subnormal arguments, when the result is tiny and inexact. This patch forces the exception in a similar way to previous fixes. Tested for x86_64, x86, mips64 and powerpc. [BZ #16517] * sysdeps/ieee754/dbl-64/s_tan.c: Include <float.h>. (tan): Force underflow exception for arguments with small absolute value. * sysdeps/ieee754/flt-32/k_tanf.c: Include <float.h>. (__kernel_tanf): Force underflow exception for arguments with small absolute value. * sysdeps/ieee754/ldbl-128/k_tanl.c: Include <float.h>. (__kernel_tanl): Force underflow exception for arguments with small absolute value. * sysdeps/ieee754/ldbl-128ibm/k_tanl.c: Include <float.h>. (__kernel_tanl): Force underflow exception for arguments with small absolute value. * sysdeps/ieee754/ldbl-96/k_tanl.c: Include <float.h>. (__kernel_tanl): Force underflow exception for arguments with small absolute value. * math/auto-libm-test-in: Add more tests of tan. * math/auto-libm-test-out: Regenerated.
* Fix sinh missing underflows (bug 16519).Joseph Myers2015-08-061-1/+8
| | | | | | | | | | | | | | | | | | | | | | | | | | | | | Similar to various other bugs in this area, some sinh implementations do not raise the underflow exception for subnormal arguments, when the result is tiny and inexact. This patch forces the exception in a similar way to previous fixes. Tested for x86_64, x86, mips64 and powerpc. [BZ #16519] * sysdeps/ieee754/dbl-64/e_sinh.c: Include <float.h>. (__ieee754_sinh): Force underflow exception for arguments with small absolute value. * sysdeps/ieee754/flt-32/e_sinhf.c: Include <float.h>. (__ieee754_sinhf): Force underflow exception for arguments with small absolute value. * sysdeps/ieee754/ldbl-128/e_sinhl.c: Include <float.h>. (__ieee754_sinhl): Force underflow exception for arguments with small absolute value. * sysdeps/ieee754/ldbl-128ibm/e_sinhl.c: Include <float.h>. (__ieee754_sinhl): Force underflow exception for arguments with small absolute value. * sysdeps/ieee754/ldbl-96/e_sinhl.c: Include <float.h>. (__ieee754_sinhl): Force underflow exception for arguments with small absolute value. * math/auto-libm-test-in: Add more tests of sinh. * math/auto-libm-test-out: Regenerated. * sysdeps/i386/fpu/libm-test-ulps: Update.
* Improve tgamma accuracy (bug 18613).Joseph Myers2015-06-291-25/+49
| | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | In non-default rounding modes, tgamma can be slightly less accurate than permitted by glibc's accuracy goals. Part of the problem is error accumulation, addressed in this patch by setting round-to-nearest for internal computations. However, there was also a bug in the code dealing with computing pow (x + n, x + n) where x + n is not exactly representable, providing another source of error even in round-to-nearest mode; it was necessary to address both bugs to get errors for all testcases within glibc's accuracy goals. Given this second fix, accuracy in round-to-nearest mode is also improved (hence regeneration of ulps for tgamma should be from scratch - truncate libm-test-ulps or at least remove existing tgamma entries - so that the expected ulps can be reduced). Some additional complications also arose. Certain tgamma tests should strictly, according to IEEE semantics, overflow or not depending on the rounding mode; this is beyond the scope of glibc's accuracy goals for any function without exactly-determined results, but gen-auto-libm-tests doesn't handle being lax there as it does for underflow. (libm-test.inc also doesn't handle being lax about whether the result in cases very close to the overflow threshold is infinity or a finite value close to overflow, but that doesn't cause problems in this case though I've seen it cause problems with random test generation for some functions.) Thus, spurious-overflow markings, with a comment, are added to auto-libm-test-in (no bug in Bugzilla because the issue is with the testsuite, not a user-visible bug in glibc). And on x86, after the patch I saw ERANGE issues as previously reported by Carlos (see my commentary in <https://sourceware.org/ml/libc-alpha/2015-01/msg00485.html>), which needed addressing by ensuring excess range and precision were eliminated at various points if FLT_EVAL_METHOD != 0. I also noticed and fixed a cosmetic issue where 1.0f was used in long double functions and should have been 1.0L. This completes the move of all functions to testing in all rounding modes with ALL_RM_TEST, so gen-libm-have-vector-test.sh is updated to remove the workaround for some functions not using ALL_RM_TEST. Tested for x86_64, x86, mips64 and powerpc. [BZ #18613] * sysdeps/ieee754/dbl-64/e_gamma_r.c (gamma_positive): Take log of X_ADJ not X when adjusting exponent. (__ieee754_gamma_r): Do intermediate computations in round-to-nearest then adjust overflowing and underflowing results as needed. * sysdeps/ieee754/flt-32/e_gammaf_r.c (gammaf_positive): Take log of X_ADJ not X when adjusting exponent. (__ieee754_gammaf_r): Do intermediate computations in round-to-nearest then adjust overflowing and underflowing results as needed. * sysdeps/ieee754/ldbl-128/e_gammal_r.c (gammal_positive): Take log of X_ADJ not X when adjusting exponent. (__ieee754_gammal_r): Do intermediate computations in round-to-nearest then adjust overflowing and underflowing results as needed. Use 1.0L not 1.0f as numerator of division. * sysdeps/ieee754/ldbl-128ibm/e_gammal_r.c (gammal_positive): Take log of X_ADJ not X when adjusting exponent. (__ieee754_gammal_r): Do intermediate computations in round-to-nearest then adjust overflowing and underflowing results as needed. Use 1.0L not 1.0f as numerator of division. * sysdeps/ieee754/ldbl-96/e_gammal_r.c (gammal_positive): Take log of X_ADJ not X when adjusting exponent. (__ieee754_gammal_r): Do intermediate computations in round-to-nearest then adjust overflowing and underflowing results as needed. Use 1.0L not 1.0f as numerator of division. * math/libm-test.inc (tgamma_test_data): Remove one test. Moved to auto-libm-test-in. (tgamma_test): Use ALL_RM_TEST. * math/auto-libm-test-in: Add one test of tgamma. Mark some other tests of tgamma with spurious-overflow. * math/auto-libm-test-out: Regenerated. * math/gen-libm-have-vector-test.sh: Do not check for START. * sysdeps/i386/fpu/libm-test-ulps: Update. * sysdeps/x86_64/fpu/libm-test-ulps: Likewise.
* Fix j1, jn missing underflows (bug 16559).Joseph Myers2015-06-291-0/+5
| | | | | | | | | | | | | | | | | | | | | | | | | | Similar to various other bugs in this area, j1 and jn implementations can fail to raise the underflow exception when the internal computation is exact although the actual function is inexact. This patch forces the exception in a similar way to other such fixes. (The ldbl-128 / ldbl-128ibm j1l implementation is different and doesn't need a change for this until spurious underflows in it are fixed.) Tested for x86_64, x86, mips64 and powerpc. [BZ #16559] * sysdeps/ieee754/dbl-64/e_j1.c: Include <float.h>. (__ieee754_j1): Force underflow exception for small results. * sysdeps/ieee754/dbl-64/e_jn.c (__ieee754_jn): Likewise. * sysdeps/ieee754/flt-32/e_j1f.c: Include <float.h>. (__ieee754_j1f): Force underflow exception for small results. * sysdeps/ieee754/flt-32/e_jnf.c (__ieee754_jnf): Likewise. * sysdeps/ieee754/ldbl-128/e_jnl.c (__ieee754_jnl): Likewise. * sysdeps/ieee754/ldbl-128ibm/e_jnl.c (__ieee754_jnl): Likewise. * sysdeps/ieee754/ldbl-96/e_j1l.c: Include <float.h>. (__ieee754_j1l): Force underflow exception for small results. * sysdeps/ieee754/ldbl-96/e_jnl.c (__ieee754_jnl): Likewise. * math/auto-libm-test-in: Add more tests of j1 and jn. * math/auto-libm-test-out: Regenerated.
* Use round-to-nearest internally in jn, test with ALL_RM_TEST (bug 18602).Joseph Myers2015-06-251-184/+190
| | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | Some existing jn tests, if run in non-default rounding modes, produce errors above those accepted in glibc, which causes problems for moving tests of jn to use ALL_RM_TEST. This patch makes jn set rounding to-nearest internally, as was done for yn some time ago, then computes the appropriate underflowing value for results that underflowed to zero in to-nearest, and moves the tests to ALL_RM_TEST. It does nothing about the general inaccuracy of Bessel function implementations in glibc, though it should make jn more accurate on average in non-default rounding modes through reduced error accumulation. The recomputation of results that underflowed to zero should as a side-effect fix some cases of bug 16559, where jn just used an exact zero, but that is *not* the goal of this patch and other cases of that bug remain unfixed. (Most of the changes in the patch are reindentation to add new scopes for SET_RESTORE_ROUND*.) Tested for x86_64, x86, powerpc and mips64. [BZ #16559] [BZ #18602] * sysdeps/ieee754/dbl-64/e_jn.c (__ieee754_jn): Set round-to-nearest internally then recompute results that underflowed to zero in the original rounding mode. * sysdeps/ieee754/flt-32/e_jnf.c (__ieee754_jnf): Likewise. * sysdeps/ieee754/ldbl-128/e_jnl.c (__ieee754_jnl): Likewise. * sysdeps/ieee754/ldbl-128ibm/e_jnl.c (__ieee754_jnl): Likewise. * sysdeps/ieee754/ldbl-96/e_jnl.c (__ieee754_jnl): Likewise * math/libm-test.inc (jn_test): Use ALL_RM_TEST. * sysdeps/i386/fpu/libm-test-ulps: Update. * sysdeps/x86_64/fpu/libm-test-ulps: Likewise.
* Fix sin, sincos missing underflows (bug 16526, bug 16538).Joseph Myers2015-06-232-7/+23
| | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | Similar to various other bugs in this area, some sin and sincos implementations do not raise the underflow exception for subnormal arguments, when the result is tiny and inexact. This patch forces the exception in a similar way to previous fixes. Tested for x86_64, x86, mips64 and powerpc. [BZ #16526] [BZ #16538] * sysdeps/ieee754/dbl-64/s_sin.c: Include <float.h>. (__sin): Force underflow exception for arguments with small absolute value. * sysdeps/ieee754/flt-32/k_sinf.c: Include <float.h>. (__kernel_sinf): Force underflow exception for arguments with small absolute value. * sysdeps/ieee754/ldbl-128/k_sincosl.c: Include <float.h>. (__kernel_sincosl): Force underflow exception for arguments with small absolute value. * sysdeps/ieee754/ldbl-128/k_sinl.c: Include <float.h>. (__kernel_sinl): Force underflow exception for arguments with small absolute value. * sysdeps/ieee754/ldbl-128ibm/k_sincosl.c: Include <float.h>. (__kernel_sincosl): Force underflow exception for arguments with small absolute value. * sysdeps/ieee754/ldbl-128ibm/k_sinl.c: Include <float.h>. (__kernel_sinl): Force underflow exception for arguments with small absolute value. * sysdeps/ieee754/ldbl-96/k_sinl.c: Include <float.h>. (__kernel_sinl): Force underflow exception for arguments with small absolute value. * sysdeps/powerpc/fpu/k_sinf.c: Include <float.h>. (__kernel_sinf): Force underflow exception for arguments with small absolute value. * math/auto-libm-test-in: Add more tests of sin and sincos. * math/auto-libm-test-out: Regenerated.
* Fix asinh missing underflows (bug 16350).Joseph Myers2015-06-181-0/+6
| | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | Similar to various other bugs in this area, some asinh implementations do not raise the underflow exception for subnormal arguments, when the result is tiny and inexact. This patch forces the exception in a similar way to previous fixes. Tested for x86_64, x86 and mips64. [BZ #16350] * sysdeps/i386/fpu/s_asinh.S (__asinh): Force underflow exception for arguments with small absolute value. * sysdeps/i386/fpu/s_asinhf.S (__asinhf): Likewise. * sysdeps/i386/fpu/s_asinhl.S (__asinhl): Likewise. * sysdeps/ieee754/dbl-64/s_asinh.c: Include <float.h>. (__asinh): Force underflow exception for arguments with small absolute value. * sysdeps/ieee754/flt-32/s_asinhf.c: Include <float.h>. (__asinhf): Force underflow exception for arguments with small absolute value. * sysdeps/ieee754/ldbl-128/s_asinhl.c: Include <float.h>. (__asinhl): Force underflow exception for arguments with small absolute value. * sysdeps/ieee754/ldbl-128ibm/s_asinhl.c: Include <float.h>. (__asinhl): Force underflow exception for arguments with small absolute value. * sysdeps/ieee754/ldbl-96/s_asinhl.c: Include <float.h>. (__asinhl): Force underflow exception for arguments with small absolute value. * math/auto-libm-test-in: Do not mark underflow exceptions as possibly missing for bug 16350. * math/auto-libm-test-out: Regenerated.
* Remove ldbl-128ibm variants of complex math functions.Joseph Myers2015-06-173-281/+0
| | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | sysdeps/ieee754/ldbl-128ibm has its own versions of cprojl, ctanhl and ctanl. Having its own versions, where otherwise the math/ copies are generally used for all floating-point formats, means they are liable to get out of sync and not benefit from bug fixes to the generic versions. The substantive differences (not arising from getting out of sync and slightly different fixes for the same issues) are: long double compat handling (also done in the ldbl-opt versions, so doesn't require special versions for ldbl-128ibm); handling of LDBL_EPSILON (conditionally undefined and redefined in other math/ implementations, so doesn't justify a special version), and: /* __gcc_qmul does not respect -0.0 so we need the following fixup. */ if ((__real__ res == 0.0L) && (__real__ x == 0.0L)) __real__ res = __real__ x; if ((__real__ res == 0.0L) && (__imag__ x == 0.0L)) __imag__ res = __imag__ x; But if that statement about __gcc_qmul was ever true for an old version of that libgcc function, it's not the case for any GCC version now supported to build glibc; there's explicit logic early in that function (and similarly in __gcc_qdiv) to return an appropriately signed zero if the product of the high parts is zero. So this patch adds the special LDBL_EPSILON handling to the generic functions and removes the ldbl-128ibm versions. Tested for powerpc32 (compared test-ldouble.out before and after the changes; there are slight changes to results for ctanl / ctanhl, arising from divergence of the implementations, but nothing that affects the overall nature of the issues shown by the testsuite, and in particular nothing related to signs of zero resutls). * math/s_ctanhl.c [LDBL_MANT_DIG == 106] (LDBL_EPSILON): Undefine and redefine. * math/s_ctanl.c [LDBL_MANT_DIG == 106] (LDBL_EPSILON): Undefine and redefine. * sysdeps/ieee754/ldbl-128ibm/s_cprojl.c: Remove file. * sysdeps/ieee754/ldbl-128ibm/s_ctanhl.c: Likewise. * sysdeps/ieee754/ldbl-128ibm/s_ctanl.c: Likewise.
* Replace finite with isfinite.Wilco Dijkstra2015-06-031-1/+1
|
* This patch renames all uses of __isinf*, __isnan*, __finite* and __signbit* ↵Wilco Dijkstra2015-06-038-10/+10
| | | | to use standard C99 macros. This has no effect on generated code.
* Fix ldbl-128 / ldbl-128ibm tanl for -Wuninitialized.Joseph Myers2015-05-221-0/+12
| | | | | | | | | | | | | | | The ldbl-128 and ldbl-128ibm implementations of tanl produce uninitialized variable warnings with -Wuninitialized because of a variable that is initialized only conditionally, then used under the same conditions under which it is set. This patch uses DIAG_* macros to suppress those warnings. Tested for powerpc and mips64. * sysdeps/ieee754/ldbl-128/k_tanl.c: Include <libc-internal.h>. (__kernel_tanl): Ignore uninitialized warnings around use of SIGN. * sysdeps/ieee754/ldbl-128ibm/k_tanl.c: Include <libc-internal.h>. (__kernel_tanl): Ignore uninitialized warnings around use of SIGN.
* Fix ldbl-128 / ldbl-128ibm erfcl for -WuninitializedJoseph Myers2015-05-221-1/+1
| | | | | | | | | | | | | | | | | | | The ldbl-128 and ldbl-128ibm implementations of erfcl produce uninitialized variable warnings with -Wuninitialized because of switch statements where in fact one of the cases will always be executed, but the compiler does not see that these cases cover all possibilities (and because the reasoning that it does involves inequalities on the representation of a floating point value leading to a set of possible values for 8.0 times that value, converted to int, it's highly nontrivial for the compiler to see that). This patch fixes those warnings by converting the last case in those switch statements to a "default" case. Tested for powerpc and mips64. * sysdeps/ieee754/ldbl-128/s_erfl.c (__erfcl): Make case 9 in switch statement into default case. * sysdeps/ieee754/ldbl-128ibm/s_erfl.c (__erfcl): Likewise.
* Fix ldbl-128 / ldbl-128ibm asinl for -Wuninitialized.Joseph Myers2015-05-221-2/+3
| | | | | | | | | | | | | | | | | | | | | | | | | | | The ldbl-128 and ldbl-128ibm implementations of asinl produce uninitialized variable warnings with -Wuninitialized because the code for small arguments in fact always returns but the compiler cannot see this and instead sees that a variable would be uninitialized if the "if (huge + x > one)" conditional used to force the "inexact" exception were false. All the code in libm trying to force "inexact" for functions that are not exactly defined is suspect and should be removed at some point given that we now have a clear definition of the accuracy goals for libm functions which, following C99/C11, does not require anything about "inexact" for most functions (likewise, the multi-precision code that tries to give correctly-rounded results, very slowly, for functions for which the goals clearly do not include correct rounding, if the faster paths are accurate enough). However, for now this patch simply changes the code to use math_force_eval, rather than "if", to ensure the evaluation of the inexact computation. Tested for powerpc and mips64. * sysdeps/ieee754/ldbl-128/e_asinl.c (__ieee754_asinl): Don't use a conditional in forcing "inexact". * sysdeps/ieee754/ldbl-128ibm/e_asinl.c (__ieee754_asinl): Likewise.
* Fix atanhl missing underflows (bug 16352).Joseph Myers2015-05-151-1/+10
| | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | Similar to various other bugs in this area, some atanh implementations do not raise the underflow exception for subnormal arguments, when the result is tiny and inexact. This patch forces the exception in a similar way to previous fixes. (No change in this regard is needed for the i386 implementation; special handling to force underflows in these cases will only be needed there when the spurious underflows, bug 18049, get fixed.) Tested for x86_64, x86, powerpc and mips64. [BZ #16352] * sysdeps/i386/fpu/e_atanh.S (dbl_min): New object. (__ieee754_atanh): Force underflow exception for results with small absolute value. * sysdeps/i386/fpu/e_atanhf.S (flt_min): New object. (__ieee754_atanhf): Force underflow exception for results with small absolute value. * sysdeps/ieee754/dbl-64/e_atanh.c: Include <float.h>. (__ieee754_atanh): Force underflow exception for results with small absolute value. * sysdeps/ieee754/flt-32/e_atanhf.c: Include <float.h>. (__ieee754_atanhf): Force underflow exception for results with small absolute value. * sysdeps/ieee754/ldbl-128/e_atanhl.c: Include <float.h>. (__ieee754_atanhl): Force underflow exception for results with small absolute value. * sysdeps/ieee754/ldbl-128ibm/e_atanhl.c: Include <float.h>. (__ieee754_atanhl): Force underflow exception for results with small absolute value. * sysdeps/ieee754/ldbl-96/e_atanhl.c: Include <float.h>. (__ieee754_atanhl): Force underflow exception for results with small absolute value. * math/auto-libm-test-in: Do not allow missing underflow exceptions from atanh. * math/auto-libm-test-out: Regenerated.
* Set errno for log1p on pole/domain error.Stefan Liebler2015-04-132-2/+23
| | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | According to bug 6792, errno is not set to ERANGE/EDOM by calling log1p/log1pf/log1pl with x = -1 or x < -1. This patch adds a wrapper which sets errno in those cases and returns the value of the existing __log1p function. The log1p is now an alias to the wrapper function instead of __log1p. The files in sysdeps are reflecting these changes. The ia64 implementation sets errno by itself, thus the wrapper-file is empty. The libm-test is adjusted for log1p-tests to check errno. [BZ #6792] * math/w_log1p.c: New file. * math/w_log1pf.c: Likewise. * math/w_log1pl.c: Likewise. * math/Makefile (libm-calls): Add w_log1p. * math/s_log1pl.c (log1pl): Remove weak_alias. * sysdeps/i386/fpu/s_log1p.S (log1p): Likewise. * sysdeps/i386/fpu/s_log1pf.S (log1pf): Likewise. * sysdeps/i386/fpu/s_log1pl.S (log1pl): Likewise. * sysdeps/x86_64/fpu/s_log1pl.S (log1pl): Likewise. * sysdeps/ieee754/dbl-64/s_log1p.c (log1p): Likewise. [NO_LONG_DOUBLE] (log1pl): Likewise. * sysdeps/ieee754/flt-32/s_log1pf.c (log1pf): Likewise. * sysdeps/ieee754/ldbl-128/s_log1pl.c (log1pl): Likewise. * sysdeps/ieee754/ldbl-64-128/s_log1pl.c (log1p): Remove long_double_symbol. * sysdeps/ieee754/ldbl-128ibm/s_log1pl.c (log1pl): Likewise. * sysdeps/ieee754/ldbl-64-128/w_log1pl.c: New file. * sysdeps/ieee754/ldbl-128ibm/w_log1pl.c: Likewise. * sysdeps/m68k/m680x0/fpu/s_log1p.c: Define empty weak_alias to remove weak_alias for corresponding log1p function. * sysdeps/m68k/m680x0/fpu/s_log1pf.c: Likewise. * sysdeps/m68k/m680x0/fpu/s_log1pl.c: Likewise. * sysdeps/ia64/fpu/w_log1p.c: New file. * sysdeps/ia64/fpu/w_log1pf.c: Likewise. * sysdeps/ia64/fpu/w_log1pl.c: Likewise. * math/libm-test.inc (log1p_test_data): Add errno expectations.
* Fix ldbl-96, ldbl-128ibm atanhl inaccuracy (bug 18046, bug 18047).Joseph Myers2015-02-271-1/+1
| | | | | | | | | | | | | | | | | | | | | | | | | The threshold in ldbl-96 atanhl for when to return the argument, 0x1p-28, is a bit too big, and that in ldbl-128ibm atanhl is much too big (the relevant condition being x^3/3 being < 0.5ulp of x), resulting in errors a bit above the limits of those considered acceptable in glibc in the ldbl-96 case, and in large errors in the ldbl-128ibm case. This patch changes those implementations to use more appropriate thresholds and adds tests around the thresholds for various formats. Tested for x86_64, x86 and powerpc. x86_64 and x86 ulps updated accordingly. [BZ #18046] [BZ #18047] * sysdeps/ieee754/ldbl-128ibm/e_atanhl.c (__ieee754_atanhl): Use 0x1p-56L as threshold for just returning the argument. * sysdeps/ieee754/ldbl-96/e_atanhl.c (__ieee754_atanhl): Use 0x1p-32L as threshold for just returning the argument. * math/auto-libm-test-in: Add more tests of atanh. * math/auto-libm-test-out: Regenerated. * sysdeps/i386/fpu/libm-test-ulps: Update. * sysdeps/x86_64/fpu/libm-test-ulp: Likewise.
* Fix ldbl-128/ldbl-128ibm acosl inaccuracy (bug 18038, bug 18039).Joseph Myers2015-02-261-1/+1
| | | | | | | | | | | | | | | | | | | | | The ldbl-128 and ldbl-128ibm implementations of acosl have similar bugs, using a threshold of 0x1p-57L to determine when they just return pi/2. Since the result pi/2 - asinl (x) is roughly pi/2 - x for small x, the relevant cut-off is actually x being < 0.5ulp of 1. This patch fixes the implementations to use that cut-off and adds tests of small acos arguments. Tested for powerpc and mips64. Also tested for x86_64 and x86; no ulps updates needed. [BZ #18038] [BZ #18039] * sysdeps/ieee754/ldbl-128/e_acosl.c (__ieee754_acosl): Only return pi/2 for arguments below 0x1p-113L. * sysdeps/ieee754/ldbl-128ibm/e_acosl.c (__ieee754_acosl): Only return pi/2 for arguments below 0x1p-106L. * math/auto-libm-test-in: Add more tests of acos. * math/auto-libm-test-out: Regenerated.
* Fix asin missing underflows (bug 16351).Joseph Myers2015-02-261-0/+6
| | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | Similar to various other bugs in this area, some asin implementations do not raise the underflow exception for subnormal arguments, when the result is tiny and inexact. This patch forces the exception in a similar way to previous fixes. Tested for x86_64, x86, powerpc and mips64. [BZ #16351] * sysdeps/i386/fpu/e_asin.S (dbl_min): New object. (MO): New macro. (__ieee754_asin): Force underflow exception for results with small absolute value. * sysdeps/i386/fpu/e_asinf.S (flt_min): New object. (MO): New macro. (__ieee754_asinf): Force underflow exception for results with small absolute value. * sysdeps/ieee754/dbl-64/e_asin.c: Include <float.h> and <math.h>. (__ieee754_asin): Force underflow exception for results with small absolute value. * sysdeps/ieee754/flt-32/e_asinf.c: Include <float.h>. (__ieee754_asinf): Force underflow exception for results with small absolute value. * sysdeps/ieee754/ldbl-128/e_asinl.c: Include <float.h>. (__ieee754_asinl): Force underflow exception for results with small absolute value. * sysdeps/ieee754/ldbl-128ibm/e_asinl.c: Include <float.h>. (__ieee754_asinl): Force underflow exception for results with small absolute value. * sysdeps/ieee754/ldbl-96/e_asinl.c: Include <float.h>. (__ieee754_asinl): Force underflow exception for results with small absolute value. * sysdeps/x86_64/fpu/multiarch/e_asin.c [HAVE_FMA4_SUPPORT]: Include <math.h>. * math/auto-libm-test-in: Do not mark underflow exceptions as possibly missing for bug 16351. * math/auto-libm-test-out: Regenerated.
* Fix ldbl-128ibm logbl near powers of 2 (bug 18030).Joseph Myers2015-02-261-3/+14
| | | | | | | | | | | | | | | | | The ldbl-128ibm implementation of logbl produces incorrect results when the high part of the argument is a power of 2 and the low part a nonzero number with the opposite sign (and so the returned exponent should be 1 less than that of the high part). For example, logbl (0x1.ffffffffffffffp1L) returns 2 but should return 1. (This is similar to (fixed) bug 16740 for frexpl, and (fixed) bug 18029 for ilogbl.) This patch adds checks for that case. Tested for powerpc. [BZ #18030] * sysdeps/ieee754/ldbl-128ibm/s_logbl.c (__logbl): Adjust exponent of power of 2 down when low part has opposite sign. * math/libm-test.inc (logb_test_data): Add more tests.
* Fix ldbl-128ibm ilogbl near powers of 2 (bug 18029).Joseph Myers2015-02-261-4/+19
| | | | | | | | | | | | | | | | | | The ldbl-128ibm implementation of ilogbl produces incorrect results when the high part of the argument is a power of 2 and the low part a nonzero number with the opposite sign (and so the returned exponent should be 1 less than that of the high part). For example, ilogbl (0x1.ffffffffffffffp1L) returns 2 but should return 1. (This is similar to (fixed) bug 16740 for frexpl, and bug 18030 for logbl.) This patch adds checks for that case. Tested for powerpc. [BZ #18029] * sysdeps/ieee754/ldbl-128ibm/e_ilogbl.c (__ieee754_ilogbl): Adjust exponent of power of 2 down when low part has opposite sign. * math/libm-test.inc (ilogb_test_data): Add more tests.
* Fix ldbl-128ibm asinhl inaccuracy (bug 18020).Joseph Myers2015-02-251-4/+4
| | | | | | | | | | | | | | | | | | | The ldbl-128ibm implementation of asinhl uses cut-offs of 0x1p28 and 0x1p-29 to determine when to use simpler formulas that avoid possible overflow / underflow. Both those cut-offs are inappropriate for this format, resulting in large errors. This patch changes the code to use more appropriate cut-offs of 0x1p56 and 0x1p-56, adding tests around the cut-offs for various floating-point formats. Tested for powerpc. Also tested for x86_64 and x86 and updated ulps. [BZ #18020] * sysdeps/ieee754/ldbl-128ibm/s_asinhl.c (__asinhl): Use 2**56 and 2**-56 not 2**28 and 2**-29 as thresholds for simpler formulas. * math/auto-libm-test-in: Add more tests of asinh. * math/auto-libm-test-out: Regenerated. * sysdeps/i386/fpu/libm-test-ulps: Update. * sysdeps/x86_64/fpu/libm-test-ulps: Likewise.
* Fix ldbl-128ibm acoshl inaccuracy (bug 18019).Joseph Myers2015-02-251-2/+2
| | | | | | | | | | | | | | | | | | The ldbl-128ibm implementation of acoshl uses a cut-off of 0x1p28 to determine when to use log(x) + log(2) as a formula. That cut-off is too small for this format, resulting in large errors. This patch changes it to a more appropriate cut-off of 0x1p56, adding tests around the cut-offs for various floating-point formats. Tested for powerpc. Also tested for x86_64 and x86 and updated ulps. [BZ #18019] * sysdeps/ieee754/ldbl-128ibm/e_acoshl.c (__ieee754_acoshl): Use 2**56 not 2**28 as threshold for log (2x) formula. * math/auto-libm-test-in: Add more tests of acosh. * math/auto-libm-test-out: Regenerated. * sysdeps/i386/fpu/libm-test-ulps: Update. * sysdeps/x86_64/fpu/libm-test-ulps: Likewise.
* Fix atan / atan2 missing underflows (bug 15319).Joseph Myers2015-02-181-0/+6
| | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | This patch fixes bug 15319, missing underflows from atan / atan2 when the result of atan is very close to its small argument (or that of atan2 is very close to the ratio of its arguments, which may be an exact division). The usual approach of doing an underflowing computation if the computed result is subnormal is followed. For 32-bit x86, there are extra complications: the inline __ieee754_atan2 in bits/mathinline.h needs to be disabled for float and double because other libm functions using it generally rely on getting proper underflow exceptions from it, while the out-of-line functions have to remove excess range and precision from the underflowing result so as to return an exact 0 in the case where errno should be set for underflow to 0. (The failures I saw without that are similar to those Carlos reported for other functions, where I haven't seen a response to <https://sourceware.org/ml/libc-alpha/2015-01/msg00485.html> confirming if my diagnosis is correct. Arguably all libm functions with float and double returns should remove excess range and precision, but that's a separate matter.) The x86_64 long double case reported in a comment in bug 15319 is not a bug (it's an argument of LDBL_MIN, and x86_64 is an after-rounding architecture so the correct IEEE result is not to raise underflow in the given rounding mode, in addition to treating the result as an exact LDBL_MIN being within the newly clarified documentation of accuracy goals). I'm presuming that the fpatan instruction can be trusted to raise appropriate exceptions when the (long double) result underflows (after rounding) and so no changes are needed for x86 / x86_64 long double functions here; empirically this is the case for the cases covered in the testsuite, on my system. Tested for x86_64, x86, powerpc and mips64. Only 32-bit x86 needs ulps updates (for the changes to inlines meaning some functions no longer get excess precision from their __ieee754_atan2* calls). [BZ #15319] * sysdeps/i386/fpu/e_atan2.S (dbl_min): New object. (MO): New macro. (__ieee754_atan2): For results with small absolute value, force underflow exception and remove excess range and precision from return value. * sysdeps/i386/fpu/e_atan2f.S (flt_min): New object. (MO): New macro. (__ieee754_atan2f): For results with small absolute value, force underflow exception and remove excess range and precision from return value. * sysdeps/i386/fpu/s_atan.S (dbl_min): New object. (MO): New macro. (__atan): For results with small absolute value, force underflow exception and remove excess range and precision from return value. * sysdeps/i386/fpu/s_atanf.S (flt_min): New object. (MO): New macro. (__atanf): For results with small absolute value, force underflow exception and remove excess range and precision from return value. * sysdeps/ieee754/dbl-64/e_atan2.c: Include <float.h> and <math.h>. (__ieee754_atan2): Force underflow exception for results with small absolute value. * sysdeps/ieee754/dbl-64/s_atan.c: Include <float.h> and <math_private.h>. (atan): Force underflow exception for results with small absolute value. * sysdeps/ieee754/flt-32/s_atanf.c: Include <float.h>. (__atanf): Force underflow exception for results with small absolute value. * sysdeps/ieee754/ldbl-128/s_atanl.c: Include <float.h> and <math.h>. (__atanl): Force underflow exception for results with small absolute value. * sysdeps/ieee754/ldbl-128ibm/s_atanl.c: Include <float.h>. (__atanl): Force underflow exception for results with small absolute value. * sysdeps/x86/fpu/bits/mathinline.h [!__SSE2_MATH__ && !__x86_64__ && __LIBC_INTERNAL_MATH_INLINES] (__ieee754_atan2): Only define inline for long double. * sysdeps/x86_64/fpu/multiarch/e_atan2.c [HAVE_FMA4_SUPPORT || HAVE_AVX_SUPPORT]: Include <math.h>. * math/auto-libm-test-in: Do not mark underflow exceptions as possibly missing for bug 15319. Add more tests of atan2. * math/auto-libm-test-out: Regenerated. * math/libm-test.inc (casin_test_data): Do not mark underflow exceptions as possibly missing for bug 15319. (casinh_test_data): Likewise. * sysdeps/i386/fpu/libm-test-ulps: Update.
* Fix sign of remquo zero remainder in round-downward mode (bug 17987).Joseph Myers2015-02-171-0/+3
| | | | | | | | | | | | | | | | | | | | | | Various remquo implementations produce a zero remainder with the wrong sign (a zero remainder should always have the sign of the first argument, as specified in IEEE 754) in round-downward mode, resulting from the sign of 0 - 0. This patch checks for zero results and fixes their sign accordingly. Tested for x86_64, x86, mips64 and powerpc. [BZ #17987] * sysdeps/ieee754/dbl-64/s_remquo.c (__remquo): Ensure sign of zero result does not depend on the sign resulting from subtraction. * sysdeps/ieee754/dbl-64/wordsize-64/s_remquo.c (__remquo): Likewise. * sysdeps/ieee754/flt-32/s_remquof.c (__remquof): Likewise. * sysdeps/ieee754/ldbl-128/s_remquol.c (__remquol): Likewise. * sysdeps/ieee754/ldbl-128ibm/s_remquol.c (__remquol): Likewise. * sysdeps/ieee754/ldbl-96/s_remquol.c (__remquol): Likewise. * math/libm-test.inc (remquo_test_data): Add more tests.
* Fix remquo spurious overflows (bug 17978).Joseph Myers2015-02-161-2/+2
| | | | | | | | | | | | | | | | | | | | | | | Various remquo implementations, when computing the last three bits of the quotient, have spurious overflows when 4 times the second argument to remquo overflows. These overflows can in turn cause bad results in rounding modes where that overflow results in a finite value. This patch adds tests to avoid the problem multiplications in cases where they would overflow, similar to those that control an earlier multiplication by 8. Tested for x86_64, x86, mips64 and powerpc. [BZ #17978] * sysdeps/ieee754/dbl-64/s_remquo.c (__remquo): Do not form products 4 * y and 2 * y where those would overflow. * sysdeps/ieee754/dbl-64/wordsize-64/s_remquo.c (__remquo): Likewise. * sysdeps/ieee754/flt-32/s_remquof.c (__remquof): Likewise. * sysdeps/ieee754/ldbl-128/s_remquol.c (__remquol): Likewise. * sysdeps/ieee754/ldbl-128ibm/s_remquol.c (__remquol): Likewise. * sysdeps/ieee754/ldbl-96/s_remquol.c (__remquol): Likewise. * math/libm-test.inc (remquo_test_data): Add more tests.
* Fix sincos errno setting (bug 15467).Joseph Myers2015-02-111-0/+3
| | | | | | | | | | | | | | | | | | | | | | | | This patch makes sincos set errno to EDOM when passed an infinity, similarly to sin and cos. Tested for x86_64, x86, powerpc and mips64. I don't know if the architecture-specific implementations for ia64 and m68k might need corresponding fixes. 2015-02-11 Joseph Myers <joseph@codesourcery.com> [BZ #15467] * sysdeps/ieee754/dbl-64/s_sincos.c: Include <errno.h>. (__sincos): Set errno to EDOM for infinite argument. * sysdeps/ieee754/flt-32/s_sincosf.c: Include <errno.h>. (SINCOSF_FUNC): Set errno to EDOM for infinite argument. * sysdeps/ieee754/ldbl-128/s_sincosl.c: Include <errno.h>. (__sincosl): Set errno to EDOM for infinite argument. * sysdeps/ieee754/ldbl-128ibm/s_sincosl.c: Include <errno.h>. (__sincosl): Set errno to EDOM for infinite argument. * sysdeps/ieee754/ldbl-96/s_sincosl.c: Include <errno.h>. (__sincosl): Set errno to EDOM for infinite argument. * math/libm-test.inc (sincos_test_data): Test errno setting.
* Update copyright dates with scripts/update-copyrights.Joseph Myers2015-01-0236-36/+36
|
* Clean up powerpc fegetround / __fegetround inlines.Joseph Myers2014-12-315-8/+8
| | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | The natural fix for some linknamespace test failures, where C90 libm functions call C99 <fenv.h> functions, is to make fe* into weak aliases for __fe* and call __fe* from within libm as needed. To do this, the __fe* names need to be available for that purpose - that is, they must not be used for something other than aliases of fe*. On powerpc, however, __fegetround is an inline function in fenv_libc.h, with no corresponding fegetround inline function; fegetround has an equivalent macro expansion in bits/fenvinline.h, but that is disabled if __NO_MATH_INLINES (which is defined for building libm). I see no need for that disabling; it's not even clear that __NO_MATH_INLINES should affect <fenv.h>, and the results of fegetround are completely defined so there is no semantic effect of that disabling at all outside glibc. The x86 inline feraiseexcept is conditioned on __USE_EXTERN_INLINES not __NO_MATH_INLINES (but that's an inline function rather than a macro). This patch removes the __NO_MATH_INLINES conditional on that fegetround macro, so resulting in it being expanded inline inside glibc. In turn, this means that direct calls to __fegetround from C99 functions in ldbl-128ibm can be changed to calls to fegetround, so that nofpu fenv_libc.h files don't need to define __fegetround at all and, by changing ldbl-128ibm files to use <fenv.h> not <fenv_libc.h>, non-e500 nofpu no longer needs an fenv_libc.h file. The other macros in fenvinline.h are left conditional on __NO_MATH_INLINES, although since the only case where this should make a difference is one involving undefined behavior (if the argument to the function is not a valid exception macro). The out-of-line definition for fegetround uses __fegetround (the inline function removed by this patch). So this continues to work, the fenvinline.h header is made to define __fegetround, and then to define fegetround to call __fegetround. Tested for powerpc32 (hard float) that installed stripped shared libraries are unchanged by this patch; also tested that powerpc-nofpu build still works. (This patch does not itself fix any bugs; it simply cleans things up in preparation for separate bug fixes.) * sysdeps/powerpc/bits/fenvinline.h (fegetround): Rename macro to __fegetround and redefine to call __fegetround. Remove condition on [!__NO_MATH_INLINES]. * sysdeps/powerpc/fpu/fenv_libc.h (__fegetround): Remove inline function. * sysdeps/powerpc/nofpu/fenv_libc.h: Remove file. * sysdeps/powerpc/powerpc32/e500/nofpu/fenv_libc.h (__fegetround): Remove macro. * sysdeps/ieee754/ldbl-128ibm/s_llrintl.c: Include <fenv.h> instead of <fenv_libc.h>. (__llrintl): Call fegetround instead of __fegetround. * sysdeps/ieee754/ldbl-128ibm/s_llroundl.c: Include <fenv.h> instead of <fenv_libc.h>. * sysdeps/ieee754/ldbl-128ibm/s_lrintl.c: Likewise. (__lrintl): Call fegetround instead of __fegetround. * sysdeps/ieee754/ldbl-128ibm/s_lroundl.c: Include <fenv.h> instead of <fenv_libc.h>. * sysdeps/ieee754/ldbl-128ibm/s_rintl.c: Likewise. (__rintl): Call fegetround instead of __fegetround.
* FIx ldbl-128ibm frexpl for 32-bit systems (bug 16619, bug 16740).Joseph Myers2014-11-261-2/+2
| | | | | | | | | | | | | | | | | | | | | | This patch fixes bugs in ldbl-128ibm frexpl for 32-bit systems shown up by warnings: ../sysdeps/ieee754/ldbl-128ibm/s_frexpl.c:82:4: warning: left shift count >= width of type ../sysdeps/ieee754/ldbl-128ibm/s_frexpl.c:129:5: warning: left shift count >= width of type This did in fact show up in test-ldouble.out (alongside all the other problems there ... maybe we should again consider running the libm tests at finer granularity from the makefiles) as already covered by the testsuite after the previous patch that fixed these bugs for 64-bit systems. The fix is simply using 1LL instead of 1L when shifting by 52. Tested for powerpc32 (soft float). [BZ #16619] [BZ #16740] * sysdeps/ieee754/ldbl-128ibm/s_frexpl.c (__frexpl): Use 1LL << 52 instead of 1L << 52.
* Remove IS_IN_libmSiddhesh Poyarekar2014-11-249-9/+9
| | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | Replace with IS_IN (libm). Generated code unchanged on x86_64. * include/math.h: Use IS_IN instead of IS_IN_libm. * sysdeps/alpha/fpu/s_copysign.c: Likewise. * sysdeps/ieee754/ldbl-128ibm/s_copysignl.c: Likewise. * sysdeps/ieee754/ldbl-128ibm/s_finitel.c: Likewise. * sysdeps/ieee754/ldbl-128ibm/s_fmal.c: Likewise. * sysdeps/ieee754/ldbl-128ibm/s_frexpl.c: Likewise. * sysdeps/ieee754/ldbl-128ibm/s_isinfl.c: Likewise. * sysdeps/ieee754/ldbl-128ibm/s_isnanl.c: Likewise. * sysdeps/ieee754/ldbl-128ibm/s_modfl.c: Likewise. * sysdeps/ieee754/ldbl-128ibm/s_scalbnl.c: Likewise. * sysdeps/ieee754/ldbl-128ibm/s_signbitl.c: Likewise. * sysdeps/ieee754/ldbl-64-128/s_copysignl.c: Likewise. * sysdeps/ieee754/ldbl-64-128/s_finitel.c: Likewise. * sysdeps/ieee754/ldbl-64-128/s_frexpl.c: Likewise. * sysdeps/ieee754/ldbl-64-128/s_isinfl.c: Likewise. * sysdeps/ieee754/ldbl-64-128/s_isnanl.c: Likewise. * sysdeps/ieee754/ldbl-64-128/s_modfl.c: Likewise. * sysdeps/ieee754/ldbl-64-128/s_scalbnl.c: Likewise. * sysdeps/ieee754/ldbl-64-128/s_signbitl.c: Likewise. * sysdeps/ieee754/ldbl-64-128/w_scalblnl.c: Likewise. * sysdeps/ieee754/ldbl-opt/s_copysign.c: Likewise. * sysdeps/ieee754/ldbl-opt/s_finite.c: Likewise. * sysdeps/ieee754/ldbl-opt/s_frexp.c: Likewise. * sysdeps/ieee754/ldbl-opt/s_isinf.c: Likewise. * sysdeps/ieee754/ldbl-opt/s_isnan.c: Likewise. * sysdeps/ieee754/ldbl-opt/s_ldexp.c: Likewise. * sysdeps/ieee754/ldbl-opt/s_ldexpl.c: Likewise. * sysdeps/ieee754/ldbl-opt/s_modf.c: Likewise. * sysdeps/ieee754/ldbl-opt/s_scalbln.c: Likewise. * sysdeps/ieee754/ldbl-opt/s_scalbn.c: Likewise. * sysdeps/powerpc/power5+/fpu/s_modf.c: Likewise. * sysdeps/powerpc/powerpc32/fpu/s_copysign.S: Likewise. * sysdeps/powerpc/powerpc32/fpu/s_copysignl.S: Likewise. * sysdeps/powerpc/powerpc32/fpu/s_isnan.S: Likewise. * sysdeps/powerpc/powerpc32/power4/fpu/multiarch/s_copysign.c: Likewise. * sysdeps/powerpc/powerpc32/power4/fpu/multiarch/s_finite.c: Likewise. * sysdeps/powerpc/powerpc32/power4/fpu/multiarch/s_isinf.c: Likewise. * sysdeps/powerpc/powerpc32/power4/fpu/multiarch/s_isnan.c: Likewise. * sysdeps/powerpc/powerpc32/power4/fpu/multiarch/s_modf.c: Likewise. * sysdeps/powerpc/powerpc32/power5/fpu/s_isnan.S: Likewise. * sysdeps/powerpc/powerpc32/power6/fpu/s_copysign.S: Likewise. * sysdeps/powerpc/powerpc32/power6/fpu/s_isnan.S: Likewise. * sysdeps/powerpc/powerpc32/power7/fpu/s_finite.S: Likewise. * sysdeps/powerpc/powerpc32/power7/fpu/s_isinf.S: Likewise. * sysdeps/powerpc/powerpc32/power7/fpu/s_isnan.S: Likewise. * sysdeps/powerpc/powerpc64/fpu/multiarch/s_copysign.c: Likewise. * sysdeps/powerpc/powerpc64/fpu/multiarch/s_finite.c: Likewise. * sysdeps/powerpc/powerpc64/fpu/multiarch/s_isinf.c: Likewise. * sysdeps/powerpc/powerpc64/fpu/multiarch/s_isnan.c: Likewise. * sysdeps/powerpc/powerpc64/fpu/multiarch/s_modf.c: Likewise. * sysdeps/powerpc/powerpc64/fpu/s_copysign.S: Likewise. * sysdeps/powerpc/powerpc64/fpu/s_copysignl.S: Likewise. * sysdeps/powerpc/powerpc64/fpu/s_isnan.S: Likewise. * sysdeps/powerpc/powerpc64/power5/fpu/s_isnan.S: Likewise. * sysdeps/powerpc/powerpc64/power6/fpu/s_copysign.S: Likewise. * sysdeps/powerpc/powerpc64/power6/fpu/s_isnan.S: Likewise. * sysdeps/powerpc/powerpc64/power6x/fpu/s_isnan.S: Likewise. * sysdeps/powerpc/powerpc64/power7/fpu/s_finite.S: Likewise. * sysdeps/powerpc/powerpc64/power7/fpu/s_isinf.S: Likewise. * sysdeps/powerpc/powerpc64/power7/fpu/s_isnan.S: Likewise. * sysdeps/powerpc/powerpc64/power8/fpu/s_finite.S: Likewise. * sysdeps/powerpc/powerpc64/power8/fpu/s_isinf.S: Likewise. * sysdeps/powerpc/powerpc64/power8/fpu/s_isnan.S: Likewise. * sysdeps/sparc/sparc32/fpu/s_signbitl.S: Likewise. * sysdeps/sparc/sparc32/sparcv9/fpu/s_isnan.S: Likewise. * sysdeps/unix/sysv/linux/alpha/fraiseexcpt.S: Likewise.
* Fix yn overflow handling in non-default rounding modes (bug 16561, bug 16562).Joseph Myers2014-06-271-66/+77
| | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | This patch fixes bugs 16561 and 16562, bad results of yn in overflow cases in non-default rounding modes, both because an intermediate overflow in the recurrence does not get detected if the result is not an infinity and because an overflowing result may occur in the wrong sign. The fix is to set FE_TONEAREST mode internally for the parts of the function where such overflows can occur (which includes the call to y1 - where yn is used to compute a Bessel function of order -1, negating the result of y1 isn't correct for overflowing results in directed rounding modes) and then compute an overflowing value in the original rounding mode if the to-nearest result was an infinity. Tested x86_64 and x86 and ulps updated accordingly. Also tested for mips64 and powerpc32 to test the ldbl-128 and ldbl-128ibm changes. (The tests for these bugs were added in my previous y1 patch, so the only thing this patch has to do with the testsuite is enable yn testing in all rounding modes.) [BZ #16561] [BZ #16562] * sysdeps/ieee754/dbl-64/e_jn.c: Include <float.h>. (__ieee754_yn): Set FE_TONEAREST mode internally and then recompute overflowing results in original rounding mode. * sysdeps/ieee754/flt-32/e_jnf.c: Include <float.h>. (__ieee754_ynf): Set FE_TONEAREST mode internally and then recompute overflowing results in original rounding mode. * sysdeps/ieee754/ldbl-128/e_jnl.c: Include <float.h>. (__ieee754_ynl): Set FE_TONEAREST mode internally and then recompute overflowing results in original rounding mode. * sysdeps/ieee754/ldbl-128ibm/e_jnl.c: Include <float.h>. (__ieee754_ynl): Set FE_TONEAREST mode internally and then recompute overflowing results in original rounding mode. * sysdeps/ieee754/ldbl-96/e_jnl.c: Include <float.h>. (__ieee754_ynl): Set FE_TONEAREST mode internally and then recompute overflowing results in original rounding mode. * sysdeps/i386/fpu/fenv_private.h [!__SSE2_MATH__] (libc_feholdsetround_ctx): New macro. * math/libm-test.inc (yn_test): Use ALL_RM_TEST. * sysdeps/i386/fpu/libm-test-ulps: Update. * sysdeps/x86_64/fpu/libm-test-ulps : Likewise.
* Fix exp10 spurious underflows (bug 16560).Joseph Myers2014-06-251-0/+2
| | | | | | | | | | | | | | | | | | | | | | | | | This patch fixes spurious underflows from exp10 for arguments near 0 (part of bug 16560; that bug also includes spurious underflows from exp2, which are not fixed by this patch). The problem is underflows in the internal computation converting the exp10 argument to arguments for exp (with extra precision), and the fix is simply to return 1 early for arguments near enough to 0 (just as arguments with large enough magnitude have their own overflow / underflow logic at the start of the function). Tested x86_64 and x86 and ulps updated accordingly; also tested for powerpc32 and mips64 to validate the ldbl-128ibm and ldbl-128 changes. [BZ #16560] * sysdeps/ieee754/dbl-64/e_exp10.c (__ieee754_exp10): Return 1 for arguments close to 0. * sysdeps/ieee754/ldbl-128/e_exp10l.c (__ieee754_exp10l): Likewise. * sysdeps/ieee754/ldbl-128ibm/e_exp10l.c (__ieee754_exp10l): Likewise. * math/auto-libm-test-in: Add more tests of exp10. * math/auto-libm-test-out: Regenerated. * sysdeps/x86_64/fpu/libm-test-ulps: Update.
* Fix cosh spurious underflows from expm1 (bug 16354), inaccurate results near ↵Joseph Myers2014-06-231-1/+1
| | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | 0 (bug 17061). This patch fixes bug 16354, spurious underflows from cosh when a tiny argument is passed to expm1 and expm1 correctly underflows although the final result of cosh should be 1. As noted in that bug, some cases are latent because of expm1 implementations not raising underflow (bug 16353), but all the implementations are fixed similarly. They already contained checks for tiny arguments, but the checks were too late to avoid underflow from expm1 (although they would avoid underflow from subsequent squaring of the result of expm1); they are moved before the expm1 calls. The thresholds used for considering arguments tiny are not particularly consistent in how they relate to the precision of the floating-point format in question. They are, however, all sufficient to ensure that the round-to-nearest result of cosh is indeed 1 below the threshold (although sometimes they are smaller than necessary). But the previous logic did not return 1, but the previously computed 1 + expm1(abs(x)) value. And the thresholds in the ldbl-128 and ldbl-128ibm code (0x1p-71L - I suspect 0x3f8b was intended in the code instead of 0x3fb8 - and (roughly) 0x1p-55L) are not sufficient for that value to be 1. So by moving the test for tiny arguments, and consequently returning 1 directly now the expm1 value hasn't been computed by that point, this patch also fixes bug 17061, the (large number of ulps) inaccuracy for small arguments in those implementations. Tests for that bug are duly added. Tested x86_64 and x86 and ulps updated accordingly. Also tested for mips64 and powerpc32 to validate the ldbl-128 and ldbl-128ibm changes. [BZ #16354] [BZ #17061] * sysdeps/ieee754/dbl-64/e_cosh.c (__ieee754_cosh): Check for small arguments before calling __expm1. * sysdeps/ieee754/flt-32/e_coshf.c (__ieee754_coshf): Check for small arguments before calling __expm1f. * sysdeps/ieee754/ldbl-128/e_coshl.c (__ieee754_coshl): Check for small arguments before calling __expm1l. * sysdeps/ieee754/ldbl-128ibm/e_coshl.c (__ieee754_coshl): Likewise. * sysdeps/ieee754/ldbl-96/e_coshl.c (__ieee754_coshl): Likewise. * math/auto-libm-test-in: Add more cosh tests. Do not allow spurious underflow for some cosh tests. * math/auto-libm-test-out: Regenerated. * sysdeps/i386/fpu/libm-test-ulps: Update.
* [BZ #6803] Set errno for scalbln, scalbnStefan Liebler2014-06-202-1/+23
| | | | | | | | | | | Errno is not set and the testcases will fail. Now the scalbln-aliases are removed in i386/m68 and the wrappers are used when calling the scalbln-functions. On ia64 only scalblnf has its own implementation. For scalbln and scalblnl the ieee754/dbl-64 and ieee754/ldbl-96 are used, thus the wrappers are needed, too.
* PowerPC: Fix nearbyintl failure for few inputsRajalakshmi Srinivasaraghavan2014-06-171-0/+5
| | | | | | | | This patch fixes few failures in nearbyintl() where the fraction part is close to 0.5.i The new tests added report few extra failures in nearbyint_downward and nearbyint_towardzero which is a known issue. Fixes #17031.
* Fix log2 (1) in round-downward mode (bug 17042).Joseph Myers2014-06-101-0/+3
| | | | | | | | | | | | | | | | | | | | | | | | | | | As with other issues of this kind, bug 17042 is log2 (1) wrongly returning -0 instead of +0 in round-downward mode because of implementations effectively in terms of log1p (x - 1). This patch fixes the issue in the same way used for log and log10. Tested x86_64 and x86 and ulps updated accordingly. Also tested for mips64 to confirm a fix was needed for ldbl-128 and to validate that fix (also applied to ldbl-128ibm since that version of log2l is essentially the same as the ldbl-128 one). [BZ #17042] * sysdeps/i386/fpu/e_log2.S (__ieee754_log2): Take absolete value when x - 1 is zero. * sysdeps/i386/fpu/e_log2f.S (__ieee754_log2f): Likewise. * sysdeps/i386/fpu/e_log2l.S (__ieee754_log2l): Likewise. * sysdeps/ieee754/ldbl-128/e_log2l.c (__ieee754_log2l): Return 0.0L for an argument of 1.0L. * sysdeps/ieee754/ldbl-128ibm/e_log2l.c (__ieee754_log2l): Likewise. * sysdeps/x86_64/fpu/e_log2l.S (__ieee754_log2l): Take absolute value when x - 1 is zero. * math/libm-test.inc (log2_test): Use ALL_RM_TEST. * sysdeps/i386/fpu/libm-test-ulps: Update. * sysdeps/x86_64/fpu/libm-test-ulps: Likewise.
* Fix log10 (1) in round-downward mode (bug 16977).Joseph Myers2014-05-231-0/+3
| | | | | | | | | | | | | | | | | | | | | | | | | | | As with various other issues of this kind, bug 16977 is log10 (1) wrongly returning -0 rather than +0 in round-downward mode because of an implementation effectively in terms of log1p (x - 1). This patch fixes the issue in the same way used for log. Tested x86_64 and x86 and ulps updated accordingly. Also tested for mips64 to confirm a fix was needed for ldbl-128 and to validate that fix (also applied to ldbl-128ibm since that version of logl is essentially the same as the ldbl-128 one). [BZ #16977] * sysdeps/i386/fpu/e_log10.S (__ieee754_log10): Take absolute value when x - 1 is zero. * sysdeps/i386/fpu/e_log10f.S (__ieee754_log10f): Likewise. * sysdeps/i386/fpu/e_log10l.S (__ieee754_log10l): Likewise. * sysdeps/ieee754/ldbl-128/e_log10l.c (__ieee754_log10l): Return 0.0L for an argument of 1.0L. * sysdeps/ieee754/ldbl-128ibm/e_log10l.c (__ieee754_log10l): Likewise. * sysdeps/x86_64/fpu/e_log10l.S (__ieee754_log10l): Take absolute value when x - 1 is zero. * math/libm-test.inc (log10_test): Use ALL_RM_TEST. * sysdeps/i386/fpu/libm-test-ulps: Update. * sysdeps/x86_64/fpu/libm-test-ulps: Likewise.
* Fix log1pl (LDBL_MAX) in FE_UPWARD mode (bug 16564).Joseph Myers2014-05-141-1/+4
| | | | | | | | | | | | | | | | | | | | | | | Bug 16564 is spurious overflow of log1pl (LDBL_MAX) in FE_UPWARD mode, resulting from log1pl adding 1 to its argument (for arguments not close to 0), which overflows in that mode. This patch fixes this by avoiding adding 1 to large arguments (precisely what counts as large depends on the floating-point format). Tested x86_64 and x86, and spot-checked log1pl tests on mips64 and powerpc64. [BZ #16564] * sysdeps/i386/fpu/s_log1pl.S (__log1pl): Do not add 1 to positive arguments with exponent 65 or above. * sysdeps/ieee754/ldbl-128/s_log1pl.c (__log1pl): Do not add 1 to arguments 0x1p113L or above. * sysdeps/ieee754/ldbl-128ibm/s_log1pl.c (__log1pl): Do not add 1 to arguments 0x1p107L or above. * sysdeps/x86_64/fpu/s_log1pl.S (__log1pl): Do not add 1 to positive arguments with exponent 65 or above. * math/auto-libm-test-in: Add more tests of log1p. * math/auto-libm-test-out: Regenerated.
* Fix erf underflow handling near 0 (bug 16516).Joseph Myers2014-05-141-5/+11
| | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | Bug 16516 reports spurious underflows from erf (for all floating-point types), when the result is close to underflowing but does not actually underflow. erf (x) is about (2/sqrt(pi))*x for x close to 0, so there are subnormal arguments for which it does not underflow. The various implementations do (x + efx*x) (for efx = 2/sqrt(pi) - 1), for greater accuracy than if just using a single multiplication by an approximation to 2/sqrt(pi) (effectively, this way there are a few more bits in the approximation to 2/sqrt(pi)). This can introduce underflows when efx*x underflows even though the final result does not, so a scaled calculation with 8*efx is done in these cases - but 8 is not a big enough scale factor to avoid all such underflows. 16 is (any underflows with a scale factor of 16 would only occur when the final result underflows), so this patch changes the code to use that factor. Rather than recomputing all the values of the efx8 variable, it is removed, leaving it to the compiler's constant folding to compute 16*efx. As such scaling can also lose underflows when the final scaling down happens to be exact, appropriate checks are added to ensure underflow exceptions occur when required in such cases. Tested x86_64 and x86; no ulps updates needed. Also spot-checked for powerpc32 and mips64 to verify the changes to the ldbl-128ibm and ldbl-128 implementations. [BZ #16516] * sysdeps/ieee754/dbl-64/s_erf.c (efx8): Remove variable. (__erf): Scale by 16 instead of 8 in potentially underflowing case. Ensure exception if result actually underflows. * sysdeps/ieee754/flt-32/s_erff.c (efx8): Remove variable. (__erff): Scale by 16 instead of 8 in potentially underflowing case. Ensure exception if result actually underflows. * sysdeps/ieee754/ldbl-128/s_erfl.c: Include <float.h>. (efx8): Remove variable. (__erfl): Scale by 16 instead of 8 in potentially underflowing case. Ensure exception if result actually underflows. * sysdeps/ieee754/ldbl-128ibm/s_erfl.c: Include <float.h>. (efx8): Remove variable. (__erfl): Scale by 16 instead of 8 in potentially underflowing case. Ensure exception if result actually underflows. * sysdeps/ieee754/ldbl-96/s_erfl.c: Include <float.h>. (efx8): Remove variable. (__erfl): Scale by 16 instead of 8 in potentially underflowing case. Ensure exception if result actually underflows. * math/auto-libm-test-in: Add more tests of erf. * math/auto-libm-test-out: Regenerated.
* Correct IBM long double frexpl.Alan Modra2014-04-161-44/+102
| | | | | | | | | | | Besides fixing the bugzilla, this also fixes corner-cases where the high and low double differ greatly in magnitude, and handles a denormal input without resorting to a fp rescale. [BZ #16740] [BZ #16619] * sysdeps/ieee754/ldbl-128ibm/s_frexpl.c (__frexpl): Rewrite. * math/libm-test.inc (frexp_test_data): Add tests.
* Correct IBM long double nextafterl.Alan Modra2014-04-021-16/+33
| | | | | | | | | | | | Fix for values near a power of two, and some tidies. [BZ #16739] * sysdeps/ieee754/ldbl-128ibm/s_nextafterl.c (__nextafterl): Correct output when value is near a power of two. Use int64_t for lx and remove casts. Use decimal rather than hex exponent constants. Don't use long double multiplication when double will suffice. * math/libm-test.inc (nextafter_test_data): Add tests. * NEWS: Add 16739 and 16786 to bug list.
* Fix nextafter overflow in non-default rounding modes (bug 16677).Joseph Myers2014-03-111-4/+10
| | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | ISO C requires the result of nextafter to be independent of the rounding mode, even when underflow or overflow occurs. This patch fixes the bug in various nextafter implementations that, having done an overflowing computation to force an overflow exception (correct), they then return the result of that computation rather than an infinity computed some other way (incorrect, when the overflowing result of arithmetic with that sign and rounding mode is finite but the correct result is infinite) - generally by falling through to existing code to return a value that in fact is correct for this case (but was computed by an integer increment and so without generating the exceptions required). Having fixed the bug, the previously deferred conversion of nextafter testing in libm-test.inc to ALL_RM_TEST is also included. Tested x86_64 and x86; also spot-checked results of nextafter tests for powerpc32 and mips64 to test the ldbl-128ibm and ldbl-128 changes. (The m68k change is untested.) [BZ #16677] * math/s_nextafter.c (__nextafter): Do not return value from overflowing computation. * sysdeps/i386/fpu/s_nextafterl.c (__nextafterl): Likewise. * sysdeps/ieee754/flt-32/s_nextafterf.c (__nextafterf): Likewise. * sysdeps/ieee754/ldbl-128/s_nextafterl.c (__nextafterl): Likewise. * sysdeps/ieee754/ldbl-128ibm/s_nextafterl.c (__nextafterl): Likewise. * sysdeps/m68k/m680x0/fpu/s_nextafterl.c (__nextafterl): Likewise. * math/libm-test.inc (nextafter_test): Use ALL_RM_TEST.
* Use glibc_likely instead __builtin_expect.Ondřej Bílka2014-02-107-19/+19
|
* PowerPC: remove wrong truncl implementation for PowerPC64Adhemerval Zanella2014-01-081-3/+0
| | | | | | | | | | | | | | | The truncl assembly implementation (sysdeps/powerpc/powerpc64/fpu/s_truncl.S) returns wrong results for some inputs where first double is a exact integer and the precision is determined by second long double. Checking on implementation comments and history, I am very confident the assembly implementation was based on a version before commit 5c68d401698a58cf7da150d9cce769fa6679ba5f that fixes BZ#2423 (Errors in long double (ldbl-128ibm) rounding functions in glibc-2.4). By just removing the implementation and make the build select sysdeps/ieee754/ldbl-128ibm/s_truncl.c instead it fixes tgammal issues regarding wrong result sign.
* Fix ldbl-128ibm expm1l on large arguments (bug 16408).Joseph Myers2014-01-081-0/+2
| | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | This patch fixes bug 16408, ldbl-128ibm expm1l returning NaN for some large arguments. The basic problem is that the approach of converting the exponent to the form n * log(2) + y, where -0.5 <= y <= 0.5, then computing 2^n * expm1(y) + (2^n - 1) falls over when 2^n overflows (starting slightly before the point where expm1 overflows, when y is negative and n is the least integer for which 2^n overflows). The ldbl-128 code, and the x86/x86_64 code, make expm1l fall back to expl for large positive arguments to avoid this issue. This patch makes the ldbl-128ibm code do the same. (The problem appears for the particular argument in the testsuite because the ldbl-128ibm code also uses an overflow threshold that's for ldbl-128 and is too big for ldbl-128ibm, but the problem described applies for large non-overflowing cases as well, although during the freeze is not a suitable time for making the expm1 tests cover cases close to overflow more thoroughly.) This leaves some code for large positive arguments in expm1l that is now dead. To keep the code for ldbl-128 and ldbl-128ibm similar, and to avoid unnecessary changes during the freeze, the patch doesn't remove it; instead I propose to file a bug in Bugzilla as a reminder that this code (for overflow, including errno setting, and for arguments of +Inf) is no longer needed and should be removed from both those expm1l implementations. Tested powerpc32. * sysdeps/ieee754/ldbl-128ibm/s_expm1l.c (__expm1l): Use __expl for large positive arguments.
* Fix ldbl-128ibm coshl spurious overflows (bug 16407).Joseph Myers2014-01-071-1/+1
| | | | | | | | | | | | | | | | | | | | | This patch fixes bug 16407, spurious overflows from ldbl-128ibm coshl. The implementation assumed that a high part (reinterpreted as an integer) of the absolute value of the argument of 0x408633ce8fb9f87dLL or more meant overflow, but the actual threshold has high part 0x408633ce8fb9f87eLL (and a negative low part). The patch adjusts the threshold accordingly. sinhl probably has the same issue, but I didn't get that far in adding tests of special cases (such as just below and above overflow) before the freeze and during the freeze is not a suitable time to add them (as they'd require ulps to be regenerated again), so I'm not changing that function for now; when I add more tests of special cases, we'll discover whether sinhl indeed has this problem. Tested powerpc32. * sysdeps/ieee754/ldbl-128ibm/e_coshl.c (__ieee754_coshl): Increase overflow threshold.