aboutsummaryrefslogtreecommitdiff
path: root/sysdeps/ieee754
AgeCommit message (Collapse)Author
2015-03-10powerpc: Fix incorrect results for pow when using FMAAdhemerval Zanella
This patch adds no FMA generation for e_pow to avoid precision issues for powerpc. This fixes BZ#18104.
2015-02-27Fix ldbl-96, ldbl-128ibm atanhl inaccuracy (bug 18046, bug 18047).Joseph Myers
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.
2015-02-27Add comment to CSTR macro in k_standard.c.Joseph Myers
* sysdeps/ieee754/k_standard.c (CSTR): Add comment.
2015-02-26Avoid -Wno-write-strings for k_standard.c.Joseph Myers
We want to avoid -Wno- options in makefiles as far as possible, by cleaning up the underlying issues if possible or failing that by using diagnostic pragmas. This patch eliminates the use of -Wno-write-strings for sysdeps/ieee754/k_standard.c by using casts in the source file to cast away const; those casts are encapsulated in a macro that also deals with the choice of strings for float / double / long double functions (for which the logic was previously replicated many times). Tested for x86_64; the only change to disassembly of installed stripped shared libraries was a line number in an assertion. * sysdeps/ieee754/k_standard.c (CSTR): New macro. (__kernel_standard): Use CSTR macro when setting exc.name. * sysdeps/ieee754/Makefile [$(subdir) = math] (CFLAGS-k_standard.c): Remove variable.
2015-02-26Avoid uninitialized warnings in Bessel functions.Joseph Myers
math/Makefile currently has: # The fdlibm code generates a lot of these warnings but is otherwise clean. override CFLAGS += -Wno-uninitialized This is of course undesirable; warnings should be disabled as narrowly as possible. To remove this override, we need to fix files that generate such warnings, or put warning-disabling pragmas in them. This patch does so for Bessel function implementations, one of the cases that have the warnings if the override is removed. The warnings arise because functions set pointer variables p and q only for certain values of the function argument, then use them unconditionally. As the static functions in question only get called for arguments that satisfy the last condition in the if/else chain, the natural fix is to change the last "else if" to just "else", which this patch does. (The ldbl-128 / ldbl-128ibm implementation of these functions is substantially different and looks like it already does use "else" in the last case in the nearest corresponding code.) Tested for x86_64 and x86. * sysdeps/ieee754/dbl-64/e_j0.c (pzero): Change last case for setting p and q from "else if" to "else". (qzero): Likewise. * sysdeps/ieee754/dbl-64/e_j1.c (pone): Likewise. (qone): Likewise. * sysdeps/ieee754/flt-32/e_j0f.c (pzerof): Likewise. (qzerof): Likewise. * sysdeps/ieee754/flt-32/e_j1f.c (ponef): Likewise. (qonef): Likewise. * sysdeps/ieee754/ldbl-96/e_j0l.c (pzero): Likewise. (qzero): Likewise. * sysdeps/ieee754/ldbl-96/e_j1l.c (pone): Likewise. (qone): Likewise.
2015-02-26Fix ldbl-128/ldbl-128ibm acosl inaccuracy (bug 18038, bug 18039).Joseph Myers
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.
2015-02-26Fix asin missing underflows (bug 16351).Joseph Myers
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.
2015-02-26Fix ldbl-128ibm logbl near powers of 2 (bug 18030).Joseph Myers
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.
2015-02-26Fix ldbl-128ibm ilogbl near powers of 2 (bug 18029).Joseph Myers
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.
2015-02-25Fix ldbl-128ibm asinhl inaccuracy (bug 18020).Joseph Myers
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.
2015-02-25Fix ldbl-128ibm acoshl inaccuracy (bug 18019).Joseph Myers
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.
2015-02-18Fix atan / atan2 missing underflows (bug 15319).Joseph Myers
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.
2015-02-17Fix sign of remquo zero remainder in round-downward mode (bug 17987).Joseph Myers
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.
2015-02-16Fix remquo spurious overflows (bug 17978).Joseph Myers
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.
2015-02-13Fix dbl-64/wordsize-64 remquo (bug 17569).Joseph Myers
The dbl-64/wordsize-64 remquo implementation follows similar logic to various other implementations, but where that logic computes some absolute values, it wrongly uses a previously computed bit-pattern for the absolute value of the first argument, where actually it needs the absolute value of the first argument mod 8 times the second. This patch fixes it to compute the correct absolute value. The integer quotient result of remquo is only specified mod 8 (including its sign); architecture-specific versions may well vary in what results they give for higher bits of that result (and indeed bug 17569 gives an example correct result from __builtin_remquo giving 9 for that result, where the particular glibc implementation used in that bug report would give 1 after this fix). Thus, this patch adapts the tests of remquo to test that result only mod 8, to allow for such variation when tests with higher quotient are included. Tested for x86_64 and x86. [BZ #17569] * sysdeps/ieee754/dbl-64/wordsize-64/s_remquo.c (__remquo): Compute absolute value of x as modified by fmod, not original value of x. * math/libm-test.inc (RUN_TEST_ffI_f1): Rename to RUN_TEST_ffI_f1_mod8. Check extra return value mod 8. (RUN_TEST_LOOP_ffI_f1): Rename to RUN_TEST_LOOP_ffI_f1_mod8. Call RUN_TEST_ffI_f1_mod8. (remquo_test_data): Add more tests.
2015-02-12Fix exp2 spurious underflows (bug 16560).Joseph Myers
This patch fixes the remaining part of bug 16560, spurious underflows from exp2 of arguments close to 0 (when the result is close to 1, so should not underflow), by just using 1+x instead of a more complicated calculation when the argument is sufficiently small. Tested for x86_64, x86 and mips64. [BZ #16560] * math/e_exp2l.c [LDBL_MANT_DIG == 106] (LDBL_EPSILON): Undefine and redefine. (__ieee754_exp2l): Do not multiply small fractional parts by M_LN2l. * sysdeps/i386/fpu/e_exp2l.S (__ieee754_exp2l): Just add 1 to small argument. * 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 (__ieee754_exp2l): Likewise. * math/auto-libm-test-in: Add more tests of exp2. * math/auto-libm-test-out: Regenerated.
2015-02-11Fix sincos errno setting (bug 15467).Joseph Myers
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.
2015-01-12Fix ldbl-96 scalblnl underflowing results (bug 17803).Joseph Myers
The ldbl-96 implementation of scalblnl (used for x86_64 and ia64) uses a condition k <= -63 to determine when a standard underflowing result tiny*__copysignl(tiny,x) should be returned. However, that condition corresponds to values with exponent -16446 or less, and in the case of -16446, the correct result for round-to-nearest depends on whether the value is exactly 0x1p-16446 (half the least subnormal) or more than that. This patch fixes the bug by changing the condition to k <= -64 and accordingly adjusting the exponent by 64 not 63 when converting to a normal value. Tested for x86_64. [BZ #17803] * sysdeps/ieee754/ldbl-96/s_scalblnl.c (twom63): Rename to twom64. Adjust value to 0x1p-64L. (__scalblnl): Only return standard underflowing result for K <= -64 not K <= -63; adjust exponent for underflowing result by 64 not 63. * math/libm-test.inc (scalbn_test_data): Add more tests. (scalbln_test_data): Likewise.
2015-01-12Fix ldbl-96 scalblnl for subnormal arguments (bug 17834).Joseph Myers
The ldbl-96 implementation of scalblnl (used for x86_64 and ia64) is incorrect for subnormal arguments (this is a separate bug from bug 17803, which is about underflowing results). There are two problems with the adjustments of subnormal arguments: the "two63" variable multiplied by is actually 0x1p52L not 0x1p63L, so is insufficient to make values normal, and then GET_LDOUBLE_EXP(es,x), used to extract the new exponent, extracts it into a variable that isn't used, while the value taken to by the new exponent is wrongly taken from the high part of the mantissa before the adjustment (hx). This patch fixes both those problems and adds appropriate tests. Tested for x86_64. [BZ #17834] * sysdeps/ieee754/ldbl-96/s_scalblnl.c (two63): Change value to 0x1p63L. (__scalblnl): Get new exponent of adjusted subnormal value from ES not HX. * math/libm-test.inc (scalbn_test_data): Add more tests. (scalbln_test_data): Likewise.
2015-01-05lround: provide cast for wordsize-64 version if neededChris Metcalf
Platforms with 64-bit registers where 32-bit values need to have the high 32 bits set in a particular way need to have an explicit cast when using the 64-bit sysdeps/ieee754/dbl-64/wordsize-64 version of llround() as lround(). This includes tilegx32, and likely MIPS. x32 does not need this, and AArch64 ILP32 will not either. Require it to be specified in sysdep.h to be explicit.
2015-01-02Fix libm fegetround namespace (bug 17748).Joseph Myers
Continuing the fixes for C90 libm functions calling C99 fe* functions, this patch fixes the case of fegetround by making it a weak alias of __fegetround and making the affected code call __fegetround. Tested for x86_64 (testsuite, and that disassembly of installed shared libraries is unchanged by the patch). Also tested for ARM (soft-float) that fegetround failures disappear from the linknamespace test failures (feholdexcept, fesetenv, fesetround and feupdateenv remain to be addressed before bug 17748 is fully fixed, although this patch may suffice to fix the failures in some cases, when the libc_fe* functions are implemented but there is no architecture-specific sqrt implementation in use so there were failures from fegetround used by sqrt but no other such failures). [BZ #17748] * include/fenv.h (__fegetround): Declare. Use libm_hidden_proto. * math/fegetround.c (fegetround): Rename to __fegetround and define as weak alias of __fegetround. Use libm_hidden_weak. * sysdeps/aarch64/fpu/fegetround.c (fegetround): Likewise. * sysdeps/alpha/fpu/fegetround.c (fegetround): Likewise. * sysdeps/arm/fegetround.c (fegetround): Likewise. * sysdeps/hppa/fpu/fegetround.c (fegetround): Likewise. * sysdeps/i386/fpu/fegetround.c (fegetround): Likewise. * sysdeps/ia64/fpu/fegetround.c (fegetround): Likewise. * sysdeps/m68k/fpu/fegetround.c (fegetround): Likewise. * sysdeps/mips/fpu/fegetround.c (fegetround): Likewise. * sysdeps/powerpc/fpu/fegetround.c (fegetround): Likewise. Undefine after rather than before function definition; use parentheses around function name in definition. (__fegetround): Also undefine macro after function definition. * sysdeps/powerpc/nofpu/fegetround.c (fegetround): Rename to __fegetround and define as weak alias of __fegetround. Use libm_hidden_weak. Do not undefine as macro. * sysdeps/powerpc/powerpc32/e500/nofpu/fegetround.c (fegetround): Likewise. * sysdeps/s390/fpu/fegetround.c (fegetround): Rename to __fegetround and define as weak alias of __fegetround. Use libm_hidden_weak. * sysdeps/sh/sh4/fpu/fegetround.c (fegetround): Likewise. * sysdeps/sparc/fpu/fegetround.c (fegetround): Likewise. * sysdeps/tile/math_private.h (__fegetround): New inline function. * sysdeps/x86_64/fpu/fegetround.c (fegetround): Rename to __fegetround and define as weak alias of __fegetround. Use libm_hidden_weak. * sysdeps/ieee754/dbl-64/e_sqrt.c (__ieee754_sqrt): Use __fegetround instead of fegetround.
2015-01-02Update copyright dates with scripts/update-copyrights.Joseph Myers
2014-12-31Clean up powerpc fegetround / __fegetround inlines.Joseph Myers
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.
2014-12-22Split __kernel_standard* functions (fixes bug 17724).Joseph Myers
Bug 17724 reports references to fesetround being brought in by ldbl-128ibm rintl via references to __rintl from __kernel_standard_l. Because all three __kernel_standard* functions are in the same file, this gets brought in even though only the long double version __kernel_standard_l needs __rintl, and the C90 functions use only __kernel_standard. This patch fixes this by splitting the three versions into separate files; it's fine for long double functions to refer to fe* functions directly, unless they get called by C90 double functions. Tested for x86_64 (testsuite; the reordering of code means disassembly of shared libraries can't usefully be compared). Tested for powerpc that the relevant issue disappears from the linknamespace test output. [BZ #17724] * sysdeps/ieee754/k_standard.c: Don't include <float.h>. (__kernel_standard_f): Remove. Moved to k_standardf.c. (__kernel_standard_l): Remove. Moved to k_standardl.c with (char *) casts added. * sysdeps/ieee754/k_standardf.c: New file. * sysdeps/ieee754/k_standardl.c: Likewise. * math/Makefile (libm-support): Remove k_standard. (libm-calls): Add k_standard.
2014-12-11powerpc: Fix lgammal_r overflow warningsAdhemerval Zanella
ldbl-128ibm uses ldbl-128 e_lgammal_r implementation as is, however some constants definitions overflows for IBM long double range. This patch suppress the compiler warnings until the ldbl-128ibm implementation is fixed.
2014-11-26FIx ldbl-128ibm frexpl for 32-bit systems (bug 16619, bug 16740).Joseph Myers
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.
2014-11-24Remove IS_IN_libmSiddhesh Poyarekar
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.
2014-11-18Fix libm mpone, mptwo namespace (bug 17616).Joseph Myers
libm uses symbols mpone and mptwo for internal purposes. This patch moves them to the implementation namespace (__mpone and __mptwo). Tested for x86_64 (testsuite, and that installed stripped shared libraries are unchanged by the patch). [BZ #17616] * sysdeps/ieee754/dbl-64/mpa.c (mpone): Rename to __mpone. (mptwo): Rename to __mptwo. (__inv): Use __mptwo instead of mptwo. * sysdeps/ieee754/dbl-64/mpa.h (mpone): Rename to __mpone. (mptwo): Rename to __mptwo. * sysdeps/ieee754/dbl-64/mpatan.c (__mpatan): Use __mpone instead of mpone and __mptwo instead of mptwo. * sysdeps/ieee754/dbl-64/mpatan2.c (__mpatan2): Use __mpone instead of mpone. * sysdeps/ieee754/dbl-64/mpexp.c (__mpexp): Likewise. * sysdeps/ieee754/dbl-64/mplog.c (__mplog): Likewise. * sysdeps/ieee754/dbl-64/sincos32.c (__c32): Use __mpone instead of mpone and __mptwo instead of mptwo. (__mpranred): Use __mpone instead of mpone. * conform/Makefile (test-xfail-ISO/math.h/linknamespace): Remove variable. (test-xfail-ISO99/complex.h/linknamespace): Likewise. (test-xfail-ISO99/math.h/linknamespace): Likewise. (test-xfail-ISO99/tgmath.h/linknamespace): Likewise. (test-xfail-ISO11/complex.h/linknamespace): Likewise. (test-xfail-ISO11/math.h/linknamespace): Likewise. (test-xfail-ISO11/tgmath.h/linknamespace): Likewise. (test-xfail-XPG3/math.h/linknamespace): Likewise. (test-xfail-XPG4/math.h/linknamespace): Likewise. (test-xfail-POSIX/math.h/linknamespace): Likewise. (test-xfail-UNIX98/math.h/linknamespace): Likewise. (test-xfail-XOPEN2K/complex.h/linknamespace): Likewise. (test-xfail-XOPEN2K/math.h/linknamespace): Likewise. (test-xfail-XOPEN2K/tgmath.h/linknamespace): Likewise. (test-xfail-POSIX2008/complex.h/linknamespace): Likewise. (test-xfail-POSIX2008/math.h/linknamespace): Likewise. (test-xfail-POSIX2008/tgmath.h/linknamespace): Likewise. (test-xfail-XOPEN2K8/complex.h/linknamespace): Likewise. (test-xfail-XOPEN2K8/math.h/linknamespace): Likewise. (test-xfail-XOPEN2K8/tgmath.h/linknamespace): Likewise.
2014-08-04Fix missing <math_private.h> in ldbl-96 fmaAndreas Schwab
2014-08-01Force eval for fma implementationsRichard Henderson
2014-06-30Fix ldbl-128 expm1l spurious underflow (bug 16539).Joseph Myers
This patch fixes spurious underflows from ldbl-128 expm1l, as reported in <https://sourceware.org/ml/libc-alpha/2014-06/msg00835.html> and exposed by the tests added for such a bug in the x86 / x86-64 version. The bug and fix are essentially the same, so no separate bug is filed in Bugzilla. Tested for mips64. [BZ #16539] * sysdeps/ieee754/ldbl-128/s_expm1l.c: Include <float.h>. (__expm1l): Return argument unchanged when small but not subnormal.
2014-06-29Fix ldbl-128 powl sign of result in overflow / underflow cases (bug 17097).Joseph Myers
This patch fixes bug 17097, ldbl-128 powl producing overflowing / underflowing results with positive sign when the result should have been negative. This was shown up by the tests in non-default rounding modes added by my patch for bug 16315, but isn't actually limited to non-default rounding modes: rather, when rounding to nearest the wrappers produced a result with the correct sign and so always hid the bug unless -lieee was used to disable the wrappers. The problem is that in the cases where Y is large enough that the result overflows or underflows for X not very close to 1, but not large enough to overflow or underflow for all X != +/- 1 (in the latter case Y is always an even integer), a positive overflowing / underflowing result is always returned, rather than one with the correct sign. This patch moves the relevant part of computation of the sign earlier and returns a result of the correct sign. Tested for mips64. [BZ #17097] * sysdeps/ieee754/ldbl-128/e_powl.c (__ieee754_powl): Return result with correct sign in case of exponents that produce overflow except for X very close to 1.
2014-06-27Fix yn overflow handling in non-default rounding modes (bug 16561, bug 16562).Joseph Myers
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.
2014-06-25Fix exp10 spurious underflows (bug 16560).Joseph Myers
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.
2014-06-24Fix ldbl-128 erfl spurious underflows (bug 16287).Joseph Myers
This patch fixes bug 16287, spurious underflows from ldbl-128 erfl arising from it calling erfcl for arguments with absolute value at least 1.0, although for large positive arguments erfcl correctly underflows but erfl shouldn't. The fix is simply to avoid calling erfcl, and just return 1, for arguments above a cut-off large enough that erfl correctly rounds to-nearest as 1 but not so large that erfcl underflows. Tested mips64. Also tested x86_64 and x86 to confirm the new tests (taken from the tests of erfc) don't cause any problems there; no ulps updates needed. [BZ #16287] * sysdeps/ieee754/ldbl-128/s_erfl.c (__erfl): Return 1 without calling __erfcl for arguments at least 16. * math/auto-libm-test-in: Add more tests of erf. * math/auto-libm-test-out: Regenerated.
2014-06-23Fix cosh spurious underflows from expm1 (bug 16354), inaccurate results near ↵Joseph Myers
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.
2014-06-23Set errno for y1 overflow (bug 17050).Joseph Myers
This patch fixes bug 17050, missing errno setting for y1 overflow (for small positive arguments). An appropriate check is added for overflow directly in the __ieee754_y1 implementation, similar to the check present for yn (doing it there rather than in the wrapper also avoids yn needing to repeat the check when called for order 1 or -1 and it uses __ieee754_y1). Tested x86_64 and x86; no ulps update needed. Also tested for mips64 to verify the ldbl-128 fix (the ldbl-128ibm code just #includes the ldbl-128 file). [BZ #17050] * sysdeps/ieee754/dbl-64/e_j1.c: Include <errno.h>. (__ieee754_y1): Set errno if return value overflows. * sysdeps/ieee754/flt-32/e_j1f.c: Include <errno.h>. (__ieee754_y1f): Set errno if return value overflows. * sysdeps/ieee754/ldbl-128/e_j1l.c: Include <errno.h>. (__ieee754_y1l): Set errno if return value overflows. * sysdeps/ieee754/ldbl-96/e_j1l.c: Include <errno.h>. (__ieee754_y1l): Set errno if return value overflows. * math/auto-libm-test-in: Add more tests of y0, y1 and yn. * math/auto-libm-test-out: Regenerated.
2014-06-23Fix pow overflow in non-default rounding modes (bug 16315).Joseph Myers
This patch fixes bug 16315, bad pow handling of overflow/underflow in non-default rounding modes. Tests of pow are duly converted to ALL_RM_TEST to run all tests in all rounding modes. There are two main issues here. First, various implementations compute a negative result by negating a positive result, but this yields inappropriate overflow / underflow values for directed rounding, so either overflow / underflow results need recomputing in the correct sign, or the relevant overflowing / underflowing operation needs to be made to have a result of the correct sign. Second, the dbl-64 implementation sets FE_TONEAREST internally; in the overflow / underflow case, the result needs recomputing in the original rounding mode. Tested x86_64 and x86 and ulps updated accordingly. [BZ #16315] * sysdeps/i386/fpu/e_pow.S (__ieee754_pow): Ensure possibly overflowing or underflowing operations take place with sign of result. * sysdeps/i386/fpu/e_powf.S (__ieee754_powf): Likewise. * sysdeps/i386/fpu/e_powl.S (__ieee754_powl): Likewise. * sysdeps/ieee754/dbl-64/e_pow.c: Include <math.h>. (__ieee754_pow): Recompute overflowing and underflowing results in original rounding mode. * sysdeps/x86/fpu/powl_helper.c: Include <stdbool.h>. (__powl_helper): Allow negative argument X and scale negated value as needed. Avoid passing value outside [-1, 1] to f2xm1. * sysdeps/x86_64/fpu/e_powl.S (__ieee754_powl): Ensure possibly overflowing or underflowing operations take place with sign of result. * sysdeps/x86_64/fpu/multiarch/e_pow.c [HAVE_FMA4_SUPPORT]: Include <math.h>. * math/auto-libm-test-in: Add more tests of pow. * math/auto-libm-test-out: Regenerated. * math/libm-test.inc (pow_test): Use ALL_RM_TEST. (pow_tonearest_test_data): Remove. (pow_test_tonearest): Likewise. (pow_towardzero_test_data): Likewise. (pow_test_towardzero): Likewise. (pow_downward_test_data): Likewise. (pow_test_downward): Likewise. (pow_upward_test_data): Likewise. (pow_test_upward): Likewise. (main): Don't call removed functions. * sysdeps/i386/fpu/libm-test-ulps: Update. * sysdeps/x86_64/fpu/libm-test-ulps: Likewise.
2014-06-20[BZ #6803] Set errno for scalbln, scalbnStefan Liebler
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.
2014-06-17PowerPC: Fix nearbyintl failure for few inputsRajalakshmi Srinivasaraghavan
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.
2014-06-10Fix log2 (1) in round-downward mode (bug 17042).Joseph Myers
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.
2014-05-31Fix excessive ULP for y1_upward (0x2p+0) in test-float and test-ifloat.David S. Miller
* sysdeps/ieee754/flt-32/e_j1f.c (__ieee754_y1f): Force computations to occur in round to nearest mode when |x| >= 2.0
2014-05-23Fix log10 (1) in round-downward mode (bug 16977).Joseph Myers
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.
2014-05-14Fix log1pl (LDBL_MAX) in FE_UPWARD mode (bug 16564).Joseph Myers
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.
2014-05-14Fix acosh (1) in round-downward mode (bug 16927).Joseph Myers
According to C99 and C11 Annex F, acosh (1) should be +0 in all rounding modes. However, some implementations in glibc wrongly return -0 in round-downward mode (which is what you get if you end up computing log1p (-0), via 1 - 1 being -0 in round-downward mode). This patch fixes the problem implementations, by correcting the test for an exact 1 value in the ldbl-96 implementation to allow for the explicit high bit of the mantissa, and by inserting fabs instructions in the i386 implementations; tests of acosh are duly converted to ALL_RM_TEST. I believe all the other sysdeps/ieee754 implementations are already OK (I haven't checked the ia64 versions, but if buggy then that will be obvious from the results of test runs after this patch is in). Tested x86_64 and x86 and ulps updated accordingly. [BZ #16927] * sysdeps/i386/fpu/e_acosh.S (__ieee754_acosh): Use fabs on x-1 value. * sysdeps/i386/fpu/e_acoshf.S (__ieee754_acoshf): Likewise. * sysdeps/i386/fpu/e_acoshl.S (__ieee754_acoshl): Likewise. * sysdeps/ieee754/ldbl-96/e_acoshl.c (__ieee754_acoshl): Correct for explicit high bit of mantissa when testing for argument equal to 1. * math/libm-test.inc (acosh_test): Use ALL_RM_TEST. * sysdeps/i386/fpu/libm-test-ulps: Update. * sysdeps/x86_64/fpu/libm-test-ulps: Likewise.
2014-05-14Fix erf underflow handling near 0 (bug 16516).Joseph Myers
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.
2014-04-29[BZ #16823] Fix log1pl returning wrong infinity signStefan Liebler
2014-04-16[BZ #16824] Fix failing y1 due to too large ulps in downward/upward rounding ↵Stefan Liebler
mode.
2014-04-16Correct IBM long double frexpl.Alan Modra
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.
2014-04-02Fix catan, catanh, __ieee754_logf in round-downward mode (bug 16799, bug 16800).Joseph Myers
This patch fixes incorrect results from catan and catanh of certain special inputs in round-downward mode (bug 16799), and incorrect results of __ieee754_logf (+/-0) in round-downward mode (bug 16800) that show up through catan/catanh when tested in all rounding modes, but not directly in the testing for logf because the bug gets hidden by the wrappers. Both bugs involve a zero that should be +0 being -0 instead: one computed as (1-x)*(1+x) in the catan/catanh case, and one as (x-x) in the logf case. The fixes ensure positive zero is used. Testing of catan and catanh in all rounding modes is duly enabled. I expect there are various other bugs in special cases in __ieee754_* functions that are normally hidden by the wrappers but would show up for testing with -lieee (or in future with -fno-math-errno if we replace -lieee and _LIB_VERSION with compile-time redirection to new *_noerrno symbol names). Tested x86_64 and x86 and ulps updated accordingly. [BZ #16799] [BZ #16800] * math/s_catan.c (__catan): Avoid passing -0 denominator to atan2 with 0 numerator. * math/s_catanf.c (__catanf): Likewise. * math/s_catanh.c (__catanh): Likewise. * math/s_catanhf.c (__catanhf): Likewise. * math/s_catanhl.c (__catanhl): Likewise. * math/s_catanl.c (__catanl): Likewise. * sysdeps/ieee754/flt-32/e_logf.c (__ieee754_logf): Always divide by positive zero when computing -Inf result. * math/libm-test.inc (catan_test): Use ALL_RM_TEST. (catanh_test): Likewise. * sysdeps/i386/fpu/libm-test-ulps: Update. * sysdeps/x86_64/fpu/libm-test-ulps: Likewise.