[libc-commits] [libc] a752460 - [libc][math] Implement exp10f function correctly rounded to all rounding modes.
Tue Ly via libc-commits
libc-commits at lists.llvm.org
Mon Sep 19 07:02:03 PDT 2022
Author: Tue Ly
Date: 2022-09-19T10:01:40-04:00
New Revision: a752460d73912b8f6759d10c78c48071dc128b8c
URL: https://github.com/llvm/llvm-project/commit/a752460d73912b8f6759d10c78c48071dc128b8c
DIFF: https://github.com/llvm/llvm-project/commit/a752460d73912b8f6759d10c78c48071dc128b8c.diff
LOG: [libc][math] Implement exp10f function correctly rounded to all rounding modes.
Implement exp10f function correctly rounded to all rounding modes.
Algorithm: perform range reduction to reduce
```
10^x = 2^(hi + mid) * 10^lo
```
where:
```
hi is an integer,
0 <= mid * 2^5 < 2^5
-log10(2) / 2^6 <= lo <= log10(2) / 2^6
```
Then `2^mid` is stored in a table of 32 entries and the product `2^hi * 2^mid` is
performed by adding `hi` into the exponent field of `2^mid`.
`10^lo` is then approximated by a degree-5 minimax polynomials generated by Sollya with:
```
> P = fpminimax((10^x - 1)/x, 4, [|D...|], [-log10(2)/64. log10(2)/64]);
```
Performance benchmark using perf tool from the CORE-MATH project on Ryzen 1700:
```
$ CORE_MATH_PERF_MODE="rdtsc" ./perf.sh exp10f
GNU libc version: 2.35
GNU libc release: stable
CORE-MATH reciprocal throughput : 10.215
System LIBC reciprocal throughput : 7.944
LIBC reciprocal throughput : 38.538
LIBC reciprocal throughput : 12.175 (with `-msse4.2` flag)
LIBC reciprocal throughput : 9.862 (with `-mfma` flag)
$ CORE_MATH_PERF_MODE="rdtsc" ./perf.sh exp10f --latency
GNU libc version: 2.35
GNU libc release: stable
CORE-MATH latency : 40.744
System LIBC latency : 37.546
BEFORE
LIBC latency : 48.989
LIBC latency : 44.486 (with `-msse4.2` flag)
LIBC latency : 40.221 (with `-mfma` flag)
```
This patch relies on https://reviews.llvm.org/D134002
Reviewed By: orex, zimmermann6
Differential Revision: https://reviews.llvm.org/D134104
Added:
libc/src/math/exp10f.h
libc/src/math/generic/exp10f.cpp
libc/test/src/math/exhaustive/exp10f_test.cpp
libc/test/src/math/exp10f_test.cpp
Modified:
libc/config/darwin/arm/entrypoints.txt
libc/config/linux/aarch64/entrypoints.txt
libc/config/linux/x86_64/entrypoints.txt
libc/config/windows/entrypoints.txt
libc/docs/math.rst
libc/spec/gnu_ext.td
libc/src/math/CMakeLists.txt
libc/src/math/generic/CMakeLists.txt
libc/src/math/generic/explogxf.h
libc/test/src/math/CMakeLists.txt
libc/test/src/math/exhaustive/CMakeLists.txt
libc/utils/MPFRWrapper/MPFRUtils.cpp
libc/utils/MPFRWrapper/MPFRUtils.h
Removed:
################################################################################
diff --git a/libc/config/darwin/arm/entrypoints.txt b/libc/config/darwin/arm/entrypoints.txt
index 02331e15d541e..cd500736ac175 100644
--- a/libc/config/darwin/arm/entrypoints.txt
+++ b/libc/config/darwin/arm/entrypoints.txt
@@ -118,6 +118,7 @@ set(TARGET_LIBM_ENTRYPOINTS
libc.src.math.coshf
libc.src.math.cosf
libc.src.math.expf
+ libc.src.math.exp10f
libc.src.math.exp2f
libc.src.math.expm1f
libc.src.math.fabs
diff --git a/libc/config/linux/aarch64/entrypoints.txt b/libc/config/linux/aarch64/entrypoints.txt
index 9fc201b389a53..b506663a04b60 100644
--- a/libc/config/linux/aarch64/entrypoints.txt
+++ b/libc/config/linux/aarch64/entrypoints.txt
@@ -172,6 +172,7 @@ set(TARGET_LIBM_ENTRYPOINTS
libc.src.math.coshf
libc.src.math.cosf
libc.src.math.expf
+ libc.src.math.exp10f
libc.src.math.exp2f
libc.src.math.expm1f
libc.src.math.fabs
diff --git a/libc/config/linux/x86_64/entrypoints.txt b/libc/config/linux/x86_64/entrypoints.txt
index 8e85b56014d61..a0313d30c1921 100644
--- a/libc/config/linux/x86_64/entrypoints.txt
+++ b/libc/config/linux/x86_64/entrypoints.txt
@@ -177,6 +177,7 @@ set(TARGET_LIBM_ENTRYPOINTS
libc.src.math.coshf
libc.src.math.cosf
libc.src.math.expf
+ libc.src.math.exp10f
libc.src.math.exp2f
libc.src.math.expm1f
libc.src.math.fabs
diff --git a/libc/config/windows/entrypoints.txt b/libc/config/windows/entrypoints.txt
index 4bd46c5a3a502..3e6dceb4622e9 100644
--- a/libc/config/windows/entrypoints.txt
+++ b/libc/config/windows/entrypoints.txt
@@ -120,6 +120,7 @@ set(TARGET_LIBM_ENTRYPOINTS
libc.src.math.cosf
libc.src.math.coshf
libc.src.math.expf
+ libc.src.math.exp10f
libc.src.math.exp2f
libc.src.math.expm1f
libc.src.math.fabs
diff --git a/libc/docs/math.rst b/libc/docs/math.rst
index e30b99c14da68..faaae72f3cc59 100644
--- a/libc/docs/math.rst
+++ b/libc/docs/math.rst
@@ -129,6 +129,7 @@ cosh |check|
erf
erfc
exp |check|
+exp10 |check|
exp2 |check|
expm1 |check|
fma |check| |check|
@@ -161,6 +162,7 @@ atanh |check|
cos |check| large
cosh |check|
exp |check|
+exp10 |check|
exp2 |check|
expm1 |check|
fma |check| |check|
@@ -219,6 +221,8 @@ Performance
+--------------+-----------+-------------------+-----------+-------------------+-------------------------------------+------------+-------------------------+--------------+---------------+
| expf | 9 | 7 | 44 | 38 | :math:`[-10, 10]` | Ryzen 1700 | Ubuntu 20.04 LTS x86_64 | Clang 12.0.0 | FMA |
+--------------+-----------+-------------------+-----------+-------------------+-------------------------------------+------------+-------------------------+--------------+---------------+
+| exp10f | 10 | 8 | 40 | 38 | :math:`[-10, 10]` | Ryzen 1700 | Ubuntu 22.04 LTS x86_64 | Clang 14.0.0 | FMA |
++--------------+-----------+-------------------+-----------+-------------------+-------------------------------------+------------+-------------------------+--------------+---------------+
| exp2f | 9 | 6 | 35 | 31 | :math:`[-10, 10]` | Ryzen 1700 | Ubuntu 22.04 LTS x86_64 | Clang 14.0.0 | FMA |
+--------------+-----------+-------------------+-----------+-------------------+-------------------------------------+------------+-------------------------+--------------+---------------+
| expm1f | 9 | 44 | 42 | 121 | :math:`[-10, 10]` | Ryzen 1700 | Ubuntu 20.04 LTS x86_64 | Clang 12.0.0 | FMA |
diff --git a/libc/spec/gnu_ext.td b/libc/spec/gnu_ext.td
index d75d6e6a2c57d..7284e4d405fb5 100644
--- a/libc/spec/gnu_ext.td
+++ b/libc/spec/gnu_ext.td
@@ -25,6 +25,7 @@ def GnuExtensions : StandardSpec<"GNUExtensions"> {
RetValSpec<VoidType>,
[ArgSpec<FloatType>, ArgSpec<FloatPtr>, ArgSpec<FloatPtr>]
>,
+ FunctionSpec<"expf", RetValSpec<FloatType>, [ArgSpec<FloatType>]>,
]
>;
diff --git a/libc/src/math/CMakeLists.txt b/libc/src/math/CMakeLists.txt
index cb2a507d9fa03..e15e4c2a48e92 100644
--- a/libc/src/math/CMakeLists.txt
+++ b/libc/src/math/CMakeLists.txt
@@ -86,6 +86,8 @@ add_math_entrypoint_object(expf)
add_math_entrypoint_object(exp2f)
+add_math_entrypoint_object(exp10f)
+
add_math_entrypoint_object(expm1f)
add_math_entrypoint_object(fabs)
diff --git a/libc/src/math/exp10f.h b/libc/src/math/exp10f.h
new file mode 100644
index 0000000000000..7f0ef25efa48f
--- /dev/null
+++ b/libc/src/math/exp10f.h
@@ -0,0 +1,18 @@
+//===-- Implementation header for exp10f ------------------------*- C++ -*-===//
+//
+// Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions.
+// See https://llvm.org/LICENSE.txt for license information.
+// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
+//
+//===----------------------------------------------------------------------===//
+
+#ifndef LLVM_LIBC_SRC_MATH_EXP10F_H
+#define LLVM_LIBC_SRC_MATH_EXP10F_H
+
+namespace __llvm_libc {
+
+float exp10f(float x);
+
+} // namespace __llvm_libc
+
+#endif // LLVM_LIBC_SRC_MATH_EXP10F_H
diff --git a/libc/src/math/generic/CMakeLists.txt b/libc/src/math/generic/CMakeLists.txt
index a10f8792fd370..48a7101c0122c 100644
--- a/libc/src/math/generic/CMakeLists.txt
+++ b/libc/src/math/generic/CMakeLists.txt
@@ -562,6 +562,26 @@ add_entrypoint_object(
-O3
)
+add_entrypoint_object(
+ exp10f
+ SRCS
+ exp10f.cpp
+ HDRS
+ ../exp10f.h
+ DEPENDS
+ .explogxf
+ libc.src.__support.FPUtil.fenv_impl
+ libc.src.__support.FPUtil.fp_bits
+ libc.src.__support.FPUtil.multiply_add
+ libc.src.__support.FPUtil.nearest_integer
+ libc.src.__support.FPUtil.polyeval
+ libc.include.errno
+ libc.src.errno.errno
+ libc.include.math
+ COMPILE_OPTIONS
+ -O3
+)
+
add_entrypoint_object(
expm1f
SRCS
diff --git a/libc/src/math/generic/exp10f.cpp b/libc/src/math/generic/exp10f.cpp
new file mode 100644
index 0000000000000..a08b57f6893d3
--- /dev/null
+++ b/libc/src/math/generic/exp10f.cpp
@@ -0,0 +1,132 @@
+//===-- Single-precision 10^x function ------------------------------------===//
+//
+// Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions.
+// See https://llvm.org/LICENSE.txt for license information.
+// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
+//
+//===----------------------------------------------------------------------===//
+
+#include "src/math/exp10f.h"
+#include "explogxf.h"
+#include "src/__support/FPUtil/BasicOperations.h"
+#include "src/__support/FPUtil/FEnvImpl.h"
+#include "src/__support/FPUtil/FPBits.h"
+#include "src/__support/FPUtil/PolyEval.h"
+#include "src/__support/FPUtil/multiply_add.h"
+#include "src/__support/FPUtil/nearest_integer.h"
+#include "src/__support/common.h"
+
+#include <errno.h>
+
+namespace __llvm_libc {
+
+LLVM_LIBC_FUNCTION(float, exp10f, (float x)) {
+ using FPBits = typename fputil::FPBits<float>;
+ FPBits xbits(x);
+
+ uint32_t x_u = xbits.uintval();
+ uint32_t x_abs = x_u & 0x7fff'ffffU;
+
+ // When |x| >= log10(2^128), or x is nan
+ if (unlikely(x_abs >= 0x421a'209bU)) {
+ // When x < log10(2^-150) or nan
+ if (x_u > 0xc234'9e35U) {
+ // exp(-Inf) = 0
+ if (xbits.is_inf())
+ return 0.0f;
+ // exp(nan) = nan
+ if (xbits.is_nan())
+ return x;
+ if (fputil::get_round() == FE_UPWARD)
+ return static_cast<float>(FPBits(FPBits::MIN_SUBNORMAL));
+ errno = ERANGE;
+ return 0.0f;
+ }
+ // x >= log10(2^128) or nan
+ if (!xbits.get_sign() && (x_u >= 0x421a'209bU)) {
+ // x is finite
+ if (x_u < 0x7f80'0000U) {
+ int rounding = fputil::get_round();
+ if (rounding == FE_DOWNWARD || rounding == FE_TOWARDZERO)
+ return static_cast<float>(FPBits(FPBits::MAX_NORMAL));
+
+ errno = ERANGE;
+ }
+ // x is +inf or nan
+ return x + static_cast<float>(FPBits::inf());
+ }
+ }
+
+ // When |x| <= log10(2)*2^-6
+ if (unlikely(x_abs <= 0x3b9a'209bU)) {
+ if (unlikely(x_u == 0xb25e'5bd9U)) { // x = -0x1.bcb7b2p-27f
+ if (fputil::get_round() == FE_TONEAREST)
+ return 0x1.fffffep-1f;
+ }
+ // |x| < 2^-25
+ // 10^x ~ 1 + log(10) * x
+ if (unlikely(x_abs <= 0x3280'0000U)) {
+ return fputil::multiply_add(x, 0x1.26bb1cp+1f, 1.0f);
+ }
+
+ return Exp10Base::powb_lo(x);
+ }
+
+ // Exceptional value.
+ if (unlikely(x_u == 0x3d14'd956U)) { // x = 0x1.29b2acp-5f
+ if (fputil::get_round() == FE_UPWARD)
+ return 0x1.1657c4p+0f;
+ }
+
+ // Exact outputs when x = 1, 2, ..., 10.
+ // Quick check mask: 0x800f'ffffU = ~(bits of 1.0f | ... | bits of 10.0f)
+ if (unlikely((x_u & 0x800f'ffffU) == 0)) {
+ switch (x_u) {
+ case 0x3f800000U: // x = 1.0f
+ return 10.0f;
+ case 0x40000000U: // x = 2.0f
+ return 100.0f;
+ case 0x40400000U: // x = 3.0f
+ return 1'000.0f;
+ case 0x40800000U: // x = 4.0f
+ return 10'000.0f;
+ case 0x40a00000U: // x = 5.0f
+ return 100'000.0f;
+ case 0x40c00000U: // x = 6.0f
+ return 1'000'000.0f;
+ case 0x40e00000U: // x = 7.0f
+ return 10'000'000.0f;
+ case 0x41000000U: // x = 8.0f
+ return 100'000'000.0f;
+ case 0x41100000U: // x = 9.0f
+ return 1'000'000'000.0f;
+ case 0x41200000U: // x = 10.0f
+ return 10'000'000'000.0f;
+ }
+ }
+
+ // Range reduction: 10^x = 2^(mid + hi) * 10^lo
+ // rr = (2^(mid + hi), lo)
+ auto rr = exp_b_range_reduc<Exp10Base>(x);
+
+ // The low part is approximated by a degree-5 minimax polynomial.
+ // 10^lo ~ 1 + COEFFS[0] * lo + ... + COEFFS[4] * lo^5
+ using fputil::multiply_add;
+ double lo2 = rr.lo * rr.lo;
+ // c0 = 1 + COEFFS[0] * lo
+ double c0 = multiply_add(rr.lo, Exp10Base::COEFFS[0], 1.0);
+ // c1 = COEFFS[1] + COEFFS[2] * lo
+ double c1 = multiply_add(rr.lo, Exp10Base::COEFFS[2], Exp10Base::COEFFS[1]);
+ // c2 = COEFFS[3] + COEFFS[4] * lo
+ double c2 = multiply_add(rr.lo, Exp10Base::COEFFS[4], Exp10Base::COEFFS[3]);
+ // p = c1 + c2 * lo^2
+ // = COEFFS[1] + COEFFS[2] * lo + COEFFS[3] * lo^2 + COEFFS[4] * lo^3
+ double p = multiply_add(lo2, c2, c1);
+ // 10^lo ~ c0 + p * lo^2
+ // 10^x = 2^(mid + hi) * 10^lo
+ // ~ mh * (c0 + p * lo^2)
+ // = (mh * c0) + p * (mh * lo^2)
+ return multiply_add(p, lo2 * rr.mh, c0 * rr.mh);
+}
+
+} // namespace __llvm_libc
diff --git a/libc/src/math/generic/explogxf.h b/libc/src/math/generic/explogxf.h
index 482f722494e3b..b249bb8ab4c71 100644
--- a/libc/src/math/generic/explogxf.h
+++ b/libc/src/math/generic/explogxf.h
@@ -51,7 +51,7 @@ struct ExpBase {
// Approximating e^dx with degree-5 minimax polynomial generated by Sollya:
// > Q = fpminimax(expm1(x)/x, 4, [|1, D...|], [-log(2)/64, log(2)/64]);
// Then:
- // e^dx ~ P(dx) = 1 + dx + COEFFS[0] * dx^2 + ... + COEFFS[4] * dx^6.
+ // e^dx ~ P(dx) = 1 + dx + COEFFS[0] * dx^2 + ... + COEFFS[3] * dx^5.
static constexpr double COEFFS[4] = {
0x1.ffffffffe5bc8p-2, 0x1.555555555cd67p-3, 0x1.5555c2a9b48b4p-5,
0x1.11112a0e34bdbp-7};
@@ -70,6 +70,46 @@ struct ExpBase {
}
};
+struct Exp10Base : public ExpBase {
+ // log2(10) * 2^5
+ static constexpr double LOG2_B = 0x1.a934f0979a371p1 * (1 << MID_BITS);
+ // High and low parts of -log10(2) * 2^(-5).
+ // Notice that since |x * log2(10)| < 150:
+ // |k| = |round(x * log2(10) * 2^5)| < 2^8 * 2^5 = 2^13
+ // So when the FMA instructions are not available, in order for the product
+ // k * M_LOGB_2_HI
+ // to be exact, we only store the high part of log10(2) up to 38 bits
+ // (= 53 - 15) of precision.
+ // It is generated by Sollya with:
+ // > round(log10(2), 44, RN);
+ static constexpr double M_LOGB_2_HI = -0x1.34413509f8p-2 / (1 << MID_BITS);
+ // > round(log10(2) - 0x1.34413509f8p-2, D, RN);
+ static constexpr double M_LOGB_2_LO = 0x1.80433b83b532ap-44 / (1 << MID_BITS);
+
+ // Approximating 10^dx with degree-5 minimax polynomial generated by Sollya:
+ // > Q = fpminimax((10^x - 1)/x, 4, [|D...|], [-log10(2)/2^6, log10(2)/2^6]);
+ // Then:
+ // 10^dx ~ P(dx) = 1 + COEFFS[0] * dx + ... + COEFFS[4] * dx^5.
+ static constexpr double COEFFS[5] = {0x1.26bb1bbb55515p1, 0x1.53524c73bd3eap1,
+ 0x1.0470591dff149p1, 0x1.2bd7c0a9fbc4dp0,
+ 0x1.1429e74a98f43p-1};
+
+ static constexpr double powb_lo(double dx) {
+ using fputil::multiply_add;
+ double dx2 = dx * dx;
+ // c0 = 1 + COEFFS[0] * dx
+ double c0 = multiply_add(dx, Exp10Base::COEFFS[0], 1.0);
+ // c1 = COEFFS[1] + COEFFS[2] * dx
+ double c1 = multiply_add(dx, Exp10Base::COEFFS[2], Exp10Base::COEFFS[1]);
+ // c2 = COEFFS[3] + COEFFS[4] * dx
+ double c2 = multiply_add(dx, Exp10Base::COEFFS[4], Exp10Base::COEFFS[3]);
+ // r = c0 + dx^2 * (c1 + c2 * dx^2)
+ // = c0 + c1 * dx^2 + c2 * dx^4
+ // = 1 + COEFFS[0] * dx + ... + COEFFS[4] * dx^5.
+ return fputil::polyeval(dx2, c0, c1, c2);
+ }
+};
+
constexpr int LOG_P1_BITS = 6;
constexpr int LOG_P1_SIZE = 1 << LOG_P1_BITS;
diff --git a/libc/test/src/math/CMakeLists.txt b/libc/test/src/math/CMakeLists.txt
index 0cce635be245b..5162b2fe68c65 100644
--- a/libc/test/src/math/CMakeLists.txt
+++ b/libc/test/src/math/CMakeLists.txt
@@ -611,6 +611,21 @@ add_fp_unittest(
libc.src.__support.FPUtil.fp_bits
)
+add_fp_unittest(
+ exp10f_test
+ NEED_MPFR
+ SUITE
+ libc_math_unittests
+ SRCS
+ exp10f_test.cpp
+ DEPENDS
+ libc.include.errno
+ libc.src.errno.errno
+ libc.include.math
+ libc.src.math.exp10f
+ libc.src.__support.FPUtil.fp_bits
+)
+
add_fp_unittest(
copysign_test
SUITE
diff --git a/libc/test/src/math/exhaustive/CMakeLists.txt b/libc/test/src/math/exhaustive/CMakeLists.txt
index 819e78922cb1a..4d812bfa4e774 100644
--- a/libc/test/src/math/exhaustive/CMakeLists.txt
+++ b/libc/test/src/math/exhaustive/CMakeLists.txt
@@ -123,6 +123,23 @@ add_fp_unittest(
-lpthread
)
+add_fp_unittest(
+ exp10f_test
+ NO_RUN_POSTBUILD
+ NEED_MPFR
+ SUITE
+ libc_math_exhaustive_tests
+ SRCS
+ exp10f_test.cpp
+ DEPENDS
+ .exhaustive_test
+ libc.include.math
+ libc.src.math.exp10f
+ libc.src.__support.FPUtil.fp_bits
+ LINK_LIBRARIES
+ -lpthread
+)
+
add_fp_unittest(
expm1f_test
NO_RUN_POSTBUILD
diff --git a/libc/test/src/math/exhaustive/exp10f_test.cpp b/libc/test/src/math/exhaustive/exp10f_test.cpp
new file mode 100644
index 0000000000000..8e37f0483d2ba
--- /dev/null
+++ b/libc/test/src/math/exhaustive/exp10f_test.cpp
@@ -0,0 +1,75 @@
+//===-- Exhaustive test for exp10f ----------------------------------------===//
+//
+// Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions.
+// See https://llvm.org/LICENSE.txt for license information.
+// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
+//
+//===----------------------------------------------------------------------===//
+
+#include "exhaustive_test.h"
+#include "src/__support/FPUtil/FPBits.h"
+#include "src/math/exp10f.h"
+#include "utils/MPFRWrapper/MPFRUtils.h"
+#include "utils/UnitTest/FPMatcher.h"
+
+#include <thread>
+
+using FPBits = __llvm_libc::fputil::FPBits<float>;
+
+namespace mpfr = __llvm_libc::testing::mpfr;
+
+struct LlvmLibcExp10fExhaustiveTest : public LlvmLibcExhaustiveTest<uint32_t> {
+ bool check(uint32_t start, uint32_t stop,
+ mpfr::RoundingMode rounding) override {
+ mpfr::ForceRoundingMode r(rounding);
+ uint32_t bits = start;
+ bool result = true;
+ do {
+ FPBits xbits(bits);
+ float x = float(xbits);
+ result &= EXPECT_MPFR_MATCH(mpfr::Operation::Exp10, x,
+ __llvm_libc::exp10f(x), 0.5, rounding);
+ } while (bits++ < stop);
+ return result;
+ }
+};
+
+// Range: [0, 89];
+static constexpr uint32_t POS_START = 0x0000'0000U;
+static constexpr uint32_t POS_STOP = 0x42b2'0000U;
+
+TEST_F(LlvmLibcExp10fExhaustiveTest, PostiveRangeRoundNearestTieToEven) {
+ test_full_range(POS_START, POS_STOP, mpfr::RoundingMode::Nearest);
+}
+
+TEST_F(LlvmLibcExp10fExhaustiveTest, PostiveRangeRoundUp) {
+ test_full_range(POS_START, POS_STOP, mpfr::RoundingMode::Upward);
+}
+
+TEST_F(LlvmLibcExp10fExhaustiveTest, PostiveRangeRoundDown) {
+ test_full_range(POS_START, POS_STOP, mpfr::RoundingMode::Downward);
+}
+
+TEST_F(LlvmLibcExp10fExhaustiveTest, PostiveRangeRoundTowardZero) {
+ test_full_range(POS_START, POS_STOP, mpfr::RoundingMode::TowardZero);
+}
+
+// Range: [-104, 0];
+static constexpr uint32_t NEG_START = 0x8000'0000U;
+static constexpr uint32_t NEG_STOP = 0xc2d0'0000U;
+
+TEST_F(LlvmLibcExp10fExhaustiveTest, NegativeRangeRoundNearestTieToEven) {
+ test_full_range(NEG_START, NEG_STOP, mpfr::RoundingMode::Nearest);
+}
+
+TEST_F(LlvmLibcExp10fExhaustiveTest, NegativeRangeRoundUp) {
+ test_full_range(NEG_START, NEG_STOP, mpfr::RoundingMode::Upward);
+}
+
+TEST_F(LlvmLibcExp10fExhaustiveTest, NegativeRangeRoundDown) {
+ test_full_range(NEG_START, NEG_STOP, mpfr::RoundingMode::Downward);
+}
+
+TEST_F(LlvmLibcExp10fExhaustiveTest, NegativeRangeRoundTowardZero) {
+ test_full_range(NEG_START, NEG_STOP, mpfr::RoundingMode::TowardZero);
+}
diff --git a/libc/test/src/math/exp10f_test.cpp b/libc/test/src/math/exp10f_test.cpp
new file mode 100644
index 0000000000000..eb17e1d3c35b7
--- /dev/null
+++ b/libc/test/src/math/exp10f_test.cpp
@@ -0,0 +1,107 @@
+//===-- Unittests for exp10f ----------------------------------------------===//
+//
+// Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions.
+// See https://llvm.org/LICENSE.txt for license information.
+// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
+//
+//===----------------------------------------------------------------------===//
+
+#include "src/__support/FPUtil/FPBits.h"
+#include "src/math/exp10f.h"
+#include "utils/MPFRWrapper/MPFRUtils.h"
+#include "utils/UnitTest/FPMatcher.h"
+#include "utils/UnitTest/Test.h"
+#include <math.h>
+
+#include <errno.h>
+#include <stdint.h>
+
+namespace mpfr = __llvm_libc::testing::mpfr;
+
+DECLARE_SPECIAL_CONSTANTS(float)
+
+TEST(LlvmLibcExp10fTest, SpecialNumbers) {
+ errno = 0;
+
+ EXPECT_FP_EQ(aNaN, __llvm_libc::exp10f(aNaN));
+ EXPECT_MATH_ERRNO(0);
+
+ EXPECT_FP_EQ(inf, __llvm_libc::exp10f(inf));
+ EXPECT_MATH_ERRNO(0);
+
+ EXPECT_FP_EQ(0.0f, __llvm_libc::exp10f(neg_inf));
+ EXPECT_MATH_ERRNO(0);
+
+ EXPECT_FP_EQ(1.0f, __llvm_libc::exp10f(0.0f));
+ EXPECT_MATH_ERRNO(0);
+
+ EXPECT_FP_EQ(1.0f, __llvm_libc::exp10f(-0.0f));
+ EXPECT_MATH_ERRNO(0);
+}
+
+TEST(LlvmLibcExp10fTest, Overflow) {
+ errno = 0;
+ EXPECT_FP_EQ(inf, __llvm_libc::exp10f(float(FPBits(0x7f7fffffU))));
+ EXPECT_MATH_ERRNO(ERANGE);
+
+ EXPECT_FP_EQ(inf, __llvm_libc::exp10f(float(FPBits(0x43000000U))));
+ EXPECT_MATH_ERRNO(ERANGE);
+
+ EXPECT_FP_EQ(inf, __llvm_libc::exp10f(float(FPBits(0x43000001U))));
+ EXPECT_MATH_ERRNO(ERANGE);
+}
+
+TEST(LlvmLibcExp10fTest, TrickyInputs) {
+ constexpr int N = 20;
+ constexpr uint32_t INPUTS[N] = {
+ 0x325e5bd8, // x = 0x1.bcb7bp-27f
+ 0x325e5bd9, // x = 0x1.bcb7b2p-27f
+ 0x325e5bda, // x = 0x1.bcb7b4p-27f
+ 0x3d14d956, // x = 0x1.29b2acp-5f
+ 0x4116498a, // x = 0x1.2c9314p3f
+ 0x4126f431, // x = 0x1.4de862p3f
+ 0x4187d13c, // x = 0x1.0fa278p4f
+ 0x4203e9da, // x = 0x1.07d3b4p5f
+ 0x420b5f5d, // x = 0x1.16bebap5f
+ 0x42349e35, // x = 0x1.693c6ap5f
+ 0x3f800000, // x = 1.0f
+ 0x40000000, // x = 2.0f
+ 0x40400000, // x = 3.0f
+ 0x40800000, // x = 4.0f
+ 0x40a00000, // x = 5.0f
+ 0x40c00000, // x = 6.0f
+ 0x40e00000, // x = 7.0f
+ 0x41000000, // x = 8.0f
+ 0x41100000, // x = 9.0f
+ 0x41200000, // x = 10.0f
+ };
+ for (int i = 0; i < N; ++i) {
+ errno = 0;
+ float x = float(FPBits(INPUTS[i]));
+ EXPECT_MPFR_MATCH_ALL_ROUNDING(mpfr::Operation::Exp10, x,
+ __llvm_libc::exp10f(x), 0.5);
+ EXPECT_MPFR_MATCH_ALL_ROUNDING(mpfr::Operation::Exp10, -x,
+ __llvm_libc::exp10f(-x), 0.5);
+ }
+}
+
+TEST(LlvmLibcExp10fTest, InFloatRange) {
+ constexpr uint32_t COUNT = 1000000;
+ constexpr uint32_t STEP = UINT32_MAX / COUNT;
+ for (uint32_t i = 0, v = 0; i <= COUNT; ++i, v += STEP) {
+ float x = float(FPBits(v));
+ if (isnan(x) || isinf(x))
+ continue;
+ errno = 0;
+ float result = __llvm_libc::exp10f(x);
+
+ // If the computation resulted in an error or did not produce valid result
+ // in the single-precision floating point range, then ignore comparing with
+ // MPFR result as MPFR can still produce valid results because of its
+ // wider precision.
+ if (isnan(result) || isinf(result) || errno != 0)
+ continue;
+ ASSERT_MPFR_MATCH_ALL_ROUNDING(mpfr::Operation::Exp10, x,
+ __llvm_libc::exp10f(x), 0.5);
+ }
+}
diff --git a/libc/utils/MPFRWrapper/MPFRUtils.cpp b/libc/utils/MPFRWrapper/MPFRUtils.cpp
index 0f71357e9302a..0b63fdd335f0c 100644
--- a/libc/utils/MPFRWrapper/MPFRUtils.cpp
+++ b/libc/utils/MPFRWrapper/MPFRUtils.cpp
@@ -238,6 +238,12 @@ class MPFRNumber {
return result;
}
+ MPFRNumber exp10() const {
+ MPFRNumber result(*this);
+ mpfr_exp10(result.value, value, mpfr_rounding);
+ return result;
+ }
+
MPFRNumber expm1() const {
MPFRNumber result(*this);
mpfr_expm1(result.value, value, mpfr_rounding);
@@ -550,6 +556,8 @@ unary_operation(Operation op, InputType input, unsigned int precision,
return mpfrInput.exp();
case Operation::Exp2:
return mpfrInput.exp2();
+ case Operation::Exp10:
+ return mpfrInput.exp10();
case Operation::Expm1:
return mpfrInput.expm1();
case Operation::Floor:
diff --git a/libc/utils/MPFRWrapper/MPFRUtils.h b/libc/utils/MPFRWrapper/MPFRUtils.h
index e8d5ecabf5ac3..7902da5e97737 100644
--- a/libc/utils/MPFRWrapper/MPFRUtils.h
+++ b/libc/utils/MPFRWrapper/MPFRUtils.h
@@ -34,6 +34,7 @@ enum class Operation : int {
Cosh,
Exp,
Exp2,
+ Exp10,
Expm1,
Floor,
Log,
More information about the libc-commits
mailing list