about summary refs log tree commit diff
path: root/math
diff options
context:
space:
mode:
authorAdhemerval Zanella Netto <adhemerval.zanella@linaro.org>2023-03-20 13:01:17 -0300
committerAdhemerval Zanella <adhemerval.zanella@linaro.org>2023-04-03 16:45:18 -0300
commitcf9cf33199fdd6550920ad43f19ad8b2435fc0c6 (patch)
treee03b3fb7d42424e4b7095a01b334e8db8296c1e8 /math
parent34b9f8bc170810c44184ad57ecf1800587e752a6 (diff)
downloadglibc-cf9cf33199fdd6550920ad43f19ad8b2435fc0c6.tar.gz
glibc-cf9cf33199fdd6550920ad43f19ad8b2435fc0c6.tar.xz
glibc-cf9cf33199fdd6550920ad43f19ad8b2435fc0c6.zip
math: Improve fmodf
This uses a new algorithm similar to already proposed earlier [1].
With x = mx * 2^ex and y = my * 2^ey (mx, my, ex, ey being integers),
the simplest implementation is:

   mx * 2^ex == 2 * mx * 2^(ex - 1)

   while (ex > ey)
     {
       mx *= 2;
       --ex;
       mx %= my;
     }

With mx/my being mantissa of double floating pointer, on each step the
argument reduction can be improved 8 (which is sizeof of uint32_t minus
MANTISSA_WIDTH plus the signal bit):

   while (ex > ey)
     {
       mx << 8;
       ex -= 8;
       mx %= my;
     }  */

The implementation uses builtin clz and ctz, along with shifts to
convert hx/hy back to doubles.  Different than the original patch,
this path assume modulo/divide operation is slow, so use multiplication
with invert values.

I see the following performance improvements using fmod benchtests
(result only show the 'mean' result):

  Architecture     | Input           | master   | patch
  -----------------|-----------------|----------|--------
  x86_64 (Ryzen 9) | subnormals      | 17.2549  | 12.0318
  x86_64 (Ryzen 9) | normal          | 85.4096  | 49.9641
  x86_64 (Ryzen 9) | close-exponents | 19.1072  | 15.8224
  aarch64 (N1)     | subnormal       | 10.2182  | 6.81778
  aarch64 (N1)     | normal          | 60.0616  | 20.3667
  aarch64 (N1)     | close-exponents | 11.5256  | 8.39685

I also see similar improvements on arm-linux-gnueabihf when running on
the N1 aarch64 chips, where it a lot of soft-fp implementation (for
modulo, and multiplication):

  Architecture     | Input           | master   | patch
  -----------------|-----------------|----------|--------
  armhf (N1)       | subnormal       | 11.6662  | 10.8955
  armhf (N1)       | normal          | 69.2759  | 34.1524
  armhf (N1)       | close-exponents | 13.6472  | 18.2131

Instead of using the math_private.h definitions, I used the
math_config.h instead which is used on newer math implementations.

Co-authored-by: kirill <kirill.okhotnikov@gmail.com>

[1] https://sourceware.org/pipermail/libc-alpha/2020-November/119794.html
Reviewed-by: Wilco Dijkstra  <Wilco.Dijkstra@arm.com>
Diffstat (limited to 'math')
-rw-r--r--math/libm-test-fmod.inc4
1 files changed, 4 insertions, 0 deletions
diff --git a/math/libm-test-fmod.inc b/math/libm-test-fmod.inc
index 8841c1332c..43376c3df2 100644
--- a/math/libm-test-fmod.inc
+++ b/math/libm-test-fmod.inc
@@ -213,6 +213,10 @@ static const struct test_ff_f_data fmod_test_data[] =
     TEST_ff_f (fmod, -0x1p127L, -0x3p-148L, -0x1p-147L, NO_INEXACT_EXCEPTION|ERRNO_UNCHANGED),
     TEST_ff_f (fmod, -0x1p127L, 0x3p-126L, -0x1p-125L, NO_INEXACT_EXCEPTION|ERRNO_UNCHANGED),
     TEST_ff_f (fmod, -0x1p127L, -0x3p-126L, -0x1p-125L, NO_INEXACT_EXCEPTION|ERRNO_UNCHANGED),
+    TEST_ff_f (fmod, 0x1.3a3e6p-127, 0x1.8b8338p-128, 0x1.d1f31p-129, NO_INEXACT_EXCEPTION|ERRNO_UNCHANGED),
+    TEST_ff_f (fmod, 0x1.3a3e6p-127, -0x1.8b8338p-128, 0x1.d1f31p-129, NO_INEXACT_EXCEPTION|ERRNO_UNCHANGED),
+    TEST_ff_f (fmod, -0x1.3a3e6p-127, 0x1.8b8338p-128, -0x1.d1f31p-129, NO_INEXACT_EXCEPTION|ERRNO_UNCHANGED),
+    TEST_ff_f (fmod, -0x1.3a3e6p-127, -0x1.8b8338p-128, -0x1.d1f31p-129, NO_INEXACT_EXCEPTION|ERRNO_UNCHANGED),
 #if !TEST_COND_binary32
     TEST_ff_f (fmod, 0x1p1023L, 0x3p-1074L, 0x1p-1073L, NO_INEXACT_EXCEPTION|ERRNO_UNCHANGED),
     TEST_ff_f (fmod, 0x1p1023L, -0x3p-1074L, 0x1p-1073L, NO_INEXACT_EXCEPTION|ERRNO_UNCHANGED),