Changes with respect to v1:
- added comment in e_j1f.c to explain the use of float is enough
Reviewed-by: Adhemerval Zanella <adhemerval.zanella@linaro.org>
Clang issues:
../sysdeps/ieee754/dbl-64/s_llround.c:83:30: error: incompatible
redeclaration of library function 'lround'
[-Werror,-Wincompatible-library-redeclaration]
libm_alias_double (__lround, lround)
^
../sysdeps/ieee754/dbl-64/s_llround.c:83:30: note: 'lround' is a builtin
with type 'long (double)'
Reviewed-by: Sam James <sam@gentoo.org>
clang issues:
../sysdeps/ieee754/dbl-64/e_lgamma_r.c:234:29: error: absolute value function 'fabsf'
given an argument of type 'double' but has parameter of type 'float' which may cause \
truncation of value [-Werror,-Wabsolute-value]
It should not matter because the value is 0.0, but using fabs is
simpler than adding a warning suppresion.
Reviewed-by: Sam James <sam@gentoo.org>
And remove some unused entries of the fallback table.
Checked on x86_64-linux-gnu and aarch64-linux-gnu.
Reviewed-by: Wilco Dijkstra <Wilco.Dijkstra@arm.com>
The fma is required only for x == -0x1.da285cp-5 in FE_TONEAREST
to provide correctly rounded results.
Checked on x86_64-linux-gnu and i686-linux-gnu.
Reviewed-by: Wilco Dijkstra <Wilco.Dijkstra@arm.com>
The fma is required only for x == +/-0x1.6371e8p-4f in FE_TOWARDZERO
to provide correctly rounded results.
Checked on x86_64-linux-gnu and aarch64-linux-gnu.
Reviewed-by: Wilco Dijkstra <Wilco.Dijkstra@arm.com>
The fma is not required to provide correctly rounded and it helps
on !__FP_FAST_FMA ISAs.
Checked on x86_64-linux-gnu and i686-linux-gnu.
Reviewed-by: Paul Zimmermann <Paul.Zimmermann@inria.fr>
The fma is required only for inputs less than 0x1.0fd288p-127. Also
only add the extra check for !__FP_FAST_FMA targets.
Checked on x86_64-linux-gnu and aarch64-linux-gnu.
Reviewed-by: Paul Zimmermann <Paul.Zimmermann@inria.fr>
The fma is not strickly required to provide correctly rounded and
it helps on !__FP_FAST_FMA ABIs.
Checked on x86_64-linux-gnu and i686-linux-gnu.
Reviewed-by: Paul Zimmermann <Paul.Zimmermann@inria.fr>
With same micro-optimization done for the double variant:
* Combine the |y| zero check.
* Rework the check to adjust result and call fmod.
* Remove one check after fmod.
* Remove float-int-float roundtrip on return.
Also use math_config.h macros and indent the code. The resulting
strategy is different in many places that I think requires a
different Copyright.
I see the following performance improvements using remainder benchtests
(using reciprocal-throughput metric):
Architecture | Input | master | patch | Improvemnt
-----------------|-----------------|----------|-----------------------
x86_64 | subnormals | 20.4176 | 19.6144 | 3.93%
x86_64 | normal | 54.0939 | 52.2343 | 3.44%
x86_64 | close-exponent | 23.9120 | 22.3768 | 6.42%
aarch64 | subnormals | 9.2423 | 8.3825 | 9.30%
aarch64 | normal | 30.5393 | 29.244 | 4.24%
aarch64 | close-exponent | 15.5405 | 13.9256 | 10.39%
The aarch64 used as Neoverse-N1, gcc 15.1.1; while the x86_64 was
a AMD Ryzen 9 5900X, gcc 15.2.1.
Checked on x86_64-linux-gnu and aarch64-linux-gnu.
Reviewed-by: Wilco Dijkstra <Wilco.Dijkstra@arm.com>
The commit 34b9f8bc17 provides an optimized fmod implementation; use
the same strategy used for remainderf and implement the double variant
on top of fmod.
I see the following performance improvements using remainder benchtests
(using reciprocal-throughput metric):
Architecture | Input | master | patch | Improvemnt
-----------------|-----------------|----------|-----------------------
x86_64 | subnormals | 76.1345 | 21.5334 | 71.72%
x86_64 | normal | 553.2670 | 426.5670 | 22.90%
x86_64 | close-exponent | 30.5111 | 22.6893 | 25.64%
aarch64 | subnormals | 26.0734 | 8.4876 | 67.45%
aarch64 | normal | 205.2590 | 200.082 | 2.52%
aarch64 | close-exponent | 13.8481 | 13.6663 | 1.31%
The aarch64 used as Neoverse-N1, gcc 15.1.1; while the x86_64 was
a AMD Ryzen 9 5900X, gcc 15.2.1.
This implementation also fixes the math/test-double-remainder issues
on alpha.
Tested on aarch64-linux-gnu and x86_64-linux-gnu.
Reviewed-by: Wilco Dijkstra <Wilco.Dijkstra@arm.com>
When compiling on x86_64 with -Wshift-overflow=2 you can see the
following warning:
../sysdeps/ieee754/flt-32/math_config.h: In function ‘is_inf’:
../sysdeps/ieee754/flt-32/math_config.h:184:37: warning: result of ‘2139095040 << 1’ requires 33 bits to represent, but ‘int’ only has 32 bits [-Wshift-overflow=]
184 | return (x << 1) == (EXPONENT_MASK << 1);
| ^~
This patch adjusts the definitions to use UINT32_C. This matches the
definitions in sysdeps/ieee754/dbl-64/math_config.h which use UINT64_C
for these definitions.
Reviewed-by: H.J. Lu <hjl.tools@gmail.com>
To avoid linknamespace issues on old standards. It is required
if the fallback fma implementation is used if/when it is also
used internally for other implementation.
Reviewed-by: DJ Delorie <dj@redhat.com>
To avoid linknamespace issues on old standards. It is required
if the fallback fma implementation is used if/when it is also
used internally for other implementation.
Reviewed-by: DJ Delorie <dj@redhat.com>
To avoid linknamespace issues on old standards. It is required
if the fallback fma implementation is used if/when it is also
used internally for other implementation.
Reviewed-by: DJ Delorie <dj@redhat.com>
To avoid linknamespace issues on old standards. It is required
if the fallback fma implementation is used if/when it is also
used internally for other implementation.
Reviewed-by: DJ Delorie <dj@redhat.com>
To avoid linknamespace issues on old standards. It is required
if the fallback fma implementation is used if/when it is also
used internally for other implementation.
Reviewed-by: DJ Delorie <dj@redhat.com>
Convert 'compare_real', 'read_real', and 'verify_input' macros to
functions so as to improve readability and avoid pitfalls.
Reviewed-by: Adhemerval Zanella <adhemerval.zanella@linaro.org>
Move the inclusion of the data class header from the individual tests to
the data-type-specific skeleton, providing for the use of the data type
under test in the data class header and reducing duplication.
Reviewed-by: Adhemerval Zanella <adhemerval.zanella@linaro.org>
Update NaN input data with 'n-char-sequence' in reference data matching
data under test, removing test failures with the M68K host.
Reviewed-by: Adhemerval Zanella <adhemerval.zanella@linaro.org>
Verify that . input is rejected by 'f' conversion (and its uppercase
counterpart). Replace 0 input with .0 rather than adding new one,
because the integral part of 0 is already covered by 0.0 data, so
there's no need to keep this duplication.
Reviewed-by: Adhemerval Zanella <adhemerval.zanella@linaro.org>
Verify that . input is rejected by 'e' conversion (and its uppercase
counterpart). Replace 0e0 input with .0e0 rather than adding new one,
because 0 significand is already covered by 0e+0 data, so there's no
need to keep this duplication.
Reviewed-by: Adhemerval Zanella <adhemerval.zanella@linaro.org>
Verify that 0x. input is rejected by 'a' and 'g' conversions (and their
uppercase counterparts). Replace 0x0p0 input with 0x.0p0 rather than
adding new one, because 0x0 significand is already covered by 0x0p+0
data, so there's no need to keep this duplication.
Reviewed-by: Adhemerval Zanella <adhemerval.zanella@linaro.org>
The generic implementation is slight more optimized than the powerpc
one, where it has a more optimized inf/nan check (by not using FP
unit checks, along with branch prediction hints), and removed one
branch by issuing trunc instead of a combination of floor/ceil (which
also generated less code).
On power10 with gcc 14.2.1:
reciprocal-throughput master patch difference
workload-0_1 1.5210 1.3942 8.34%
workload-1_maxint 2.0926 1.3940 33.38%
workload-maxint_maxfloat 1.7851 1.3940 21.91%
workload-integral 1.5216 1.3941 8.37%
latency master patch difference
workload-0_1 1.5928 2.6337 -65.35%
workload-1_maxint 3.2929 2.6337 20.02%
workload-maxint_maxfloat 1.9697 2.6341 -33.73%
workload-integral 2.0597 2.6337 -27.87%
Checked on powerpc64le-linux-gnu.
Reviewed-by: Sachin Monga <smonga@linux.ibm.com>
Refactor the generic implementation to use math_config.h definitions,
and add an alternative one if the ABI supports truncf instructions
(gated through math-use-builtins-trunc.h).
The generic implementation generates similar code on x86_64, while
the optimization one for aarch64 (where truncf is supported as a
builtin by through frintz), the improvements are:
reciprocal-throughput master patch difference
workload-0_1 3.0595 3.0698 -0.34%
workload-1_maxint 5.1747 3.0542 40.98%
workload-maxint_maxfloat 3.4391 3.0349 11.75%
workload-integral 3.2732 3.0293 7.45%
latency master patch difference
workload-0_1 3.5267 4.7107 -33.57%
workload-1_maxint 6.9074 4.7282 31.55%
workload-maxint_maxfloat 3.7210 4.7506 -27.67%
workload-integral 3.8634 4.8137 -24.60%
Checked on aarch64-linux-gnu and x86_64-linux-gnu.
Reviewed-by: Wilco Dijkstra <Wilco.Dijkstra@arm.com>
Refactor the generic implementation to use math_config.h definitions,
and add an alternative one if the ABI supports truncf instructions
(gated through math-use-builtins-trunc.h).
The generic implementation generates similar code for x86_64, while
the optimization path aarch64 (where truncf is supported as a builtin)
through frintz), the improvements are:
reciprocal-throughput master patch difference
workload-0_1 3.0740 3.0326 1.35%
workload-1_maxint 5.2231 3.0436 41.73%
workload-maxint_maxfloat 4.0962 3.0551 25.42%
workload-integral 3.7093 3.0612 17.47%
latency master patch difference
workload-0_1 3.5521 4.7313 -33.20%
workload-1_maxint 6.7148 4.7314 29.54%
workload-maxint_maxfloat 4.0458 4.7518 -17.45%
workload-integral 3.9719 4.7427 -19.40%
Checked on aarch64-linux-gnu and x86_64-linux-gnu.
Reviewed-by: Wilco Dijkstra <Wilco.Dijkstra@arm.com>
It removes the wrapper by moving the error/EDOM handling to an
out-of-line implementation (__math_invalidf_i/__math_invalidf_li).
Also, __glibc_unlikely is used on errors case since it helps
code generation on recent gcc.
The code now builds to with gcc-14 on aarch64:
0000000000000000 <__ilogbf>:
0: 1e260000 fmov w0, s0
4: d3577801 ubfx x1, x0, #23, #8
8: 340000e1 cbz w1, 24 <__ilogbf+0x24>
c: 5101fc20 sub w0, w1, #0x7f
10: 7103fc3f cmp w1, #0xff
14: 54000040 b.eq 1c <__ilogbf+0x1c> // b.none
18: d65f03c0 ret
1c: 12b00000 mov w0, #0x7fffffff // #2147483647
20: 14000000 b 0 <__math_invalidf_i>
24: 53175800 lsl w0, w0, #9
28: 340000a0 cbz w0, 3c <__ilogbf+0x3c>
2c: 5ac01000 clz w0, w0
30: 12800fc1 mov w1, #0xffffff81 // #-127
34: 4b000020 sub w0, w1, w0
38: d65f03c0 ret
3c: 320107e0 mov w0, #0x80000001 // #-2147483647
40: 14000000 b 0 <__math_invalidf_i>
Some ABI requires additional adjustments:
* i386 and m68k requires to use the template version, since
both provide __ieee754_ilogb implementatations.
* loongarch uses a custom implementation as well.
* powerpc64le also has a custom implementation for POWER9, which
is also used for float and float128 version. The generic
e_ilogb.c implementation is moved on powerpc to keep the
current code as-is.
Checked on aarch64-linux-gnu and x86_64-linux-gnu.
Reviewed-by: Wilco Dijkstra <Wilco.Dijkstra@arm.com>
It removes the wrapper by moving the error/EDOM handling to an
out-of-line implementation (__math_invalid_i/__math_invalid_li).
Also, __glibc_unlikely is used on errors case since it helps
code generation on recent gcc.
The code now builds to with gcc-14 on aarch64:
0000000000000000 <__ilogb>:
0: 9e660000 fmov x0, d0
4: d374f801 ubfx x1, x0, #52, #11
8: 340000e1 cbz w1, 24 <__ilogb+0x24>
c: 510ffc20 sub w0, w1, #0x3ff
10: 711ffc3f cmp w1, #0x7ff
14: 54000040 b.eq 1c <__ilogb+0x1c> // b.none
18: d65f03c0 ret
1c: 12b00000 mov w0, #0x7fffffff // #2147483647
20: 14000000 b 0 <__math_invalid_i>
24: d374cc00 lsl x0, x0, #12
28: b40000a0 cbz x0, 3c <__ilogb+0x3c>
2c: dac01000 clz x0, x0
30: 12807fc1 mov w1, #0xfffffc01 // #-1023
34: 4b000020 sub w0, w1, w0
38: d65f03c0 ret
3c: 320107e0 mov w0, #0x80000001 // #-2147483647
40: 14000000 b 0 <__math_invalid_i>
Some ABI requires additional adjustments:
* i386 and m68k requires to use the template version, since
both provide __ieee754_ilogb implementatations.
* loongarch uses a custom implementation as well.
* powerpc64le also has a custom implementation for POWER9, which
is also used for float and float128 version. The generic
e_ilogb.c implementation is moved on powerpc to keep the
current code as-is.
Checked on aarch64-linux-gnu and x86_64-linux-gnu.
Reviewed-by: Wilco Dijkstra <Wilco.Dijkstra@arm.com>
C23 adds various <math.h> function families originally defined in TS
18661-4. Add the rootn functions, which compute the Yth root of X for
integer Y (with a domain error if Y is 0, even if X is a NaN). The
integer exponent has type long long int in C23; it was intmax_t in TS
18661-4, and as with other interfaces changed after their initial
appearance in the TS, I don't think we need to support the original
version of the interface.
As with pown and compoundn, I strongly encourage searching for worst
cases for ulps error for these implementations (necessarily
non-exhaustively, given the size of the input space). I also expect a
custom implementation for a given format could be much faster as well
as more accurate, although the implementation is simpler than those
for pown and compoundn.
This completes adding to glibc those TS 18661-4 functions (ignoring
DFP) that are included in C23. See
https://gcc.gnu.org/bugzilla/show_bug.cgi?id=118592 regarding the C23
mathematical functions (not just the TS 18661-4 ones) missing built-in
functions in GCC, where such functions might usefully be added.
Tested for x86_64 and x86, and with build-many-glibcs.py.
C23 adds various <math.h> function families originally defined in TS
18661-4. Add the compoundn functions, which compute (1+X) to the
power Y for integer Y (and X at least -1). The integer exponent has
type long long int in C23; it was intmax_t in TS 18661-4, and as with
other interfaces changed after their initial appearance in the TS, I
don't think we need to support the original version of the interface.
Note that these functions are "compoundn" with a trailing "n", *not*
"compound" (CORE-MATH has the wrong name, for example).
As with pown, I strongly encourage searching for worst cases for ulps
error for these implementations (necessarily non-exhaustively, given
the size of the input space). I also expect a custom implementation
for a given format could be much faster as well as more accurate (I
haven't tested or benchmarked the CORE-MATH implementation for
binary32); this is one of the more complicated and less efficient
functions to implement in a type-generic way.
As with exp2m1 and exp10m1, this showed up places where the
powerpc64le IFUNC setup is not as self-contained as one might hope (in
this case, without the changes specific to powerpc64le, there were
undefined references to __GI___expf128).
Tested for x86_64 and x86, and with build-many-glibcs.py.
The left shift overflows for 'int', use uint32_t instead. It syncs
with CORE-MATH commit bbfabd99.
Checked on aarch64-linux-gnu, x86_64-linux-gnu, and i686-linux-gnu.
Reviewed-by: Carlos O'Donell <carlos@redhat.com>
The left shift overflows for 'int', use uint64_t instead. It syncs
with CORE-MATH commit d0a2be200cbc1344d800d9ef0ebee9ad67dd3ad8.
Checked on aarch64-linux-gnu, x86_64-linux-gnu, and i686-linux-gnu.
Reviewed-by: Carlos O'Donell <carlos@redhat.com>
The left shift overflows for 'int', use uint32_t instead. It syncs
with CORE-MATH commit bbfabd993a71b049c210b0febfd06d18369fadc1.
Checked on aarch64-linux-gnu, x86_64-linux-gnu, and i686-linux-gnu.
Reviewed-by: Carlos O'Donell <carlos@redhat.com>
The left shift overflows for 'int64_t', use unsigned instead. It syncs
with CORE-MATH commit f7c7408d1749ec2859ea249495af699359ae559b.
Checked on aarch64-linux-gnu, x86_64-linux-gnu, and i686-linux-gnu.
Reviewed-by: Carlos O'Donell <carlos@redhat.com>
The left shift overflows for 'int', use uint64_t instead. It syncs
with CORE-MATH commit bbfabd99.
Checked on aarch64-linux-gnu, x86_64-linux-gnu, and i686-linux-gnu.
Reviewed-by: Carlos O'Donell <carlos@redhat.com>
The left shift overflows for 'int', use a literal instead. It syncs
with OPTIMIZED-ROUTINES commit 0f87f607b976820ef41fe64d004fe67dc7af8236.
Checked on aarch64-linux-gnu, x86_64-linux-gnu, and i686-linux-gnu.
Reviewed-by: Carlos O'Donell <carlos@redhat.com>
The left shift overflows for 'int', use uint64_t instead. It syncs
with CORE-MATH commit 4d6192d2.
Checked on aarch64-linux-gnu, x86_64-linux-gnu, and i686-linux-gnu.
Reviewed-by: Carlos O'Donell <carlos@redhat.com>
The left shift overflows for 'int', use unsigned instead. It syncs
with CORE-MATH commit 4d6192d2.
Checked on aarch64-linux-gnu, x86_64-linux-gnu, and i686-linux-gnu.
Reviewed-by: Carlos O'Donell <carlos@redhat.com>
As mentioned by the reporter in a pull request against gcc-mirror,
the THREEp96 constant in e_expl.c is incorrect, it is actually 0x3.p+94f128
rather than 0x3.p+96f128.
The algorithm uses that to compute the t2 integer (tval2), by whose
delta it adjusts the x+xl pair and then in the result uses the precomputed
exp value for that entry.
Using 0x3.p+94f128 rather than 0x3.p+96f128 results in tval2 sometimes
being one smaller, sometimes one larger than the desired value, thus can mean
the x+xl pair after adjustment will be larger in absolute value than it
should be.
DesWursters created a test program for this
https://github.com/DesWurstes/comparefloats
and his results were
total: 1135000000 not_equal: 4322 earlier_score: 674 later_score: 3648
I've modified this so with
https://sourceware.org/bugzilla/show_bug.cgi?id=32411#c3
so that it actually tests pseudo-random _Float128 values with range
(-16384.,16384) with strong bias on values larger than 0.0002 in absolute
value (so that tval1/tval2 aren't zero most of the time) and that gave
total: 10000000000 not_equal: 29861 earlier_score: 4606 later_score: 25255
So, in both cases, in most cases the change doesn't result in any differences,
and in those rare cases where does, about 85% have smaller ulp than without
the patch.
Additionally I've tried
https://sourceware.org/bugzilla/show_bug.cgi?id=32411#c4
and in 2 billion iterations it didn't find any case where x+xl after the
adjustments without this change would be smaller in absolute value compared
to x+xl after the adjustments with this change.
Reviewed-by: Joseph Myers <josmyers@redhat.com>
