Commit Graph

1205 Commits

Author SHA1 Message Date
H.J. Lu
d3e4f5a101 s_sincosf.h: Change pio4 type to float [BZ #28713]
s_cosf.c and s_sinf.c have

  if (abstop12 (y) < abstop12 (pio4))

where abstop12 takes a float argument, but pio4 is static const double.
pio4 is used only in calls to abstop12 and never in arithmetic.  Apply

-static const double pio4 = 0x1.921FB54442D18p-1;
+static const float pio4 = 0x1.921FB6p-1f;

to fix:

FAIL: math/test-float-cos
FAIL: math/test-float-sin
FAIL: math/test-float-sincos
FAIL: math/test-float32-cos
FAIL: math/test-float32-sin
FAIL: math/test-float32-sincos

when compiling with GCC 12.

Reviewed-by: Paul Zimmermann <Paul.Zimmermann@inria.fr>
2021-12-21 08:56:12 -08:00
Akila Welihinda
3b1402b3fc sysdeps: Simplify sin Taylor Series calculation
The macro TAYLOR_SIN adds the term `-0.5*da*a^2 + da` in hopes
of regaining some precision as a function of da. However the
comment says we add the term `-0.5*da*a^2 + 0.5*da` which is
different. This fix updates the comment to reflect the
code and also simplifies the calculation by replacing `a` with `x`
because they always have the same value.

Signed-off-by: Akila Welihinda <akilawelihinda@ucla.edu>
Reviewed-by: Paul Zimmermann <Paul.Zimmermann@inria.fr>
2021-12-13 15:31:05 +01:00
Adhemerval Zanella
104d2005d5 math: Remove the error handling wrapper from hypot and hypotf
The error handling is moved to sysdeps/ieee754 version with no SVID
support.  The compatibility symbol versions still use the wrapper with
SVID error handling around the new code.  There is no new symbol version
nor compatibility code on !LIBM_SVID_COMPAT targets (e.g. riscv).

Only ia64 is unchanged, since it still uses the arch specific
__libm_error_region on its implementation.

Checked on x86_64-linux-gnu, i686-linux-gnu, and aarch64-linux-gnu.
2021-12-13 10:08:46 -03:00
Wilco Dijkstra
2f44eef584 math: Use fmin/fmax on hypot
It optimizes for architectures that provides fast builtins.

Checked on aarch64-linux-gnu.
2021-12-13 10:08:46 -03:00
Adhemerval Zanella
c212d6397e math: Use an improved algorithm for hypotl (ldbl-128)
This implementation is based on 'An Improved Algorithm for hypot(a,b)'
by Carlos F. Borges [1] using the MyHypot3 with the following changes:

  - Handle qNaN and sNaN.
  - Tune the 'widely varying operands' to avoid spurious underflow
    due the multiplication and fix the return value for upwards
    rounding mode.
  - Handle required underflow exception for subnormal results.

The main advantage of the new algorithm is its precision.  With a
random 1e9 input pairs in the range of [LDBL_MIN, LDBL_MAX], glibc
current implementation shows around 0.05% results with an error of
1 ulp (453266 results) while the new implementation only shows
0.0001% of total (1280).

Checked on aarch64-linux-gnu and x86_64-linux-gnu.

[1] https://arxiv.org/pdf/1904.09481.pdf
2021-12-13 09:02:34 -03:00
Adhemerval Zanella
aa9c28cde3 math: Use an improved algorithm for hypotl (ldbl-96)
This implementation is based on 'An Improved Algorithm for hypot(a,b)'
by Carlos F. Borges [1] using the MyHypot3 with the following changes:

 - Handle qNaN and sNaN.
 - Tune the 'widely varying operands' to avoid spurious underflow
   due the multiplication and fix the return value for upwards
   rounding mode.
 - Handle required underflow exception for subnormal results.

The main advantage of the new algorithm is its precision.  With a
random 1e8 input pairs in the range of [LDBL_MIN, LDBL_MAX], glibc
current implementation shows around 0.02% results with an error of
1 ulp (23158 results) while the new implementation only shows
0.0001% of total (111).

[1] https://arxiv.org/pdf/1904.09481.pdf
2021-12-13 09:02:34 -03:00
Wilco Dijkstra
ccfa865a82 math: Improve hypot performance with FMA
Improve hypot performance significantly by using fma when available. The
fma version has twice the throughput of the previous version and 70% of
the latency.  The non-fma version has 30% higher throughput and 10%
higher latency.

Max ULP error is 0.949 with fma and 0.792 without fma.

Passes GLIBC testsuite.
2021-12-13 09:02:34 -03:00
Wilco Dijkstra
6c848d7038 math: Use an improved algorithm for hypot (dbl-64)
This implementation is based on the 'An Improved Algorithm for
hypot(a,b)' by Carlos F. Borges [1] using the MyHypot3 with the
following changes:

 - Handle qNaN and sNaN.
 - Tune the 'widely varying operands' to avoid spurious underflow
   due the multiplication and fix the return value for upwards
   rounding mode.
 - Handle required underflow exception for denormal results.

The main advantage of the new algorithm is its precision: with a
random 1e9 input pairs in the range of [DBL_MIN, DBL_MAX], glibc
current implementation shows around 0.34% results with an error of
1 ulp (3424869 results) while the new implementation only shows
0.002% of total (18851).

The performance result are also only slight worse than current
implementation.  On x86_64 (Ryzen 5900X) with gcc 12:

Before:

  "hypot": {
   "workload-random": {
    "duration": 3.73319e+09,
    "iterations": 1.12e+08,
    "reciprocal-throughput": 22.8737,
    "latency": 43.7904,
    "max-throughput": 4.37184e+07,
    "min-throughput": 2.28361e+07
   }
  }

After:

  "hypot": {
   "workload-random": {
    "duration": 3.7597e+09,
    "iterations": 9.8e+07,
    "reciprocal-throughput": 23.7547,
    "latency": 52.9739,
    "max-throughput": 4.2097e+07,
    "min-throughput": 1.88772e+07
   }
  }

Co-Authored-By: Adhemerval Zanella  <adhemerval.zanella@linaro.org>

Checked on x86_64-linux-gnu and aarch64-linux-gnu.

[1] https://arxiv.org/pdf/1904.09481.pdf
2021-12-13 09:02:34 -03:00
Adhemerval Zanella
7fe0ace3e2 math: Simplify hypotf implementation
Use a more optimized comparison for check for NaN and infinite and
add an inlined issignaling implementation for float.  With gcc it
results in 2 FP comparisons.

The file Copyright is also changed to use  GPL, the implementation was
completely changed by 7c10fd3515 to use double precision instead of
scaling and this change removes all the GET_FLOAT_WORD usage.

Checked on x86_64-linux-gnu.
2021-12-13 09:02:30 -03:00
Paul Zimmermann
6bbf729832 Fixed inaccuracy of j0f (BZ #28185)
The largest errors over the full binary32 range are after this
patch (on x86_64):

RNDN: libm wrong by up to 9.00e+00 ulp(s) [9] for x=0x1.04c39cp+6
RNDZ: libm wrong by up to 9.00e+00 ulp(s) [9] for x=0x1.04c39cp+6
RNDU: libm wrong by up to 9.00e+00 ulp(s) [9] for x=0x1.04c39cp+6
RNDD: libm wrong by up to 8.98e+00 ulp(s) [9] for x=0x1.4b7066p+7

Inputs that were yielding huge errors have been added to "make check".
Reviewed-by: Adhemeral Zanella  <adhemerval.zanella@linaro.org>
2021-10-05 13:45:37 +02:00
Joseph Myers
90f0ac10a7 Add fmaximum, fminimum functions
C2X adds new <math.h> functions for floating-point maximum and
minimum, corresponding to the new operations that were added in IEEE
754-2019 because of concerns about the old operations not being
associative in the presence of signaling NaNs.  fmaximum and fminimum
handle NaNs like most <math.h> functions (any NaN argument means the
result is a quiet NaN).  fmaximum_num and fminimum_num handle both
quiet and signaling NaNs the way fmax and fmin handle quiet NaNs (if
one argument is a number and the other is a NaN, return the number),
but still raise "invalid" for a signaling NaN argument, making them
exceptions to the normal rule that a function with a floating-point
result raising "invalid" also returns a quiet NaN.  fmaximum_mag,
fminimum_mag, fmaximum_mag_num and fminimum_mag_num are corresponding
functions returning the argument with greatest or least absolute
value.  All these functions also treat +0 as greater than -0.  There
are also corresponding <tgmath.h> type-generic macros.

Add these functions to glibc.  The implementations use type-generic
templates based on those for fmax, fmin, fmaxmag and fminmag, and test
inputs are based on those for those functions with appropriate
adjustments to the expected results.  The RISC-V maintainers might
wish to add optimized versions of fmaximum_num and fminimum_num (for
float and double), since RISC-V (F extension version 2.2 and later)
provides instructions corresponding to those functions - though it
might be at least as useful to add architecture-independent built-in
functions to GCC and teach the RISC-V back end to expand those
functions inline, which is what you generally want for functions that
can be implemented with a single instruction.

Tested for x86_64 and x86, and with build-many-glibcs.py.
2021-09-28 23:31:35 +00:00
Tulio Magno Quites Machado Filho
54ff4f1e39 powerpc64le: Avoid conflicting types for f64xfmaf128 when IFUNC is not used
Avoid defining f64xfmaf128 twice when building s_fmaf128.c.
This can be reproduced on powerpc64le whenever f128 functions do not
have IFUNC enabled, e.g. using "--with-cpu=power8 --disable-multi-arch", or
when using "-with-cpu=power9".

Fixes: b3f27d8150 ("Add narrowing fma functions")
2021-09-23 19:29:54 -03:00
Joseph Myers
b3f27d8150 Add narrowing fma functions
This patch adds the narrowing fused multiply-add functions from TS
18661-1 / TS 18661-3 / C2X to glibc's libm: ffma, ffmal, dfmal,
f32fmaf64, f32fmaf32x, f32xfmaf64 for all configurations; f32fmaf64x,
f32fmaf128, f64fmaf64x, f64fmaf128, f32xfmaf64x, f32xfmaf128,
f64xfmaf128 for configurations with _Float64x and _Float128;
__f32fmaieee128 and __f64fmaieee128 aliases in the powerpc64le case
(for calls to ffmal and dfmal when long double is IEEE binary128).
Corresponding tgmath.h macro support is also added.

The changes are mostly similar to those for the other narrowing
functions previously added, especially that for sqrt, so the
description of those generally applies to this patch as well.  As with
sqrt, I reused the same test inputs in auto-libm-test-in as for
non-narrowing fma rather than adding extra or separate inputs for
narrowing fma.  The tests in libm-test-narrow-fma.inc also follow
those for non-narrowing fma.

The non-narrowing fma has a known bug (bug 6801) that it does not set
errno on errors (overflow, underflow, Inf * 0, Inf - Inf).  Rather
than fixing this or having narrowing fma check for errors when
non-narrowing does not (complicating the cases when narrowing fma can
otherwise be an alias for a non-narrowing function), this patch does
not attempt to check for errors from narrowing fma and set errno; the
CHECK_NARROW_FMA macro is still present, but as a placeholder that
does nothing, and this missing errno setting is considered to be
covered by the existing bug rather than needing a separate open bug.
missing-errno annotations are duly added to many of the
auto-libm-test-in test inputs for fma.

This completes adding all the new functions from TS 18661-1 to glibc,
so will be followed by corresponding stdc-predef.h changes to define
__STDC_IEC_60559_BFP__ and __STDC_IEC_60559_COMPLEX__, as the support
for TS 18661-1 will be at a similar level to that for C standard
floating-point facilities up to C11 (pragmas not implemented, but
library functions done).  (There are still further changes to be done
to implement changes to the types of fromfp functions from N2548.)

Tested as followed: natively with the full glibc testsuite for x86_64
(GCC 11, 7, 6) and x86 (GCC 11); with build-many-glibcs.py with GCC
11, 7 and 6; cross testing of math/ tests for powerpc64le, powerpc32
hard float, mips64 (all three ABIs, both hard and soft float).  The
different GCC versions are to cover the different cases in tgmath.h
and tgmath.h tests properly (GCC 6 has _Float* only as typedefs in
glibc headers, GCC 7 has proper _Float* support, GCC 8 adds
__builtin_tgmath).
2021-09-22 21:25:31 +00:00
Joseph Myers
1356f38df5 Fix f64xdivf128, f64xmulf128 spurious underflows (bug 28358)
As described in bug 28358, the round-to-odd computations used in the
libm functions that round their results to a narrower format can yield
spurious underflow exceptions in the following circumstances: the
narrowing only narrows the precision of the type and not the exponent
range (i.e., it's narrowing _Float128 to _Float64x on x86_64, x86 or
ia64), the architecture does after-rounding tininess detection (which
applies to all those architectures), the result is inexact, tiny
before rounding but not tiny after rounding (with the chosen rounding
mode) for _Float64x (which is possible for narrowing mul, div and fma,
not for narrowing add, sub or sqrt), so the underflow exception
resulting from the toward-zero computation in _Float128 is spurious
for _Float64x.

Fixed by making ROUND_TO_ODD call feclearexcept (FE_UNDERFLOW) in the
problem cases (as indicated by an extra argument to the macro); there
is never any need to preserve underflow exceptions from this part of
the computation, because the conversion of the round-to-odd value to
the narrower type will underflow in exactly the cases in which the
function should raise that exception, but it may be more efficient to
avoid the extra manipulation of the floating-point environment when
not needed.

Tested for x86_64 and x86, and with build-many-glibcs.py.
2021-09-21 21:54:37 +00:00
Joseph Myers
4b6574a6f6 Redirect fma calls to __fma in libm
include/math.h has a mechanism to redirect internal calls to various
libm functions, that can often be inlined by the compiler, to call
non-exported __* names for those functions in the case when the calls
aren't inlined, with the redirection being disabled when
NO_MATH_REDIRECT.  Add fma to the functions to which this mechanism is
applied.

At present, libm-internal fma calls (generally to __builtin_fma*
functions) are only done when it's known the call will be inlined,
with alternative code not relying on an fma operation being used in
the caller otherwise.  This patch is in preparation for adding the TS
18661 / C2X narrowing fma functions to glibc; it will be natural for
the narrowing function implementations to call the underlying fma
functions unconditionally, with this either being inlined or resulting
in an __fma* call.  (Using two levels of round-to-odd computation like
that, in the case where there isn't an fma hardware instruction, isn't
optimal but is certainly a lot simpler for the initial implementation
than writing different narrowing fma implementations for all the
various pairs of formats.)

Tested with build-many-glibcs.py that installed stripped shared
libraries are unchanged by the patch (using
<https://sourceware.org/pipermail/libc-alpha/2021-September/130991.html>
to fix installed library stripping in build-many-glibcs.py).  Also
tested for x86_64.
2021-09-15 22:57:35 +00:00
Joseph Myers
abd383584b Add narrowing square root functions
This patch adds the narrowing square root functions from TS 18661-1 /
TS 18661-3 / C2X to glibc's libm: fsqrt, fsqrtl, dsqrtl, f32sqrtf64,
f32sqrtf32x, f32xsqrtf64 for all configurations; f32sqrtf64x,
f32sqrtf128, f64sqrtf64x, f64sqrtf128, f32xsqrtf64x, f32xsqrtf128,
f64xsqrtf128 for configurations with _Float64x and _Float128;
__f32sqrtieee128 and __f64sqrtieee128 aliases in the powerpc64le case
(for calls to fsqrtl and dsqrtl when long double is IEEE binary128).
Corresponding tgmath.h macro support is also added.

The changes are mostly similar to those for the other narrowing
functions previously added, so the description of those generally
applies to this patch as well.  However, the not-actually-narrowing
cases (where the two types involved in the function have the same
floating-point format) are aliased to sqrt, sqrtl or sqrtf128 rather
than needing a separately built not-actually-narrowing function such
as was needed for add / sub / mul / div.  Thus, there is no
__nldbl_dsqrtl name for ldbl-opt because no such name was needed
(whereas the other functions needed such a name since the only other
name for that entry point was e.g. f32xaddf64, not reserved by TS
18661-1); the headers are made to arrange for sqrt to be called in
that case instead.

The DIAG_* calls in sysdeps/ieee754/soft-fp/s_dsqrtl.c are because
they were observed to be needed in GCC 7 testing of
riscv32-linux-gnu-rv32imac-ilp32.  The other sysdeps/ieee754/soft-fp/
files added didn't need such DIAG_* in any configuration I tested with
build-many-glibcs.py, but if they do turn out to be needed in more
files with some other configuration / GCC version, they can always be
added there.

I reused the same test inputs in auto-libm-test-in as for
non-narrowing sqrt rather than adding extra or separate inputs for
narrowing sqrt.  The tests in libm-test-narrow-sqrt.inc also follow
those for non-narrowing sqrt.

Tested as followed: natively with the full glibc testsuite for x86_64
(GCC 11, 7, 6) and x86 (GCC 11); with build-many-glibcs.py with GCC
11, 7 and 6; cross testing of math/ tests for powerpc64le, powerpc32
hard float, mips64 (all three ABIs, both hard and soft float).  The
different GCC versions are to cover the different cases in tgmath.h
and tgmath.h tests properly (GCC 6 has _Float* only as typedefs in
glibc headers, GCC 7 has proper _Float* support, GCC 8 adds
__builtin_tgmath).
2021-09-10 20:56:22 +00:00
Siddhesh Poyarekar
30891f35fa Remove "Contributed by" lines
We stopped adding "Contributed by" or similar lines in sources in 2012
in favour of git logs and keeping the Contributors section of the
glibc manual up to date.  Removing these lines makes the license
header a bit more consistent across files and also removes the
possibility of error in attribution when license blocks or files are
copied across since the contributed-by lines don't actually reflect
reality in those cases.

Move all "Contributed by" and similar lines (Written by, Test by,
etc.) into a new file CONTRIBUTED-BY to retain record of these
contributions.  These contributors are also mentioned in
manual/contrib.texi, so we just maintain this additional record as a
courtesy to the earlier developers.

The following scripts were used to filter a list of files to edit in
place and to clean up the CONTRIBUTED-BY file respectively.  These
were not added to the glibc sources because they're not expected to be
of any use in future given that this is a one time task:

https://gist.github.com/siddhesh/b5ecac94eabfd72ed2916d6d8157e7dc
https://gist.github.com/siddhesh/15ea1f5e435ace9774f485030695ee02

Reviewed-by: Carlos O'Donell <carlos@redhat.com>
2021-09-03 22:06:44 +05:30
H.J. Lu
3213ed770c Update math: redirect roundeven function
Redirect target specific roundeven functions for aarch64, ldbl-128ibm
and riscv.
2021-06-27 07:56:57 -07:00
Shen-Ta Hsieh
eb9066203f Use GCC builtins for roundeven functions if desired.
This patch is using the corresponding GCC builtin for roundevenf,
roundeven and roundevenl if the USE_FUNCTION_BUILTIN macros are defined
to one in math-use-builtins.h.

These builtin functions is supported since GCC 10.

The code of the generic implementation is not changed.

Signed-off-by: Shen-Ta Hsieh <ibmibmibm.tw@gmail.com>
Reviewed-by: H.J. Lu <hjl.tools@gmail.com>
2021-06-27 07:56:57 -07:00
Shen-Ta Hsieh
447954a206 math: redirect roundeven function
This patch redirect roundeven function for futhermore changes.

Signed-off-by: Shen-Ta Hsieh <ibmibmibm.tw@gmail.com>
Reviewed-by: H.J. Lu <hjl.tools@gmail.com>
2021-06-27 07:56:57 -07:00
Naohiro Tamura
b190bccc8a configure: Replaced obsolete AC_TRY_COMPILE
This patch replaced obsolete AC_TRY_COMPILE to AC_COMPILE_IFELSE or
AC_PREPROC_IFELSE.
It has been confirmed that GNU 'autoconf' 2.69 suppressed obsolete
warnings, updated the following files:
  - configure
  - sysdeps/mach/configure
  - sysdeps/mach/hurd/configure
  - sysdeps/s390/configure
  - sysdeps/unix/sysv/linux/configure
and didn't change the following files:
  - sysdeps/ieee754/ldbl-opt/configure
  - sysdeps/unix/sysv/linux/powerpc/configure

Reviewed-by: Adhemerval Zanella  <adhemerval.zanella@linaro.org>
2021-06-04 10:16:00 -03:00
Florian Weimer
c8a11c5867 stdio-common: Remove _IO_vfwscanf
The symbol has never been exported, so no compatibility symbol is
needed.  Removing this file prevents ld from creation an exported
symbol in case GLIBC_2_0 expands to a symbol version which
does not have a local: *; directive in the symbol version map file.

Reviewed-by: Adhemerval Zanella  <adhemerval.zanella@linaro.org>
2021-06-01 16:00:52 +02:00
Paul Zimmermann
8d0985b055 add workload traces for cbrtl
These workload traces cover the whole "long double" range.
This patch was prepared with the help of Adhemerval Zanella.
Reviewed-by: Carlos O'Donell <carlos@redhat.com>
2021-05-10 18:45:34 +02:00
Paul Zimmermann
43576de04a Improve the accuracy of tgamma (BZ #26983)
With this patch, the maximal known error for tgamma is now reduced to 9 ulps
for dbl-64, for all rounding modes. Since exhaustive testing is not possible
for dbl-64, it might be that there are still cases with an error larger than
9 ulps, but all known cases are fixed (intensive tests were done to find cases
with large errors).

Tested on x86_64 and powerpc (and by Adhemerval Zanella on aarch64, arm,
s390x, sparc, and i686).
Reviewed-by: Adhemerval Zanella  <adhemerval.zanella@linaro.org>
2021-04-07 13:23:39 +02:00
Paul Zimmermann
9acda61d94 Fix the inaccuracy of j0f/j1f/y0f/y1f [BZ #14469, #14470, #14471, #14472]
For j0f/j1f/y0f/y1f, the largest error for all binary32
inputs is reduced to at most 9 ulps for all rounding modes.

The new code is enabled only when there is a cancellation at the very end of
the j0f/j1f/y0f/y1f computation, or for very large inputs, thus should not
give any visible slowdown on average.  Two different algorithms are used:

* around the first 64 zeros of j0/j1/y0/y1, approximation polynomials of
  degree 3 are used, computed using the Sollya tool (https://www.sollya.org/)

* for large inputs, an asymptotic formula from [1] is used

[1] Fast and Accurate Bessel Function Computation,
    John Harrison, Proceedings of Arith 19, 2009.

Inputs yielding the new largest errors are added to auto-libm-test-in,
and ulps are regenerated for various targets (thanks Adhemerval Zanella).

Tested on x86_64 with --disable-multi-arch and on powerpc64le-linux-gnu.
Reviewed-by: Adhemerval Zanella  <adhemerval.zanella@linaro.org>
2021-04-02 06:15:48 +02:00
Siddhesh Poyarekar
abadbef5c8 Move __isnanf128 to libc.so
All of the isnan functions are in libc.so due to printf_fp, so move
__isnanf128 there too for consistency.

Reviewed-by: Tulio Magno Quites Machado Filho <tuliom@ascii.art.br>
Reviewed-by: Florian Weimer <fweimer@redhat.com>
2021-03-30 14:58:19 +05:30
Wilco Dijkstra
92cfc9ad82 math: Remove mpa files (part 2) [BZ #15267]
Previous commit was missing deleted files in sysdeps/ieee754/dbl-64.

Finally remove all mpa related files, headers, declarations, probes, unused
tables and update makefiles.

Reviewed-By: Paul Zimmermann <Paul.Zimmermann@inria.fr>
2021-03-11 15:45:19 +00:00
Wilco Dijkstra
47ad14d789 math: Remove mpa files [BZ #15267]
Finally remove all mpa related files, headers, declarations, probes, unused
tables and update makefiles.

Reviewed-By: Paul Zimmermann <Paul.Zimmermann@inria.fr>
2021-03-11 14:26:36 +00:00
Wilco Dijkstra
4e1a870b9a math: Remove slow paths from atan2 [BZ #15267]
Remove slow paths from atan2. Add ULP annotations.

Reviewed-By: Paul Zimmermann <Paul.Zimmermann@inria.fr>
2021-03-11 14:26:36 +00:00
Wilco Dijkstra
e898cd1593 math: Remove slow paths from atan [BZ #15267]
Remove slow paths from atan. Add ULP annotations.

Reviewed-By: Paul Zimmermann <Paul.Zimmermann@inria.fr>
2021-03-11 14:26:36 +00:00
Wilco Dijkstra
476d692e8a math: Remove slow paths in tan [BZ #15267]
Remove slow paths in tan. Add ULP annotations. Merge 'number' into 'mynumber'.
Remove unused entries from tan constants.

Reviewed-By: Paul Zimmermann <Paul.Zimmermann@inria.fr>
2021-03-11 14:26:36 +00:00
Wilco Dijkstra
db3f7bb558 math: Remove slow paths from asin and acos [BZ #15267]
This patch series removes all remaining slow paths and related code.
First asin/acos, tan, atan, atan2 implementations are updated, and the final
patch removes the unused mpa files, headers and probes. Passes buildmanyglibc.

Remove slow paths from asin/acos. Add ULP annotations based on previous slow
path checks (which are approximate). Update AArch64 and x86_64 libm-test-ulps.

Reviewed-By: Paul Zimmermann <Paul.Zimmermann@inria.fr>
2021-03-11 14:26:36 +00:00
Adhemerval Zanella
bf7db6d369 math: Add BZ#18980 fix back on dbl-64 cosh
It is regression from 9e97f239ea (Remove dbl-64/wordsize-64
(part 2)) where is missed to add the BZ#18980 fix (9e97f239ea).

Checked on i686-linux-gnu.
2021-01-11 16:56:33 -03:00
Wilco Dijkstra
9e97f239ea Remove dbl-64/wordsize-64 (part 2)
Remove the wordsize-64 implementations by merging them into the main dbl-64
directory.  The second patch just moves all wordsize-64 files and removes a
few wordsize-64 uses in comments and Implies files.

Reviewed-by: Adhemerval Zanella  <adhemerval.zanella@linaro.org>
2021-01-07 15:26:26 +00:00
Wilco Dijkstra
caa884dda7 Remove dbl-64/wordsize-64
Remove the wordsize-64 implementations by merging them into the main dbl-64
directory.  The first patch adds special cases needed for 32-bit targets
(FIX_INT_FP_CONVERT_ZERO and FIX_DBL_LONG_CONVERT_OVERFLOW) to the
wordsize-64 versions.  This has no effect on 64-bit targets since they don't
define these macros.

Reviewed-by: Adhemerval Zanella  <adhemerval.zanella@linaro.org>
2021-01-07 15:02:51 +00:00
Siddhesh Poyarekar
8cc1e39a36 Drop nan-pseudo-number.h usage from tests
Make the tests use TEST_COND_intel96 to decide on whether to build the
unnormal tests instead of the macro in nan-pseudo-number.h and then
drop the header inclusion.  This unbreaks test runs on all
architectures that do not have ldbl-96.

Also drop the HANDLE_PSEUDO_NUMBERS macro since it is not used
anywhere.
2021-01-04 20:49:56 +05:30
Siddhesh Poyarekar
fee3b889d8 Move generic nan-pseudo-number.h to ldbl-96
The concept of pseudo number formats only exists in the realm of the
96 bit long double format.
2021-01-04 14:51:52 +05:30
Paul Eggert
2b778ceb40 Update copyright dates with scripts/update-copyrights
I used these shell commands:

../glibc/scripts/update-copyrights $PWD/../gnulib/build-aux/update-copyright
(cd ../glibc && git commit -am"[this commit message]")

and then ignored the output, which consisted lines saying "FOO: warning:
copyright statement not found" for each of 6694 files FOO.
I then removed trailing white space from benchtests/bench-pthread-locks.c
and iconvdata/tst-iconv-big5-hkscs-to-2ucs4.c, to work around this
diagnostic from Savannah:
remote: *** pre-commit check failed ...
remote: *** error: lines with trailing whitespace found
remote: error: hook declined to update refs/heads/master
2021-01-02 12:17:34 -08:00
Siddhesh Poyarekar
7525c1c71d x86 long double: Consider pseudo numbers as signaling
Add support to treat pseudo-numbers specially and implement x86
version to consider all of them as signaling.

Reviewed-by: Adhemerval Zanella <adhemerval.zanella@linaro.org>
2020-12-30 10:52:45 +05:30
Anssi Hannula
69a7ca7705 ieee754: Remove unused __sin32 and __cos32
The __sin32 and __cos32 functions were only used in the now removed slow
path of asin and acos.
2020-12-18 12:10:31 +05:30
Anssi Hannula
f67f9c9af2 ieee754: Remove slow paths from asin and acos
asin and acos have slow paths for rounding the last bit that cause some
calls to be 500-1500x slower than average calls.

These slow paths are rare, a test of a trillion (1.000.000.000.000)
random inputs between -1 and 1 showed 32870 slow calls for acos and 4473
for asin, with most occurrences between -1.0 .. -0.9 and 0.9 .. 1.0.

The slow paths claim correct rounding and use __sin32() and __cos32()
(which compare two result candidates and return the closest one) as the
final step, with the second result candidate (res1) having a small offset
applied from res. This suggests that res and res1 are intended to be 1
ULP apart (which makes sense for rounding), barring bugs, allowing us to
pick either one and still remain within 1 ULP of the exact result.

Remove the slow paths as the accuracy is better than 1 ULP even without
them, which is enough for glibc.

Also remove code comments claiming correctly rounded results.

After slow path removal, checking the accuracy of 14.400.000.000 random
asin() and acos() inputs showed only three incorrectly rounded
(error > 0.5 ULP) results:
- asin(-0x1.ee2b43286db75p-1) (0.500002 ULP, same as before)
- asin(-0x1.f692ba202abcp-4)  (0.500003 ULP, same as before)
- asin(-0x1.9915e876fc062p-1) (0.50000000001 ULP, previously exact)
The first two had the same error even before this commit, and they did
not use the slow path at all.

Checking 4934 known randomly found previously-slow-path asin inputs
shows 25 calls with incorrectly rounded results, with a maximum error of
0.500000002 ULP (for 0x1.fcd5742999ab8p-1). The previous slow-path code
rounded all these inputs correctly (error < 0.5 ULP).
The observed average speed increase was 130x.

Checking 36240 known randomly found previously-slow-path acos inputs
shows 42 calls with incorrectly rounded results, with a maximum error of
0.500000008 ULP (for 0x1.f63845056f35ep-1). The previous "exact"
slow-path code showed 34 calls with incorrectly rounded results, with the
same maximum error of 0.500000008 ULP (for 0x1.f63845056f35ep-1).
The observed average speed increase was 130x.

The functions could likely be trimmed more while keeping acceptable
accuracy, but this at least gets rid of the egregiously slow cases.

Tested on x86_64.
2020-12-18 12:09:23 +05:30
Adhemerval Zanella
01bd62517c Remove tls.h inclusion from internal errno.h
The tls.h inclusion is not really required and limits possible
definition on more arch specific headers.

This is a cleanup to allow inline functions on sysdep.h, more
specifically on i386 and ia64 which requires to access some tls
definitions its own.

No semantic changes expected, checked with a build against all
affected ABIs.
2020-11-13 12:59:19 -03:00
Paul Zimmermann
b7dd366dbe math: Fix inaccuracy of j0f for x >= 2^127 when sin(x)+cos(x) is tiny
Checked on x86_64-linux-gnu and i686-linux-gnu.
2020-08-07 16:33:13 -03:00
Joseph Myers
6c010c5dde Use C2x return value from getpayload of non-NaN (bug 26073).
In TS 18661-1, getpayload had an unspecified return value for a
non-NaN argument, while C2x requires the return value -1 in that case.

This patch implements the return value of -1.  I don't think this is
worth having a new symbol version that's an alias of the old one,
although occasionally we do that in such cases where the new function
semantics are a refinement of the old ones (to avoid programs relying
on the new semantics running on older glibc versions but not behaving
as intended).

Tested for x86_64 and x86; also ran math/ tests for aarch64 and
powerpc.
2020-07-06 16:18:02 +00:00
Adhemerval Zanella
be668a8d78 New exp10f version without SVID compat wrapper
This patch changes the exp10f error handling semantics to only set
errno according to POSIX rules.  New symbol version is introduced at
GLIBC_2.32.  The old wrappers are kept for compat symbols.

There are some outliers that need special handling:

  - ia64 provides an optimized implementation of exp10f that uses ia64
    specific routines to set SVID compatibility.  The new symbol version
    is aliased to the exp10f one.

  - m68k also provides an optimized implementation, and the new version
    uses it instead of the sysdeps/ieee754/flt32 one.

  - riscv and csky uses the generic template implementation that
    does not provide SVID support.  For both cases a new exp10f
    version is not added, but rather the symbols version of the
    generic sysdeps/ieee754/flt32 is adjusted instead.

Checked on aarch64-linux-gnu, x86_64-linux-gnu, i686-linux-gnu,
powerpc64le-linux-gnu.
2020-06-19 12:08:47 -03:00
Paul Zimmermann
6e98983c09 math: Optimized generic exp10f with wrappers
It is inspired by expf and reuses its tables and internal functions.
The error checks are inlined and errno setting is in separate tail
called functions, but the wrappers are kept in this patch to handle
the _LIB_VERSION==_SVID_ case.

Double precision arithmetics is used which is expected to be faster on
most targets (including soft-float) than using single precision and it
is easier to get good precision result with it.

Result for x86_64 (i7-4790K CPU @ 4.00GHz) are:

Before new code:
  "exp10f": {
   "workload-spec2017.wrf (adapted)": {
    "duration": 4.0414e+09,
    "iterations": 1.00128e+08,
    "reciprocal-throughput": 26.6818,
    "latency": 54.043,
    "max-throughput": 3.74787e+07,
    "min-throughput": 1.85038e+07
   }

With new code:
  "exp10f": {
   "workload-spec2017.wrf (adapted)": {
    "duration": 4.11951e+09,
    "iterations": 1.23968e+08,
    "reciprocal-throughput": 21.0581,
    "latency": 45.4028,
    "max-throughput": 4.74876e+07,
    "min-throughput": 2.20251e+07
   }

Result for aarch64 (A72 @ 2GHz) are:

Before new code:
  "exp10f": {
   "workload-spec2017.wrf (adapted)": {
    "duration": 4.62362e+09,
    "iterations": 3.3376e+07,
    "reciprocal-throughput": 127.698,
    "latency": 149.365,
    "max-throughput": 7.831e+06,
    "min-throughput": 6.69501e+06
   }

With new code:
  "exp10f": {
   "workload-spec2017.wrf (adapted)": {
    "duration": 4.29108e+09,
    "iterations": 6.6752e+07,
    "reciprocal-throughput": 51.2111,
    "latency": 77.3568,
    "max-throughput": 1.9527e+07,
    "min-throughput": 1.29271e+07
   }

Checked on x86_64-linux-gnu, powerpc64le-linux-gnu, aarch64-linux-gnu,
and sparc64-linux-gnu.
2020-06-19 10:48:15 -03:00
Vineet Gupta
e93c264336 ieee754/dbl-64: Reduce the scope of temporary storage variables
This came to light when adding hard-flaot support to ARC glibc port
without hardware sqrt support causing glibc build to fail:

| ../sysdeps/ieee754/dbl-64/e_sqrt.c: In function '__ieee754_sqrt':
| ../sysdeps/ieee754/dbl-64/e_sqrt.c:58:54: error: unused variable 'ty' [-Werror=unused-variable]
|   double y, t, del, res, res1, hy, z, zz, p, hx, tx, ty, s;

The reason being EMULV() macro uses the hardware provided
__builtin_fma() variant, leaving temporary variables 'p, hx, tx, hy, ty'
unused hence compiler warning and ensuing error.

The intent of the patch was to fix that error, but EMULV is pervasive
and used fair bit indirectly via othe rmacros, hence this patch.
Functionally it should not result in code gen changes and if at all
those would be better since the scope of those temporaries is greatly
reduced now

Built tested with aarch64-linux-gnu arm-linux-gnueabi arm-linux-gnueabihf hppa-linux-gnu x86_64-linux-gnu arm-linux-gnueabihf riscv64-linux-gnu-rv64imac-lp64 riscv64-linux-gnu-rv64imafdc-lp64 powerpc-linux-gnu microblaze-linux-gnu nios2-linux-gnu hppa-linux-gnu

Also as suggested by Joseph [1] used --strip and compared the libs with
and w/o patch and they are byte-for-byte unchanged (with gcc 9).

| for i in `find . -name libm-2.31.9000.so`;
| do
|    echo $i; diff $i /SCRATCH/vgupta/gnu2/install/glibcs/$i ; echo $?;
| done

| ./aarch64-linux-gnu/lib64/libm-2.31.9000.so
| 0
| ./arm-linux-gnueabi/lib/libm-2.31.9000.so
| 0
| ./x86_64-linux-gnu/lib64/libm-2.31.9000.so
| 0
| ./arm-linux-gnueabihf/lib/libm-2.31.9000.so
| 0
| ./riscv64-linux-gnu-rv64imac-lp64/lib64/lp64/libm-2.31.9000.so
| 0
| ./riscv64-linux-gnu-rv64imafdc-lp64/lib64/lp64/libm-2.31.9000.so
| 0
| ./powerpc-linux-gnu/lib/libm-2.31.9000.so
| 0
| ./microblaze-linux-gnu/lib/libm-2.31.9000.so
| 0
| ./nios2-linux-gnu/lib/libm-2.31.9000.so
| 0
| ./hppa-linux-gnu/lib/libm-2.31.9000.so
| 0
| ./s390x-linux-gnu/lib64/libm-2.31.9000.so

[1] https://sourceware.org/pipermail/libc-alpha/2019-November/108267.html
2020-06-15 13:09:21 -07:00
Vineet Gupta
628d90c5f9 ieee754: provide gcc builtins based generic fma functions
Reviewed-by: Adhemerval Zanella  <adhemerval.zanella@linaro.org>
2020-06-03 10:23:28 -07:00
Vineet Gupta
3374868668 ieee754: provide gcc builtins based generic sqrt functions
Reviewed-by: Adhemerval Zanella  <adhemerval.zanella@linaro.org>
2020-06-03 10:23:22 -07:00
Paul E. Murphy
a49e56a945 float128: use builtin_signbitf128 always
The minimum GCC version has been raised to 6.2 for building
glibc.  Therefore, follow the advice inside the implementation
and remove the GCC < 6 codepath.

Likewise, remove the hidden_proto as all internal usages should
inline now.
2020-05-04 13:18:45 -05:00
Paul E. Murphy
5c7ccc2983 powerpc64le: blacklist broken GCC compilers (e.g GCC 7.5.0)
GCC 7.5.0 (PR94200) will refuse to compile if both -mabi=% and
-mlong-double-128 are passed on the command line.  Surprisingly,
it will work happily if the latter is not.  For the sake of
maintaining status quo, test for and blacklist such compilers.

Tested with a GCC 8.3.1 and GCC 7.5.0 compiler for ppc64le.

Reviewed-by: Tulio Magno Quites Machado Filho <tuliom@linux.ibm.com>
2020-04-30 08:52:08 -05:00
Paul E. Murphy
e2239af353 Rename __LONG_DOUBLE_USES_FLOAT128 to __LDOUBLE_REDIRECTS_TO_FLOAT128_ABI
Improve the commentary to aid future developers who will stumble
upon this novel, yet not always perfect, mechanism to support
alternative formats for long double.

Likewise, rename __LONG_DOUBLE_USES_FLOAT128 to
__LDOUBLE_REDIRECTS_TO_FLOAT128_ABI now that development work
has settled down.  The command used was

git grep -l __LONG_DOUBLE_USES_FLOAT128 ':!./ChangeLog*' | \
  xargs sed -i 's/__LONG_DOUBLE_USES_FLOAT128/__LDOUBLE_REDIRECTS_TO_FLOAT128_ABI/g'

Reviewed-by: Tulio Magno Quites Machado Filho <tuliom@linux.ibm.com>
2020-04-30 08:52:08 -05:00
Paul E. Murphy
bd98471eb2 ldbl-128ibm: simplify iscanonical.h
The test for enabling _Float128 or IEEE 128 long double can be
greatly simplified knowing that there is no ibm128, thus we require
no special cases, and everything is canonical.

This reverts the changes to ldbl-128ibm iscanonical.h from commit
8dbfea3a20 and extends the check
for __NO_LONG_DOUBLE_MATH to include a check for float128 redirects
to long double.

Reviewed-by: Tulio Magno Quites Machado Filho <tuliom@linux.ibm.com>
2020-04-06 10:23:58 -05:00
Paul E. Murphy
3618e5fece Makeconfig: sandwich gnulib-tests between libc/ld linking of tests
This better resembles the default linking process with the gnulibs,
and also resolves the increasingly difficult to maintain
f128-loader-link usage on powerpc64le as some libgcc symbols are
dependent on those found in the loader (ld).
2020-03-25 14:34:23 -05:00
Paul E. Murphy
45ae17dd7e ldbl-128ibm-compat: PLT redirects for using ldbl redirects internally
Tweak the PLT bypass magic when building glibc with long double
redirects.  This is made more difficult by the fact we only get
one chance to redirect functions.  This happens via the public
headers.

There are roughly three classes of redirect we need to attend to
today:

 1. Simple redirects, redirected via cdef macro overrides and
    and new libc_hidden_ldbl_proto macro.
 2. Internal usage of internal API, e.g __snprintf, which has
    no direct analogue.  This is bypassed directly on case-by-
    case basis.
 3. Double redirects, e.g sscanf and related.  These require
    a heavier handed approach of macro renaming to existing
    symbols.

Most simple redirects are handled via 1.  Ideally, the libc_*
macro would live in libc-symbols.h, but in practice the macros
needed for it to do anything useful live in cdefs.h, so they
are defined in the local override.

Notably, the internal name of the asprintf generated for ieee ldbl
redirects is renamed to work with internal prefixed usage.

This resolves the local plt usage introduced when building glibc
with ldbl == ieee128 on ppc64le.

Reviewed-by: Tulio Magno Quites Machado Filho <tuliom@linux.ibm.com>
2020-03-25 14:34:23 -05:00
Adhemerval Zanella
1c15464ca0 math: Remove inline math tests
With mathinline removal there is no need to keep building and testing
inline math tests.

The gen-libm-tests.py support to generate ULP_I_* is removed and all
libm-test-ulps files are updated to longer have the
i{float,double,ldouble} entries.  The support for no-test-inline is
also removed from both gen-auto-libm-tests and the
auto-libm-test-out-* were regenerated.

Checked on x86_64-linux-gnu and i686-linux-gnu.
2020-03-19 11:45:44 -03:00
Tulio Magno Quites Machado Filho
d2b2acefa6 ldbl-128ibm: Let long double files have specific compiler flags
Soon, powerpc64le will need to provide extra compiler flags to the long
double files in order to continue to build using the IBM 128-bit
extended floating point type as long double.
2020-03-06 09:17:32 -06:00
Rajalakshmi Srinivasaraghavan
0059122aa0 ldbl-128ibm-compat: Add tests for IBM long double functions
This patch creates test-ibm128* tests from the long double function tests.
In order to explicitly test IBM long double functions -mabi=ibmlongdouble is
added to CFLAGS.

Likewise, update the test headers to correct choose ULPs when redirects
are enabled.

Co-authored-by: Tulio Magno Quites Machado Filho  <tuliom@linux.ibm.com>
Co-authored-by: Paul E. Murphy  <murphyp@linux.vnet.ibm.com>
2020-03-06 09:17:32 -06:00
Paul E. Murphy
4068b54c70 ldbl-128ibm-compat: link tst-ldbl-efgcvt against loader too
This also requires the linker workaround to ensure everything links
correctly.  See comment in sysdeps/powerpc/powerpc64/le/Makefile
for details.
2020-02-28 08:20:02 -06:00
Paul E. Murphy
77d5bddd8f ldbl-128ibm-compat: enforce ibm128 on compat tests
For lack of a more comprehensive solution, tack on the ibm128 ABI
compiler options for the totalorder{,mag}l compat tests which exist
prior to enabling this feature.
2020-02-28 08:20:02 -06:00
Gabriel F. T. Gomes
d967237443 ldbl-128ibm-compat: Provide nexttoward functions
The functions in the nexttoward family are special, in the sense that
they always have a long double argument, regardless of their suffix
(i.e.: nexttowardf and nexttoward have a long double argument, besides
the float and double arguments).

On top of that, they are also special because nexttoward functions are
not part of the _FloatN API, hence __nexttowardf128 do not exist.

This patch adds 4 new function implementations for the new long double
format:

  __nexttoward_to_ieee128
  __nexttowardf_to_ieee128
  __nexttowardieee128 (as an alias to __nextafterieee128)

Likewise, rename "long double" "_Float128" in shared ldbl-128
files to ensure correct type is used irrespective of ABI
switches.

Thank you to those who helped out with this patch:

Co-Authored-By: Tulio Magno Quites Machado Filho  <tuliom@linux.ibm.com>
2020-02-28 08:20:02 -06:00
Tulio Magno Quites Machado Filho
f1a166b0d7 ldbl-128ibm-compat: Provide a significand implementation
Reuse the template in order to provide the global symbol
__significandieee128.
2020-02-28 08:20:02 -06:00
Tulio Magno Quites Machado Filho
8dbfea3a20 ldbl-128ibm-compat: Redirect long double functions to f128/ieee128 functions
Modify the headers to redirect long double functions to global __*f128
symbols or to __*ieee128 otherwise.

Most of the functions in math.h benefit from the infrastructure already
available for __LDBL_COMPAT.  The only exceptions are nexttowardf and
nexttoward that need especial treatment.

Both math/bits/mathcalls-helper-functions.h and math/bits/mathcalls.h
were modified in order to provide alternative redirection destinations
that are essential to support functions that should not be redirected to
the same name pattern of the rest of the functions, i.e.: __fpclassify,
__signbit, __iseqsig, __issignaling, isinf, finite and isnan, which will
be redirected to __*f128 instead of __*ieee128 used for the rest.
2020-02-28 08:20:02 -06:00
Paul E. Murphy
a92ac9f494 ldbl-128ibm: make ieee754.h work with IEEE 128 long double
Instead of attempting something more creative, just copy
the small struct from ldbl-128 and enable it when IEEE
long double is present, and update the ibm long double
variant if supported.

Likewise, provide a shadow copy of math_ldbl.h to prevent
the ibm128 specific long double header from poisoning
unrelated files due to it's usage in math_private.h.

Reviewed-by: Tulio Magno Quites Machado Filho <tuliom@linux.ibm.com>
2020-02-21 08:24:36 -06:00
Paul E. Murphy
bb05281822 ldbl-128ibm-compat: fixup subdir location of several funcs
We want to ensure that if a second file is built to support
ieee128 long double, we built its companion implementation
with ibm128 long double.  The shared object versions of these
files build correctly because the aliasing is sufficiently
complex to prevent the redirects from applying when defining
them.

However, this does not prevent the static object variants
from becoming quietly broken due to redirects.  This is
intentionally avoided by marking such objects to be built
with -mabi=ibmlongdouble.

Shuffle the misplaced routines to build against the subdir
which defines the needed symbols.
2020-02-21 08:24:15 -06:00
Paul E. Murphy
74f614b70a ldbl-128ibm-compat: enforce correct abi flags on internal file
A number of utility files and helper objects should also be
explicitly configured to build with the ibm128 ABI to prevent
gremlins when enabling IEEE long double.
2020-02-21 08:24:04 -06:00
Tulio Magno Quites Machado Filho
218dad29e8 ldbl-128ibm-compat: Provide ieee128 symbols to narrow functions
Move the narrow math aliasing macros into a new sysdep header file
math-narrow-alias-float128.h.  Then, provide an override header
to supply the necessary changes to supply the *ieee128 aliases of
these symbols.

This adds ieee128 aliases for faddl, fdivl, fmull, fsubl, daddl, ddivl,
dmull, dsubl.
2020-02-20 17:12:14 -06:00
Tulio Magno Quites Machado Filho
77ad97356c Undefine redirections after long double definition on __LDBL_COMPAT [BZ #23294]
After defining the long double redirections to double, __MATHDECL_1 has
to be redefined to its previous state in order to avoid redirecting all
subsequent types.
2020-02-20 17:11:06 -06:00
Tulio Magno Quites Machado Filho
2b7dc93f82 ldbl-128ibm-compat: Provide a scalb implementation
Reuse the template in order to provide the redirect for
scalbl to __scalbieee128, but avoid any extra aliasing
as this is intended to support long double redirects only.
2020-02-14 08:24:56 -06:00
Tulio Magno Quites Machado Filho
c624d23260 Add a generic scalb implementation
This is a preparatory patch to enable building a _Float128
variant to ease reuse when building a _Float128 variant to
alias this long double only symbol.

Notably, stubs are added where missing to the native _Float128
sysdep dir to prevent building these newly templated variants
created inside the build directories.

Also noteworthy are the changes around LIBM_SVID_COMPAT.  These
changes are not intuitive.  The templated version is only
enabled when !LIBM_SVID_COMPAT, and the compat version is
predicated entirely on LIBM_SVID_COMPAT.  Thus, exactly one is
stubbed out entirely when building.  The nldbl scalb compat
files are updated to account for this.

Likewise, fixup the reuse of m68k's e_scalb{f,l}.c to include
it's override of e_scalb.c.  Otherwise, the search path finds
the templated copy in the build directory.  This could be
futher simplified by providing an overridden template, but I
lack the hardware to verify.
2020-02-14 08:24:56 -06:00
Joseph Myers
ad180676b8 Adjust thresholds in Bessel function implementations (bug 14469).
A recent discussion in bug 14469 notes that a threshold in float
Bessel function implementations, used to determine when to use a
simpler implementation approach, results in substantially inaccurate
results.

As I discussed in
<https://sourceware.org/ml/libc-alpha/2013-03/msg00345.html>, a
heuristic argument suggests 2^(S+P) as the right order of magnitude
for a suitable threshold, where S is the number of significand bits in
the floating-point type and P is the number of significant bits in the
representation of the floating-point type, and the float and ldbl-96
implementations use thresholds that are too small.  Some threshold
does need using, there or elsewhere in the implementation, to avoid
spurious underflow and overflow for large arguments.

This patch sets the thresholds in the affected implementations to more
heuristically justifiable values.  Results will still be inaccurate
close to zeroes of the functions (thus this patch does *not* fix any
of the bugs for Bessel function inaccuracy); fixing that would require
a different implementation approach, likely along the lines described
in <http://www.cl.cam.ac.uk/~jrh13/papers/bessel.ps.gz>.

So the justification for a change such as this would be statistical
rather than based on particular tests that had excessive errors and no
longer do so (no doubt such tests could be found, but would probably
be too fragile to add to the testsuite, as liable to give large errors
again from very small implementation changes or even from compiler
changes).  See
<https://sourceware.org/ml/libc-alpha/2020-02/msg00638.html> for such
statistics of the resulting improvements for float functions.

Tested (glibc testsuite) for x86_64.
2020-02-14 14:16:25 +00:00
Florian Weimer
c10acd4026 math/test-sinl-pseudo: Use stack protector only if available
This fixes commit 9333498794 ("Avoid ldbl-96 stack
corruption from range reduction of pseudo-zero (bug 25487).").
2020-02-13 17:01:15 +01:00
Joseph Myers
9333498794 Avoid ldbl-96 stack corruption from range reduction of pseudo-zero (bug 25487).
Bug 25487 reports stack corruption in ldbl-96 sinl on a pseudo-zero
argument (an representation where all the significand bits, including
the explicit high bit, are zero, but the exponent is not zero, which
is not a valid representation for the long double type).

Although this is not a valid long double representation, existing
practice in this area (see bug 4586, originally marked invalid but
subsequently fixed) is that we still seek to avoid invalid memory
accesses as a result, in case of programs that treat arbitrary binary
data as long double representations, although the invalid
representations of the ldbl-96 format do not need to be consistently
handled the same as any particular valid representation.

This patch makes the range reduction detect pseudo-zero and unnormal
representations that would otherwise go to __kernel_rem_pio2, and
returns a NaN for them instead of continuing with the range reduction
process.  (Pseudo-zero and unnormal representations whose unbiased
exponent is less than -1 have already been safely returned from the
function before this point without going through the rest of range
reduction.)  Pseudo-zero representations would previously result in
the value passed to __kernel_rem_pio2 being all-zero, which is
definitely unsafe; unnormal representations would previously result in
a value passed whose high bit is zero, which might well be unsafe
since that is not a form of input expected by __kernel_rem_pio2.

Tested for x86_64.
2020-02-12 23:31:56 +00:00
Paul E. Murphy
cd7965bd97 ldbl-128ibm-compat: set PRINTF_CHK flag in {,v}sprintf_chk
This should be unconditionally set to match the common implementation,
and fixes multiple test failures related to sprintf.

Reviewed-by: Tulio Magno Quites Machado Filho <tuliom@linux.ibm.com>
2020-02-10 16:51:03 -06:00
Wilco Dijkstra
814309f0c2 Remove a comment claiming that sin/cos round correctly. 2020-02-07 17:15:37 +00:00
Matheus Castanho
9f8b135f76 Fix maybe-uninitialized error on powerpc
The build has been failing on powerpc64le-linux-gnu with GCC 10
due to a maybe-uninitialized error:

../sysdeps/ieee754/dbl-64/mpa.c:875:6: error: ‘w.e’ may be used
uninitialized in this function [-Werror=maybe-uninitialized]
  875 |   EY -= EX;
      |      ^~

The warning is thrown because when __inv is called by __dvd *y is not
initialized and if t == 0 before calling __dbl_mp, EY will stay
uninitialized, as the function does not touch it in this case.

However, since t will be set to 1/t before calling __dbl_mp, t == 0 will
never happen, so we can instruct the compiler to ignore this case, which
suppresses the warning.

Tested on powerpc64le.

Reviewed-by: Tulio Magno Quites Machado Filho <tuliom@linux.ibm.com>
2020-01-17 21:02:13 -03:00
Wilco Dijkstra
220622dde5 Add libm_alias_finite for _finite symbols
This patch adds a new macro, libm_alias_finite, to define all _finite
symbol.  It sets all _finite symbol as compat symbol based on its first
version (obtained from the definition at built generated first-versions.h).

The <fn>f128_finite symbols were introduced in GLIBC 2.26 and so need
special treatment in code that is shared between long double and float128.
It is done by adding a list, similar to internal symbol redifinition,
on sysdeps/ieee754/float128/float128_private.h.

Alpha also needs some tricky changes to ensure we still emit 2 compat
symbols for sqrt(f).

Passes buildmanyglibc.

Co-authored-by: Adhemerval Zanella <adhemerval.zanella@linaro.org>
Reviewed-by: Siddhesh Poyarekar <siddhesh@sourceware.org>
2020-01-03 10:02:04 -03:00
Joseph Myers
d614a75396 Update copyright dates with scripts/update-copyrights. 2020-01-01 00:14:33 +00:00
Gabriel F. T. Gomes
9ae967bf45 ldbl-128ibm-compat: Do not mix -mabi=*longdouble and -mlong-double-128
Some compiler versions, e.g. GCC 7, complain when -mlong-double-128 is
used together with -mabi=ibmlongdouble or -mabi=ieeelongdouble,
producing the following error message:

  cc1: error: ‘-mabi=ibmlongdouble’ requires ‘-mlong-double-128’

This patch removes -mlong-double-128 from the compilation lines that
explicitly request -mabi=*longdouble.

Tested for powerpc64le.

Reviewed-by: Tulio Magno Quites Machado Filho <tuliom@linux.ibm.com>
2019-12-27 15:02:10 -03:00
Tulio Magno Quites Machado Filho
5d73c96f64 ldbl-128ibm-compat: Compiler flags for stdio functions
Some of the files that provide stdio.h and wchar.h functions have a
filename prefixed with 'io', such as 'iovsprintf.c'.  On platforms that
imply ldbl-128ibm-compat, these files must be compiled with the flag
-mabi=ibmlongdouble.  This patch adds this flag to their compilation.

Notice that this is not required for the other files that provide
similar functions, because filenames that are not prefixed with 'io'
have ldbl-128ibm-compat counterparts in the Makefile, which already adds
-mabi=ibmlongdouble to them.

Reviewed-by: Gabriel F. T. Gomes <gabrielftg@linux.ibm.com>
2019-12-27 15:02:10 -03:00
Tulio Magno Quites Machado Filho
1ef9b6e0bf Do not redirect calls to __GI_* symbols, when redirecting to *ieee128
On platforms where long double has IEEE binary128 format as a third
option (initially, only powerpc64le), many exported functions are
redirected to their __*ieee128 equivalents.  This redirection is
provided by installed headers such as stdio-ldbl.h, and is supposed to
work correctly with user code.

However, during the build of glibc, similar redirections are employed,
in internal headers, such as include/stdio.h, in order to avoid extra
PLT entries.  These redirections conflict with the redirections to
__*ieee128, and must be avoided during the build.  This patch protects
the second redirections with a test for __LONG_DOUBLE_USES_FLOAT128, a
new macro that is defined to 1 when functions that deal with long double
typed values reuses the _Float128 implementation (this is currently only
true for powerpc64le).

Tested for powerpc64le, x86_64, and with build-many-glibcs.py.

Co-authored-by: Gabriel F. T. Gomes <gabrielftg@linux.ibm.com>
Reviewed-by: Florian Weimer <fweimer@redhat.com>
2019-12-27 15:02:10 -03:00
Gabriel F. T. Gomes
f8cd102081 Avoid compat symbols for totalorder in powerpc64le IEEE long double
On powerpc64le, the libm_alias_float128_other_r_ldbl macro is
used to create an alias between totalorderf128 and __totalorderlieee128,
as well as between the totalordermagf128 and __totalordermaglieee128.

However, the totalorder* and totalordermag* functions changed their
parameter type since commit ID 42760d7646 and got compat symbols for
their old versions.  With this change, the aforementioned macro would
create two conflicting aliases for __totalorderlieee128 and
__totalordermaglieee128.

This patch avoids the creation of the alias between the IEEE long double
symbols (__totalorderl*ieee128) and the compat symbols, because the IEEE
long double functions have never been exported thus don't need such
compat symbol.

Tested for powerpc64le.

Reviewed-by: Joseph Myers <joseph@codesourcery.com>
2019-12-23 16:32:20 -03:00
Gabriel F. T. Gomes
3021e78178 ldbl-128ibm-compat: Add *cvt functions
This patch adds IEEE long double versions of q*cvt* functions for
powerpc64le.  Unlike all other long double to/from string conversion
functions, these do not rely on internal functions that can take
floating-point numbers with different formats and act on them
accordingly, instead, the related files are rebuilt with the
-mabi=ieeelongdouble compiler flag set.

Having -mabi=ieeelongdouble passed to the compiler causes the object
files to be marked with a .gnu_attribute that is incompatible with the
.gnu_attribute in files built with -mabi=ibmlongdouble (the default).
The difference causes error messages similar to the following:

  ld: libc_pic.a(s_isinfl.os) uses IBM long double,
      libc_pic.a(ieee128-qefgcvt_r.os) uses IEEE long double.
  collect2: error: ld returned 1 exit status
  make[2]: *** [../Makerules:649: libc_pic.os] Error 1

Although this warning is useful in other situations, the library
actually needs to have functions with different long double formats, so
.gnu_attribute generation is explicitly disabled for these files with
the use of -mno-gnu-attribute.

Tested for powerpc64le on the branch that actually enables the
sysdeps/ieee754/ldbl-128ibm-compat for powerpc64le.

Reviewed-by: Paul E. Murphy <murphyp@linux.ibm.com>
2019-12-23 16:32:20 -03:00
Gabriel F. T. Gomes
f1a0eb5b67 ldbl-128ibm-compat: Add ISO C99 versions of scanf functions
In the format string for *scanf functions, the '%as', '%aS', and '%a[]'
modifiers behave differently depending on ISO C99 compatibility.  When
_GNU_SOURCE is defined and -std=c89 is passed to the compiler, these
functions behave like ascanf, and the modifiers allocate memory for the
output.  Otherwise, the ISO C99 compliant version of these functions is
used, and the modifiers consume a floating-point argument.  This patch
adds the IEEE binary128 variant of ISO C99 compliant functions for the
third long double format on powerpc64le.

Tested for powerpc64le.

Reviewed-by: Paul E. Murphy <murphyp@linux.ibm.com>
2019-12-13 18:02:09 -03:00
Gabriel F. T. Gomes
348787f069 ldbl-128ibm-compat: Fix selection of GNU and ISO C99 scanf
Since commit

commit 03992356e6
Author: Zack Weinberg <zackw@panix.com>
Date:   Sat Feb 10 11:58:35 2018 -0500

    Use C99-compliant scanf under _GNU_SOURCE with modern compilers.

the selection of the GNU versions of scanf functions requires both
_GNU_SOURCE and -std=c89.  This patch changes the tests in
ldbl-128ibm-compat so that they actually test the GNU versions (without
this change, the redirection to the ISO C99 version always happens, so
GNU versions of the new implementation (e.g. __scanfieee128) were left
untested).

Tested for powerpc64le.

Reviewed-by: Paul E. Murphy <murphyp@linux.ibm.com>
2019-12-13 18:01:25 -03:00
Stefan Liebler
1902d5d5ff Adjust s_copysignl.c regarding code style.
This patch just adjusts the generic implementation regarding code style.
No functional change.

Reviewed-by: Adhemerval Zanella <adhemerval.zanella@linaro.org>
2019-12-11 15:12:21 +01:00
Stefan Liebler
171d23d7cb Adjust s_ceilf.c and s_ceill.c regarding code style.
This patch just adjusts the generic implementation regarding code style.
No functional change.

Reviewed-by: Adhemerval Zanella <adhemerval.zanella@linaro.org>
2019-12-11 15:12:21 +01:00
Stefan Liebler
d3a0409ab6 Adjust s_floorf.c and s_floorl.c regarding code style.
This patch just adjusts the generic implementation regarding code style.
No functional change.

Reviewed-by: Adhemerval Zanella <adhemerval.zanella@linaro.org>
2019-12-11 15:12:20 +01:00
Stefan Liebler
99b39a83e7 Adjust s_rintf.c and s_rintl.c regarding code style.
This patch just adjusts the generic implementation regarding code style.
No functional change.

Reviewed-by: Adhemerval Zanella <adhemerval.zanella@linaro.org>
2019-12-11 15:12:20 +01:00
Stefan Liebler
6a3866dae9 Adjust s_nearbyintf.c and s_nearbyintl.c regarding code style.
This patch just adjusts the generic implementation regarding code style.
No functional change.

Reviewed-by: Adhemerval Zanella <adhemerval.zanella@linaro.org>
2019-12-11 15:12:19 +01:00
Stefan Liebler
f818afdd3b Use GCC builtins for copysign functions if desired.
This patch is always using the corresponding GCC builtin for copysignf, copysign,
and is using the builtin for copysignl, copysignf128 if the USE_FUNCTION_BUILTIN
macros are defined to one in math-use-builtins.h.

Altough the long double version is enabled by default we still need
the macro and the alternative implementation as the _Float128 version
of the builtin is not available with all supported GCC versions.

Reviewed-by: Adhemerval Zanella <adhemerval.zanella@linaro.org>
2019-12-11 15:12:19 +01:00
Stefan Liebler
f82996f815 Use GCC builtins for round functions if desired.
This patch is using the corresponding GCC builtin for roundf, round,
roundl and roundf128 if the USE_FUNCTION_BUILTIN macros are defined to one
in math-use-builtins.h.

This is the case for s390 if build with at least --march=z196 --mzarch.
Otherwise the generic implementation is used.  The code of the generic
implementation is not changed.

Reviewed-by: Adhemerval Zanella <adhemerval.zanella@linaro.org>
2019-12-11 15:12:18 +01:00
Stefan Liebler
1ac9c1cf87 Use GCC builtins for trunc functions if desired.
This patch is using the corresponding GCC builtin for truncf, trunc,
truncl and truncf128 if the USE_FUNCTION_BUILTIN macros are defined to one
in math-use-builtins.h.

This is the case for s390 if build with at least --march=z196 --mzarch.
Otherwise the generic implementation is used.  The code of the generic
implementation is not changed.

Reviewed-by: Adhemerval Zanella <adhemerval.zanella@linaro.org>
2019-12-11 15:12:17 +01:00
Stefan Liebler
62560ee840 Use GCC builtins for ceil functions if desired.
This patch is using the corresponding GCC builtin for ceilf, ceil,
ceill and ceilf128 if the USE_FUNCTION_BUILTIN macros are defined to one
in math-use-builtins.h.

This is the case for s390 if build with at least --march=z196 --mzarch.
Otherwise the generic implementation is used.  The code of the generic
implementation is not changed.

Reviewed-by: Adhemerval Zanella <adhemerval.zanella@linaro.org>
2019-12-11 15:12:17 +01:00
Stefan Liebler
6c1b6a5e8c Use GCC builtins for floor functions if desired.
This patch is using the corresponding GCC builtin for floorf, floor,
floorl and floorf128 if the USE_FUNCTION_BUILTIN macros are defined to one
in math-use-builtins.h.

This is the case for s390 if build with at least --march=z196 --mzarch.
Otherwise the generic implementation is used.  The code of the generic
implementation is not changed.

Reviewed-by: Adhemerval Zanella <adhemerval.zanella@linaro.org>
2019-12-11 15:12:16 +01:00
Stefan Liebler
a2a9b00429 Use GCC builtins for rint functions if desired.
This patch is using the corresponding GCC builtin for rintf, rint,
rintl and rintf128 if the USE_FUNCTION_BUILTIN macros are defined to one
in math-use-builtins.h.

This is the case for s390 if build with at least --march=z196 --mzarch.
Otherwise the generic implementation is used.  The code of the generic
implementation is not changed.

Reviewed-by: Adhemerval Zanella <adhemerval.zanella@linaro.org>
2019-12-11 15:12:16 +01:00
Stefan Liebler
ae3577f607 Use GCC builtins for nearbyint functions if desired.
This patch is using the corresponding GCC builtin for nearbyintf, nearbyint,
nearbintl and nearbyintf128 if the USE_FUNCTION_BUILTIN macros are defined to one
in math-use-builtins.h.

This is the case for s390 if build with at least --march=z196 --mzarch.
Otherwise the generic implementation is used.  The code of the generic
implementation is not changed.

Reviewed-by: Adhemerval Zanella <adhemerval.zanella@linaro.org>
2019-12-11 15:12:15 +01:00
Stefan Liebler
36e9acbd5c Always use wordsize-64 version of s_round.c.
This patch replaces s_round.c in sysdeps/dbl-64 with the one in
sysdeps/dbl-64/wordsize-64 and removes the latter one.

Reviewed-by: Adhemerval Zanella <adhemerval.zanella@linaro.org>
2019-12-11 15:12:14 +01:00
Stefan Liebler
1c94bf0f0a Always use wordsize-64 version of s_trunc.c.
This patch replaces s_trunc.c in sysdeps/dbl-64 with the one in
sysdeps/dbl-64/wordsize-64 and removes the latter one.
The code is not changed except changes in code style.

Also adjusted the include path in x86_64 and sparc64 files.

Reviewed-by: Adhemerval Zanella <adhemerval.zanella@linaro.org>
2019-12-11 15:12:14 +01:00
Stefan Liebler
9f234eafe8 Always use wordsize-64 version of s_ceil.c.
This patch replaces s_ceil.c in sysdeps/dbl-64 with the one in
sysdeps/dbl-64/wordsize-64 and removes the latter one.
The code is not changed except changes in code style.

Also adjusted the include path in x86_64 and sparc64 files.

Reviewed-by: Adhemerval Zanella <adhemerval.zanella@linaro.org>
2019-12-11 15:12:13 +01:00