about summary refs log tree commit diff
path: root/sysdeps/ieee754
Commit message (Collapse)AuthorAgeFilesLines
* Fix tgamma missing underflows (bug 18951).Joseph Myers2015-09-175-0/+25
| | | | | | | | | | | | | | | | | | | | | | | | Similar to various other bugs in this area, tgamma functions can fail to raise the underflow exception when the result is tiny and inexact but one or more low bits of the intermediate result that is scaled down are zero. This patch forces the exception in a similar way to previous fixes. Tested for x86_64, x86, mips64 and powerpc. [BZ #18951] * sysdeps/ieee754/dbl-64/e_gamma_r.c (__ieee754_gamma_r): Force underflow exception for small results. * sysdeps/ieee754/flt-32/e_gammaf_r.c (__ieee754_gammaf_r): Likewise. * sysdeps/ieee754/ldbl-128/e_gammal_r.c (__ieee754_gammal_r): Likewise. * sysdeps/ieee754/ldbl-128ibm/e_gammal_r.c (__ieee754_gammal_r): Likewise. * sysdeps/ieee754/ldbl-96/e_gammal_r.c (__ieee754_gammal_r): Likewise. * math/auto-libm-test-in: Add more tests of tgamma. * math/auto-libm-test-out: Regenerated.
* Make scalbn set errno (bug 6803).Joseph Myers2015-09-167-16/+3
| | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | As noted in bug 6803, scalbn fails to set errno on overflow and underflow. This patch fixes this by making scalbn an alias of ldexp, which has exactly the same semantics (for floating-point types with radix 2) and already has wrappers that deal with setting errno, instead of an alias of the internal __scalbn (which ldexp calls). Notes: * Where compat symbols were defined for scalbn functions, I didn't change what they point to (to keep the patch minimal), so such compat symbols continue to go directly to the non-errno-setting functions. * Mike, I didn't do anything with the IA64 versions of these functions, where I think both the ldexp and scalbn functions already deal with setting errno. As a cleanup (not needed to fix this bug) however you might want to make those functions into aliases for IA64; there is no need for them to be separate function implementations at all. * This concludes the fix for bug 6803 since the scalb and scalbln cases of that bug were fixed some time ago. Tested for x86_64, x86, mips64 and powerpc. [BZ #6803] * math/s_ldexp.c (scalbn): Define as weak alias of __ldexp. [NO_LONG_DOUBLE] (scalbnl): Define as weak alias of __ldexp. * math/s_ldexpf.c (scalbnf): Define as weak alias of __ldexpf. * math/s_ldexpl.c (scalbnl): Define as weak alias of __ldexpl. * sysdeps/i386/fpu/s_scalbn.S (scalbn): Remove alias. * sysdeps/i386/fpu/s_scalbnf.S (scalbnf): Likewise. * sysdeps/i386/fpu/s_scalbnl.S (scalbnl): Likewise. * sysdeps/ieee754/dbl-64/s_scalbn.c (scalbn): Likewise. [NO_LONG_DOUBLE] (scalbnl): Likewise. * sysdeps/ieee754/dbl-64/wordsize-64/s_scalbn.c (scalbn): Likewise. [NO_LONG_DOUBLE] (scalbnl): Likewise. * sysdeps/ieee754/flt-32/s_scalbnf.c (scalbnf): Likewise. * sysdeps/ieee754/ldbl-128/s_scalbnl.c (scalbnl): Likewise. * sysdeps/ieee754/ldbl-128ibm/s_scalbnl.c (scalbnl): Remove long_double_symbol calls. * sysdeps/ieee754/ldbl-64-128/s_scalbnl.c (scalbnl): Likewise. * sysdeps/ieee754/ldbl-opt/s_ldexpl.c (__ldexpl_2): Define as strong alias of __ldexpl. (scalbnl): Define using long_double_symbol. * sysdeps/m68k/m680x0/fpu/s_scalbn.c (__CONCATX(scalbn,suffix)): Remove alias. * sysdeps/sparc/sparc64/soft-fp/s_scalbnl.c (scalbnl): Likewise. * sysdeps/x86_64/fpu/s_scalbnl.S (scalbnl): Likewise. * math/libm-test.inc (scalbn_test_data): Add errno expectations. (scalbln_test_data): Add more errno expectations.
* Clean up ldbl-128 / ldbl-128ibm expm1l dead code (bug 16415).Joseph Myers2015-09-162-32/+4
| | | | | | | | | | | | | | | | The ldbl-128 and ldbl-128ibm expm1l implementations have code to handle +Inf and finite arguments above an overflow threshold. Since they now use __expl for large positive arguments to fix other problems, this code is unreachable; this patch removes it. Tested for mips64 and powerpc. [BZ #16415] * sysdeps/ieee754/ldbl-128/s_expm1l.c (maxlog): Remove variable. (__expm1l): Remove code to handle positive infinity and overflow. * sysdeps/ieee754/ldbl-128ibm/s_expm1l.c (maxlog): Remove variable. (__expm1l): Remove code to handle positive infinity and overflow.
* Fix ldbl-128ibm nearbyintl use of signaling comparisons on NaNs (bug 18857).Joseph Myers2015-09-151-1/+3
| | | | | | | | | | | | | | | | | | | | | | The ldbl-128ibm implementation of nearbyintl wrongly uses signaling comparisons such as "if (fabs (u.d[0].d) < TWO52)" on arguments that might be NaNs, when "invalid" exceptions should not be raised. (For hard float, this issue may be hidden by <https://gcc.gnu.org/bugzilla/show_bug.cgi?id=58684>, powerpc GCC wrongly only using unordered comparison instructions.) This patch fixes this by just returning the argument if it is not finite (because of the arbitrary value of the low part of a NaN in IBM long double, there are quite a lot of comparisons that could end up involving a NaN when the argument to nearbyintl is a NaN, so excluding NaN arguments at the start is the simplest and safest fix). Tested for powerpc-nofpu, where it removes failures for spurious "invalid" exceptions from nearbyintl. [BZ #18857] * sysdeps/ieee754/ldbl-128ibm/s_nearbyintl.c (__nearbyintl): Just return non-finite argument without doing ordered comparisons on it.
* Simplify hypotf infinity handling (bug 15918).Joseph Myers2015-09-151-10/+2
| | | | | | | | | | | | | | Bug 15918 points out that the handling of infinities in hypotf can be simplified: it's enough to return the absolute value of the infinite argument without first comparing it to the other argument and possibly returning that other argument's absolute value. This patch makes that cleanup (which should not change how hypotf behaves on any input). Tested for x86_64. [BZ #15918] * sysdeps/ieee754/flt-32/e_hypotf.c (__ieee754_hypotf): Simplify handling of cases where one argument is an infinity.
* Fix exp2 missing underflows (bug 16521).Joseph Myers2015-09-142-2/+18
| | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | Various exp2 implementations in glibc can miss underflow exceptions when the scaling down part of the calculation is exact (or, in the x86 case, when the conversion from extended precision to the target precision is exact). This patch forces the exception in a similar way to previous fixes. The x86 exp2f changes may in fact not be needed for this purpose - it's likely to be the case that no argument of type float has an exp2 result so close to an exact subnormal float value that it equals that value when rounded to 64 bits (even taking account of variation between different x86 implementations). However, they are included for consistency with the changes to exp2 and so as to fix the exp2f part of bug 18875 by ensuring that excess range and precision is removed from underflowing return values. Tested for x86_64, x86 and mips64. [BZ #16521] [BZ #18875] * math/e_exp2l.c (__ieee754_exp2l): Force underflow exception for small results. * sysdeps/i386/fpu/e_exp2.S (dbl_min): New object. (MO): New macro. (__ieee754_exp2): For small results, force underflow exception and remove excess range and precision from return value. * sysdeps/i386/fpu/e_exp2f.S (flt_min): New object. (MO): New macro. (__ieee754_exp2f): For small results, force underflow exception and remove excess range and precision from return value. * sysdeps/i386/fpu/e_exp2l.S (ldbl_min): New object. (MO): New macro. (__ieee754_exp2l): Force underflow exception for small results. * sysdeps/ieee754/dbl-64/e_exp2.c (__ieee754_exp2): Likewise. * sysdeps/ieee754/flt-32/e_exp2f.c (__ieee754_exp2f): Likewise. * sysdeps/x86_64/fpu/e_exp2l.S (ldbl_min): New object. (MO): New macro. (__ieee754_exp2l): Force underflow exception for small results. * math/auto-libm-test-in: Add more tests or exp2. * math/auto-libm-test-out: Regenerated.
* Fix ldbl-128/ldbl-128ibm lgamma spurious "invalid", incorrect signgam (bug ↵Joseph Myers2015-09-111-2/+2
| | | | | | | | | | | | | | | | | | | | 18952). The ldbl-128 / ldbl-128ibm implementation of lgammal converts (the floor of minus) non-integer negative arguments to int to determine the value of signgam. When those values are outside the range of int, this produces spurious "invalid" exceptions and incorrect values of signgam. This patch fixes this by instead determining signgam through comparing half the integer in question to floor of half the integer. Tested for mips64, x86_64 and x86. [BZ #18952] * sysdeps/ieee754/ldbl-128/e_lgammal_r.c (__ieee754_lgammal_r): Do not convert non-integer negative arguments to int to determine the value of signgam. * math/auto-libm-test-in: Add more tests of lgamma. * math/auto-libm-test-out: Regenerated.
* Fix lgamma (negative) inaccuracy (bug 2542, bug 2543, bug 2558).Joseph Myers2015-09-1015-0/+2519
| | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | The existing implementations of lgamma functions (except for the ia64 versions) use the reflection formula for negative arguments. This suffers large inaccuracy from cancellation near zeros of lgamma (near where the gamma function is +/- 1). This patch fixes this inaccuracy. For arguments above -2, there are no zeros and no large cancellation, while for sufficiently large negative arguments the zeros are so close to integers that even for integers +/- 1ulp the log(gamma(1-x)) term dominates and cancellation is not significant. Thus, it is only necessary to take special care about cancellation for arguments around a limited number of zeros. Accordingly, this patch uses precomputed tables of relevant zeros, expressed as the sum of two floating-point values. The log of the ratio of two sines can be computed accurately using log1p in cases where log would lose accuracy. The log of the ratio of two gamma(1-x) values can be computed using Stirling's approximation (the difference between two values of that approximation to lgamma being computable without computing the two values and then subtracting), with appropriate adjustments (which don't reduce accuracy too much) in cases where 1-x is too small to use Stirling's approximation directly. In the interval from -3 to -2, using the ratios of sines and of gamma(1-x) can still produce too much cancellation between those two parts of the computation (and that interval is also the worst interval for computing the ratio between gamma(1-x) values, which computation becomes more accurate, while being less critical for the final result, for larger 1-x). Because this can result in errors slightly above those accepted in glibc, this interval is instead dealt with by polynomial approximations. Separate polynomial approximations to (|gamma(x)|-1)(x-n)/(x-x0) are used for each interval of length 1/8 from -3 to -2, where n (-3 or -2) is the nearest integer to the 1/8-interval and x0 is the zero of lgamma in the relevant half-integer interval (-3 to -2.5 or -2.5 to -2). Together, the two approaches are intended to give sufficient accuracy for all negative arguments in the problem range. Outside that range, the previous implementation continues to be used. Tested for x86_64, x86, mips64 and powerpc. The mips64 and powerpc testing shows up pre-existing problems for ldbl-128 and ldbl-128ibm with large negative arguments giving spurious "invalid" exceptions (exposed by newly added tests for cases this patch doesn't affect the logic for); I'll address those problems separately. [BZ #2542] [BZ #2543] [BZ #2558] * sysdeps/ieee754/dbl-64/e_lgamma_r.c (__ieee754_lgamma_r): Call __lgamma_neg for arguments from -28.0 to -2.0. * sysdeps/ieee754/flt-32/e_lgammaf_r.c (__ieee754_lgammaf_r): Call __lgamma_negf for arguments from -15.0 to -2.0. * sysdeps/ieee754/ldbl-128/e_lgammal_r.c (__ieee754_lgammal_r): Call __lgamma_negl for arguments from -48.0 or -50.0 to -2.0. * sysdeps/ieee754/ldbl-96/e_lgammal_r.c (__ieee754_lgammal_r): Call __lgamma_negl for arguments from -33.0 to -2.0. * sysdeps/ieee754/dbl-64/lgamma_neg.c: New file. * sysdeps/ieee754/dbl-64/lgamma_product.c: Likewise. * sysdeps/ieee754/flt-32/lgamma_negf.c: Likewise. * sysdeps/ieee754/flt-32/lgamma_productf.c: Likewise. * sysdeps/ieee754/ldbl-128/lgamma_negl.c: Likewise. * sysdeps/ieee754/ldbl-128/lgamma_productl.c: Likewise. * sysdeps/ieee754/ldbl-128ibm/lgamma_negl.c: Likewise. * sysdeps/ieee754/ldbl-128ibm/lgamma_productl.c: Likewise. * sysdeps/ieee754/ldbl-96/lgamma_negl.c: Likewise. * sysdeps/ieee754/ldbl-96/lgamma_product.c: Likewise. * sysdeps/ieee754/ldbl-96/lgamma_productl.c: Likewise. * sysdeps/generic/math_private.h (__lgamma_negf): New prototype. (__lgamma_neg): Likewise. (__lgamma_negl): Likewise. (__lgamma_product): Likewise. (__lgamma_productl): Likewise. * math/Makefile (libm-calls): Add lgamma_neg and lgamma_product. * math/auto-libm-test-in: Add more tests of lgamma. * math/auto-libm-test-out: Regenerated. * sysdeps/i386/fpu/libm-test-ulps: Update. * sysdeps/x86_64/fpu/libm-test-ulps: Likewise.
* Move bits/libc-lock.h and bits/libc-lockP.h out of bits/ (bug 14912).Joseph Myers2015-09-081-1/+1
| | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | It was noted in <https://sourceware.org/ml/libc-alpha/2012-09/msg00305.html> that the bits/*.h naming scheme should only be used for installed headers. This patch renames bits/libc-lock.h to plain libc-lock.h and bits/libc-lockP.h to plain libc-lockP.h to follow that convention. Note that I don't know where libc-lockP.h comes from for Hurd (the Hurd libc-lock.h includes libc-lockP.h, but the only libc-lockP.h in the glibc source tree is for NPTL) - some unmerged patch? - but I updated the #include in the Hurd libc-lock.h anyway. Tested for x86_64 (testsuite, and that installed stripped shared libraries are unchanged by the patch). [BZ #14912] * bits/libc-lock.h: Move to ... * sysdeps/generic/libc-lock.h: ...here. (_BITS_LIBC_LOCK_H): Rename macro to _LIBC_LOCK_H. * sysdeps/mach/hurd/bits/libc-lock.h: Move to ... * sysdeps/mach/hurd/libc-lock.h: ...here. (_BITS_LIBC_LOCK_H): Rename macro to _LIBC_LOCK_H. [_LIBC]: Include <libc-lockP.h> instead of <bits/libc-lockP.h>. * sysdeps/mach/bits/libc-lock.h: Move to ... * sysdeps/mach/libc-lock.h: ...here. (_BITS_LIBC_LOCK_H): Rename macro to _LIBC_LOCK_H. * sysdeps/nptl/bits/libc-lock.h: Move to ... * sysdeps/nptl/libc-lock.h: ...here. (_BITS_LIBC_LOCK_H): Rename macro to _LIBC_LOCK_H. * sysdeps/nptl/bits/libc-lockP.h: Move to ... * sysdeps/nptl/libc-lockP.h: ...here. (_BITS_LIBC_LOCKP_H): Rename macro to _LIBC_LOCKP_H. * crypt/crypt_util.c: Include <libc-lock.h> instead of <bits/libc-lock.h>. * dirent/scandir-tail.c: Likewise. * dlfcn/dlerror.c: Likewise. * elf/dl-close.c: Likewise. * elf/dl-iteratephdr.c: Likewise. * elf/dl-lookup.c: Likewise. * elf/dl-open.c: Likewise. * elf/dl-support.c: Likewise. * elf/dl-writev.h: Likewise. * elf/rtld.c: Likewise. * grp/fgetgrent.c: Likewise. * gshadow/fgetsgent.c: Likewise. * gshadow/sgetsgent.c: Likewise. * iconv/gconv_conf.c: Likewise. * iconv/gconv_db.c: Likewise. * iconv/gconv_dl.c: Likewise. * iconv/gconv_int.h: Likewise. * iconv/gconv_trans.c: Likewise. * include/link.h: Likewise. * inet/getnameinfo.c: Likewise. * inet/getnetgrent.c: Likewise. * inet/getnetgrent_r.c: Likewise. * intl/bindtextdom.c: Likewise. * intl/dcigettext.c: Likewise. * intl/finddomain.c: Likewise. * intl/gettextP.h: Likewise. * intl/loadmsgcat.c: Likewise. * intl/localealias.c: Likewise. * intl/textdomain.c: Likewise. * libidn/idn-stub.c: Likewise. * libio/libioP.h: Likewise. * locale/duplocale.c: Likewise. * locale/freelocale.c: Likewise. * locale/newlocale.c: Likewise. * locale/setlocale.c: Likewise. * login/getutent_r.c: Likewise. * login/getutid_r.c: Likewise. * login/getutline_r.c: Likewise. * login/utmp-private.h: Likewise. * login/utmpname.c: Likewise. * malloc/mtrace.c: Likewise. * misc/efgcvt.c: Likewise. * misc/error.c: Likewise. * misc/fstab.c: Likewise. * misc/getpass.c: Likewise. * misc/mntent.c: Likewise. * misc/syslog.c: Likewise. * nis/nis_call.c: Likewise. * nis/nis_callback.c: Likewise. * nis/nss-default.c: Likewise. * nis/nss_compat/compat-grp.c: Likewise. * nis/nss_compat/compat-initgroups.c: Likewise. * nis/nss_compat/compat-pwd.c: Likewise. * nis/nss_compat/compat-spwd.c: Likewise. * nis/nss_nis/nis-alias.c: Likewise. * nis/nss_nis/nis-ethers.c: Likewise. * nis/nss_nis/nis-grp.c: Likewise. * nis/nss_nis/nis-hosts.c: Likewise. * nis/nss_nis/nis-network.c: Likewise. * nis/nss_nis/nis-proto.c: Likewise. * nis/nss_nis/nis-pwd.c: Likewise. * nis/nss_nis/nis-rpc.c: Likewise. * nis/nss_nis/nis-service.c: Likewise. * nis/nss_nis/nis-spwd.c: Likewise. * nis/nss_nisplus/nisplus-alias.c: Likewise. * nis/nss_nisplus/nisplus-ethers.c: Likewise. * nis/nss_nisplus/nisplus-grp.c: Likewise. * nis/nss_nisplus/nisplus-hosts.c: Likewise. * nis/nss_nisplus/nisplus-initgroups.c: Likewise. * nis/nss_nisplus/nisplus-network.c: Likewise. * nis/nss_nisplus/nisplus-proto.c: Likewise. * nis/nss_nisplus/nisplus-pwd.c: Likewise. * nis/nss_nisplus/nisplus-rpc.c: Likewise. * nis/nss_nisplus/nisplus-service.c: Likewise. * nis/nss_nisplus/nisplus-spwd.c: Likewise. * nis/ypclnt.c: Likewise. * nptl/libc_pthread_init.c: Likewise. * nss/getXXbyYY.c: Likewise. * nss/getXXent.c: Likewise. * nss/getXXent_r.c: Likewise. * nss/nss_db/db-XXX.c: Likewise. * nss/nss_db/db-netgrp.c: Likewise. * nss/nss_db/nss_db.h: Likewise. * nss/nss_files/files-XXX.c: Likewise. * nss/nss_files/files-alias.c: Likewise. * nss/nsswitch.c: Likewise. * posix/regex_internal.h: Likewise. * posix/wordexp.c: Likewise. * pwd/fgetpwent.c: Likewise. * resolv/res_hconf.c: Likewise. * resolv/res_libc.c: Likewise. * shadow/fgetspent.c: Likewise. * shadow/lckpwdf.c: Likewise. * shadow/sgetspent.c: Likewise. * socket/opensock.c: Likewise. * stdio-common/reg-modifier.c: Likewise. * stdio-common/reg-printf.c: Likewise. * stdio-common/reg-type.c: Likewise. * stdio-common/vfprintf.c: Likewise. * stdio-common/vfscanf.c: Likewise. * stdlib/abort.c: Likewise. * stdlib/cxa_atexit.c: Likewise. * stdlib/fmtmsg.c: Likewise. * stdlib/random.c: Likewise. * stdlib/setenv.c: Likewise. * string/strsignal.c: Likewise. * sunrpc/auth_none.c: Likewise. * sunrpc/bindrsvprt.c: Likewise. * sunrpc/create_xid.c: Likewise. * sunrpc/key_call.c: Likewise. * sunrpc/rpc_thread.c: Likewise. * sysdeps/arm/backtrace.c: Likewise. * sysdeps/generic/ldsodefs.h: Likewise. * sysdeps/generic/stdio-lock.h: Likewise. * sysdeps/generic/unwind-dw2-fde.c: Likewise. * sysdeps/i386/backtrace.c: Likewise. * sysdeps/ieee754/ldbl-opt/nldbl-compat.c: Likewise. * sysdeps/m68k/backtrace.c: Likewise. * sysdeps/mach/hurd/cthreads.c: Likewise. * sysdeps/mach/hurd/dirstream.h: Likewise. * sysdeps/mach/hurd/malloc-machine.h: Likewise. * sysdeps/nptl/malloc-machine.h: Likewise. * sysdeps/nptl/stdio-lock.h: Likewise. * sysdeps/posix/dirstream.h: Likewise. * sysdeps/posix/getaddrinfo.c: Likewise. * sysdeps/posix/system.c: Likewise. * sysdeps/pthread/aio_suspend.c: Likewise. * sysdeps/s390/s390-32/backtrace.c: Likewise. * sysdeps/s390/s390-64/backtrace.c: Likewise. * sysdeps/unix/sysv/linux/check_pf.c: Likewise. * sysdeps/unix/sysv/linux/if_index.c: Likewise. * sysdeps/unix/sysv/linux/s390/s390-32/getutent_r.c: Likewise. * sysdeps/unix/sysv/linux/s390/s390-32/getutid_r.c: Likewise. * sysdeps/unix/sysv/linux/s390/s390-32/getutline_r.c: Likewise. * sysdeps/unix/sysv/linux/shm-directory.c: Likewise. * sysdeps/unix/sysv/linux/system.c: Likewise. * sysdeps/x86_64/backtrace.c: Likewise. * time/alt_digit.c: Likewise. * time/era.c: Likewise. * time/tzset.c: Likewise. * wcsmbs/wcsmbsload.c: Likewise. * nptl/tst-initializers1.c (do_test): Refer to <libc-lock.h> instead of <bits/libc-lock.h> in comment.
* Disable uninitialized warning with GCC 4.8Stan Shebs2015-08-261-0/+12
| | | | | | | | | As with other spots in the code, GCC 4.8 unnecessarily complains about an uninitialized variable in tanl calcs, so this patch disables. With it, the library and sees the usual set of test passes. * sysdeps/ieee754/ldbl-96/k_tanl.c: Include <libc-internal.h>. (__kernel_tanl): Ignore uninitialized warnings around use of SIGN.
* 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 fma spurious underflows (bug 18824).Joseph Myers2015-08-143-3/+3
| | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | Various fma implementations have logic that, when computing fma (x, y, z) where z is large (so care needs taking to avoid internal overflow) but x * y is small, scale x * y up instead of down to avoid internal underflows resulting from scaling down. (In these cases, x * y is small enough that only its sign actually matters rather than the exact value.) The threshold for scaling up instead of down was correct for "if the unscaled values were multiplied, the low part of the multiplication could underflow", and the scaling was sufficient to ensure that the low part of the multiplication did not underflow (given that cases of very small x * y - less than half the least subnormal - were previously dealt with). However, the choice in the functions wasn't between scaling up or no scaling, but between scaling up and scaling down (scaling down actually being needed when x * y isn't so small compared to z and so the exact value does matter). Thus a larger threshold is needed to ensure that scaling down doesn't produce values the multiplication of whose low parts underflows. This patch increases the thresholds accordingly. Tested for x86_64, x86 and mips64 (with the MIPS version of s_fmal.c removed so that the ldbl-128 version gets tested instead of the soft-fp one). [BZ #18824] * sysdeps/ieee754/dbl-64/s_fma.c (__fma): Increase threshold for scaling x * y up instead of down. * sysdeps/ieee754/ldbl-128/s_fmal.c (__fmal): Likewise. * sysdeps/ieee754/ldbl-96/s_fmal.c (__fmal): Likewise. * math/auto-libm-test-in: Add more tests of fma. * math/auto-libm-test-out: Regenerated.
* Fix tanh missing underflows (bug 16520).Joseph Myers2015-08-135-2/+42
| | | | | | | | | | | | | | | | | | | | | | | | | | | | | 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-075-4/+51
| | | | | | | | | | | | | | | | | | | | | | | | | | | | 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-065-7/+43
| | | | | | | | | | | | | | | | | | | | | | | | | | | | | 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.
* Fix powf (close to -1, large) (bug 18647).Joseph Myers2015-08-051-1/+1
| | | | | | | | | | | | | | | | The flt-32 implementation of powf wrongly uses x-1 instead of |x|-1 when computing log (x) for the case where |x| is close to 1 and y is large. This patch fixes the logic accordingly. Relevant tests existed for x close to 1, and corresponding tests are added for x close to -1, as well as for some new variant cases. Tested for x86_64 and x86. [BZ #18647] * sysdeps/ieee754/flt-32/e_powf.c (__ieee754_powf): For large y and |x| close to 1, use absolute value of x when computing log. * math/auto-libm-test-in: Add more tests of pow. * math/auto-libm-test-out: Regenerated.
* Remove unused file sysdeps/ieee754/support.cWilco Dijkstra2015-07-061-520/+0
|
* Fix ldbl-128 expm1l (-min_subnorm) result sign (bug 18619).Joseph Myers2015-07-011-3/+12
| | | | | | | | | | | | | In the ldbl-128 implementation of expm1l, when expm1l's result should underflow to 0 (argument minus the least subnormal, in some rounding modes), it can be a zero of the wrong sign. This patch fixes this in the same way previously used for the x86 / x86_64 versions. Tested for mips64. [BZ #18619] * sysdeps/ieee754/ldbl-128/s_expm1l.c (__expm1l): Force underflow and return argument in case of subnormal argument.
* Improve tgamma accuracy (bug 18613).Joseph Myers2015-06-295-126/+271
| | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | 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 ldbl-128 j1l spurious underflows (bug 18612).Joseph Myers2015-06-291-0/+10
| | | | | | | | | | | | | | | | | The ldbl-128 implementation of j1l produces spurious underflow exceptions for some small arguments, as a result of squaring the argument. This patch fixes it just to use a linear approximation for sufficiently small arguments, and then to force an underflow exception only in the cases where it is required. Tested for mips64. [BZ #18612] * sysdeps/ieee754/ldbl-128/e_j1l.c (__ieee754_j1l): For small arguments, just return 0.5 times the argument, with underflow forced as needed. * math/auto-libm-test-in: Add more tests of j1. * math/auto-libm-test-out: Regenerated.
* Fix j1, jn missing underflows (bug 16559).Joseph Myers2015-06-298-5/+55
| | | | | | | | | | | | | | | | | | | | | | | | | | 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-255-710/+741
| | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | 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 ldbl-128 expl missing underflows (bug 18586).Joseph Myers2015-06-241-1/+9
| | | | | | | | | | | | | | | Similar to various other bugs in this area, the ldbl-128 expl implementation does not raise the underflow exception for all subnormal results, if the scaling down is exact although the actual result is inexact. This patch fixes this by forcing the exception in this case (the tests that failed before and pass after the test are already in the testsuite). Tested for mips64. [BZ #18586] * sysdeps/ieee754/ldbl-128/e_expl.c (__ieee754_expl): Force underflow exception for small results.
* Fix sin, sincos missing underflows (bug 16526, bug 16538).Joseph Myers2015-06-237-17/+74
| | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | 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 spurious "inexact" exceptions from __kernel_standard_l (bug 18245, bug ↵Joseph Myers2015-06-231-24/+8
| | | | | | | | | | | | | | | | | | | | | | | | | | | | | 18583). __kernel_standard_l converts long double arguments to double for use in SVID "struct exception". This has special-case handling for when that conversion would overflow or underflow but the original long double function wouldn't. However, it turns out that "inexact" exceptions can be spurious here as well, when the function is exactly determined and __kernel_standard_l is being called for a domain error. This patch fixes this by using feholdexcept / fesetenv to avoid exceptions from the conversion, replacing the previous special-case logic for overflow and underflow (this covers all functions using __kernel_standard_l, not just those that actually need a change, since there doesn't seem to be much point in restricting things just to the functions that mustn't get "inexact" here). Tested for x86_64 and x86. [BZ #18245] [BZ #18583] * sysdeps/ieee754/k_standardl.c: Include <fenv.h>. (__kernel_standard_l): Use feholdexcept and fesetenv around conversion to double instead of special-casing overflow and underflow. * math/libm-test.inc (fmod_test_data): Add more tests. (remainder_test_data): Likewise. (sqrt_test_data): Likewise.
* Fix exp2, exp2f spurious underflows (bug 18219).Joseph Myers2015-06-232-2/+6
| | | | | | | | | | | | | | | | | | | | | The dbl-64 and flt-32 implementations of exp2 functions produce spurious underflow exceptions. The underlying reason is the same in both cases: the computation works as (2^a - 1)*2^b + 2^b for suitably chosen a and b, where a has small magnitude so 2^a - 1 can be computed with a low-degree polynomial approximation, and (2^a - 1)*2^b can underflow even when the final result does not. This patch fixes this by adjusting the threshold for when scaling is used to avoid intermediate underflow so it works for any possible value of a where the final result would not underflow. Tested for x86_64 and x86. [BZ #18219] * sysdeps/ieee754/dbl-64/e_exp2.c (__ieee754_exp2): Reduce threshold on absolute value of exponent for which scaling is used. * sysdeps/ieee754/flt-32/e_exp2f.c (__ieee754_exp2f): Likewise. * math/auto-libm-test-in: Add more tests of exp2. * math/auto-libm-test-out: Regenerated.
* Fix expm1 missing underflows (bug 16353).Joseph Myers2015-06-223-1/+14
| | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | Similar to various other bugs in this area, some expm1 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. (The issue does not apply to the ldbl-* implementations or to those for x86 / x86_64 long double. The change to sysdeps/ieee754/dbl-64/wordsize-64/e_cosh.c is one I missed when previously fixing bug 16354; the bug in that implementation was previously latent, but the expm1 fixes stopped it being latent and so required it to be fixed to avoid spurious underflows from cosh.) Tested for x86_64 and x86. [BZ #16353] * sysdeps/i386/fpu/s_expm1.S (dbl_min): New object. (__expm1): Force underflow exception for arguments with small absolute value. * sysdeps/i386/fpu/s_expm1f.S (flt_min): New object. (__expm1f): Force underflow exception for arguments with small absolute value. * sysdeps/ieee754/dbl-64/s_expm1.c: Include <float.h>. (__expm1): Force underflow exception for arguments with small absolute value. * sysdeps/ieee754/flt-32/s_expm1f.c: Include <float.h>. (__expm1f): Force underflow exception for arguments with small absolute value. * sysdeps/ieee754/dbl-64/wordsize-64/e_cosh.c (__ieee754_cosh): Check for small arguments before calling __expm1. * math/auto-libm-test-in: Do not mark underflow exceptions as possibly missing for bug 16353. * math/auto-libm-test-out: Regenerated.
* Fixed powerpc64 build.Andrew Senkevich2015-06-191-0/+4
| | | | | | * sysdeps/ieee754/ldbl-opt/s_sin.c (__DECL_SIMD_sincos_disable, __DECL_SIMD_sincos_disablef, __DECL_SIMD_sincos_disablel): Added empty definitions for proper unfolding of __MATHDECL_VEC.
* Fix asinh missing underflows (bug 16350).Joseph Myers2015-06-185-0/+30
| | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | 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-033-3/+3
|
* This patch renames all uses of __isinf*, __isnan*, __finite* and __signbit* ↵Wilco Dijkstra2015-06-0334-53/+53
| | | | to use standard C99 macros. This has no effect on generated code.
* 2015-05-28 Wilco Dijkstra <wdijkstr@arm.com>Wilco Dijkstra2015-05-282-10/+2
| | | | | * sysdeps/ieee754/dbl-64/s_fabs.c: (__fabs): Call __builtin_fabs. * sysdeps/ieee754/flt-32/s_fabsf.c: (__fabsf): Likewise.
* Fix ldbl-128 / ldbl-128ibm tanl for -Wuninitialized.Joseph Myers2015-05-222-0/+24
| | | | | | | | | | | | | | | 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-222-2/+2
| | | | | | | | | | | | | | | | | | | 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-222-4/+6
| | | | | | | | | | | | | | | | | | | | | | | | | | | 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 lgamma implementations for -Wuninitialized.Joseph Myers2015-05-213-0/+35
| | | | | | | | | | | | | | | | | | | | | | | | | | | | | | If you remove the "override CFLAGS += -Wno-uninitialized" in math/Makefile, you get errors from lgamma implementations of the form: ../sysdeps/ieee754/dbl-64/e_lgamma_r.c: In function '__ieee754_lgamma_r': ../sysdeps/ieee754/dbl-64/e_lgamma_r.c:297:13: error: 'nadj' may be used uninitialized in this function [-Werror=maybe-uninitialized] if(hx<0) r = nadj - r; This is one of the standard kinds of false positive uninitialized warnings: nadj is set under a certain condition, and then later used under the same condition. This patch uses DIAG_* macros to suppress the warning on the use of nadj. The ldbl-128 / ldbl-128ibm implementation has a substantially different structure that avoids this issue. Tested for x86_64. (In fact this patch eliminates the need for that -Wno-uninitialized on x86_64, but I want to test on more architectures before removing it.) * sysdeps/ieee754/dbl-64/e_lgamma_r.c: Include <libc-internal.h>. (__ieee754_lgamma_r): Ignore uninitialized warnings around use of NADJ. * sysdeps/ieee754/flt-32/e_lgammaf_r.c: Include <libc-internal.h>. (__ieee754_lgammaf_r): Ignore uninitialized warnings around use of NADJ. * sysdeps/ieee754/ldbl-96/e_lgammal_r.c: Include <libc-internal.h>. (__ieee754_lgammal_r): Ignore uninitialized warnings around use of NADJ.
* Fix sysdeps/ieee754/dbl-64/mpa.c for -Wuninitialized.Joseph Myers2015-05-211-2/+3
| | | | | | | | | | | | | | | | | | | | If you remove the "override CFLAGS += -Wno-uninitialized" in math/Makefile, one of the errors you get is: ../sysdeps/ieee754/dbl-64/mpa.c: In function '__mp_dbl.part.0': ../sysdeps/ieee754/dbl-64/mpa.c:183:5: error: 'c' may be used uninitialized in this function [-Werror=maybe-uninitialized] c *= X[0]; The problem is that the p < 5 case initializes c if p is 1, 2, 3 or 4 but not otherwise, and in fact p is positive for all calls to this function so the uninitialized case can't actually occur. This patch replaces the "if (p == 4)" last case with a comment so the compiler can see that all paths do initialize c. Tested for x86_64. * sysdeps/ieee754/dbl-64/mpa.c (norm): Remove if condition on (p == 4) case.
* Fix ldbl-96 remquol (finite, Inf) (bug 18244).Joseph Myers2015-05-191-1/+1
| | | | | | | | | | | | | | | | | | ldbl-96 remquol wrongly handles the case where the first argument is finite and the second infinite, because the check for the second argument being a NaN fails to disregard the explicit high mantissa bit and so wrongly interprets an infinity as being a NaN. This patch fixes this by masking off that bit, and improves test coverage for both remainder and remquo (various cases were missing tests, or, as in the case of the bug, were tested only for one of the two functions). Tested for x86_64 and x86. [BZ #18244] * sysdeps/ieee754/ldbl-96/s_remquol.c (__remquol): Ignore explicit high mantissa bit when testing whether P is a NaN. * math/libm-test.inc (remainder_test_data): Add more tests. (remquo_test_data): Likewise.
* Fix atanhl missing underflows (bug 16352).Joseph Myers2015-05-155-2/+38
| | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | 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.
* Fix tanf spurious underflows (bug 18221).Joseph Myers2015-05-151-1/+1
| | | | | | | | | | | | | | | | The flt-32 implementation of tanf produces spurious underflow exceptions for some small arguments, through computing values on the order of x^5. This patch fixes this by adjusting the threshold for returning x (or, as applicable, +/- 1/x) to 2**-13 (the next term in the power series being x^3/3). Tested for x86_64 and x86. [BZ #18221] * sysdeps/ieee754/flt-32/k_tanf.c (__kernel_tanf): Use 2**-13 not 2**-28 as threshold for returning x or +/- 1/x. * math/auto-libm-test-in: Add more tests of tan. * math/auto-libm-test-out: Regenerated.
* Fix lgammaf spurious underflows (bug 18220).Joseph Myers2015-05-151-3/+3
| | | | | | | | | | | | | | | | | The flt-32 implementation of lgammaf produces spurious underflow exceptions for some large arguments, because of calculations involving x^-2 multiplied by small constants. This patch fixes this by adjusting the threshold for a simpler computation to 2**26 (the error in the simpler computation is on the order of 0.5 * log (x), for a result on the order of x * log (x)). Tested for x86_64 and x86. [BZ #18220] * sysdeps/ieee754/flt-32/e_lgammaf_r.c (__ieee754_lgammaf_r): Use 2**26 not 2**58 as threshold for returning x * (log (x) - 1). * math/auto-libm-test-in: Add another test of lgamma. * math/auto-libm-test-out: Regenerated.
* Remove various ABS macros and replace uses with fabs (or in one case abs)Wilco Dijkstra2015-05-1511-85/+90
| | | | which is more efficient on all targets.
* Fix erfcf spurious underflows (bug 18217).Joseph Myers2015-05-151-1/+1
| | | | | | | | | | | | | | | | | | The flt-32 implementation of erfcf produces spurious underflow exceptions for some arguments close to 0, because of calculations squaring the argument and then multiplying by small constants. This patch fixes this by adjusting the threshold for arguments for which the result is so close to 1 that 1 - x will give the right result from 2**-56 to 2**-26. (If 1 - x * 2/sqrt(pi) were used, the errors would be on the order of x^3 and a much larger threshold could be used.) Tested for x86_64 and x86. [BZ #18217] * sysdeps/ieee754/flt-32/s_erff.c (__erfcf): Use 2**-26 not 2**-56 as threshold for returning 1 - x. * math/auto-libm-test-in: Add more tests of erfc. * math/auto-libm-test-out: Regenerated.
* Fix atanf spurious underflows (bug 18196).Joseph Myers2015-05-141-1/+1
| | | | | | | | | | | | | | | | The sysdeps/ieee754/flt-32 version of atanf produces spurious underflow exceptions for some large arguments, because of computations that compute x^-4. This patch fixes this by adjusting the threshold for large arguments (for which +/- pi/2 can just be returned, the correct result being roughly +/- pi/2 - 1/x) from 2^34 to 2^25. Tested for x86_64 and x86. [BZ #18196] * sysdeps/ieee754/flt-32/s_atanf.c (__atanf): Use 2^25 not 2^34 as threshold for large arguments. * math/auto-libm-test-in: Add another test of atan. * math/auto-libm-test-out: Regenerated.
* Fix log1p missing underflows (bug 16339).Joseph Myers2015-05-143-1/+23
| | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | | Similar to various other bugs in this area, some log1p 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. (The ldbl-128ibm implementation doesn't currently need any change as it already generates this exception, albeit through code that would generate spurious exceptions in other cases; special code for this issue will only be needed there when fixing the spurious exceptions.) Tested for x86_64, x86, powerpc and mips64. [BZ #16339] * sysdeps/i386/fpu/s_log1p.S (dbl_min): New object. (__log1p): Force underflow exception for results with small absolute value. * sysdeps/i386/fpu/s_log1pf.S (flt_min): New object. (__log1pf): Force underflow exception for results with small absolute value. * sysdeps/ieee754/dbl-64/s_log1p.c: Include <float.h>. (__log1p): Force underflow exception for results with small absolute value. * sysdeps/ieee754/flt-32/s_log1pf.c: Include <float.h>. (__log1pf): Force underflow exception for results with small absolute value. * sysdeps/ieee754/ldbl-128/s_log1pl.c: Include <float.h>. (__log1pl): Force underflow exception for results with small absolute value. * math/auto-libm-test-in: Do not allow missing underflow exceptions from log1p. * math/auto-libm-test-out: Regenerated.
* Fix ldbl-128 roundl for exponents in [31, 47] (bug 18346).Joseph Myers2015-04-281-1/+1
| | | | | | | | | | | | | | | | | | | | | | | | | | | The implementation of roundl for ldbl-128 involves undefined behavior for arguments with exponents from 31 to 47 inclusive, from the shift: u_int64_t i = -1ULL >> (j0 - 48); For example, on mips64, this means roundl (0xffffffffffff.8p0L) wrongly returns its argument, which is not an integer. A condition checking for exponents < 31 should actually be checking for exponents < 48, and this patch makes it do so. (That condition is for whether the bit representing 0.5 is in the high 64-bit half of the floating-point number. The value 31 might have arisen from an incorrect conversion of the ldbl-96 version to handle ldbl-128.) This was originally reported as a GCC libquadmath bug <https://gcc.gnu.org/bugzilla/show_bug.cgi?id=65757>. Tested for mips64; also tested for x86_64 and x86 to make sure the new tests pass there. [BZ #18346] * sysdeps/ieee754/ldbl-128/s_roundl.c (__roundl): Handle all exponents less than 48 as cases where high part of mantissa needs examining to determine whether argument is integral. * math/libm-test.inc (round_test_data): Add more tests.
* Use __copysign rather than copysign.Wilco Dijkstra2015-04-222-2/+2
|
generated by cgit-pink 1.4.1 (git 2.36.1) at 2024-11-07 18:18:49 +0000