aboutsummaryrefslogtreecommitdiff
path: root/libc
diff options
context:
space:
mode:
Diffstat (limited to 'libc')
-rw-r--r--libc/shared/math.h1
-rw-r--r--libc/shared/math/exp10m1f16.h29
-rw-r--r--libc/src/__support/math/CMakeLists.txt17
-rw-r--r--libc/src/__support/math/exp10m1f16.h185
-rw-r--r--libc/src/math/generic/CMakeLists.txt13
-rw-r--r--libc/src/math/generic/exp10m1f16.cpp158
-rw-r--r--libc/test/shared/CMakeLists.txt1
-rw-r--r--libc/test/shared/shared_math_test.cpp1
8 files changed, 237 insertions, 168 deletions
diff --git a/libc/shared/math.h b/libc/shared/math.h
index cccd6a3..4b2a0d8 100644
--- a/libc/shared/math.h
+++ b/libc/shared/math.h
@@ -46,6 +46,7 @@
#include "math/exp10f.h"
#include "math/exp10f16.h"
#include "math/exp10m1f.h"
+#include "math/exp10m1f16.h"
#include "math/expf.h"
#include "math/expf16.h"
#include "math/frexpf.h"
diff --git a/libc/shared/math/exp10m1f16.h b/libc/shared/math/exp10m1f16.h
new file mode 100644
index 0000000..5f18f29
--- /dev/null
+++ b/libc/shared/math/exp10m1f16.h
@@ -0,0 +1,29 @@
+//===-- Shared exp10m1f16 function ------------------------------*- 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_SHARED_MATH_EXP10M1F16_H
+#define LLVM_LIBC_SHARED_MATH_EXP10M1F16_H
+
+#include "include/llvm-libc-macros/float16-macros.h"
+#include "shared/libc_common.h"
+
+#ifdef LIBC_TYPES_HAS_FLOAT16
+
+#include "src/__support/math/exp10m1f16.h"
+
+namespace LIBC_NAMESPACE_DECL {
+namespace shared {
+
+using math::exp10m1f16;
+
+} // namespace shared
+} // namespace LIBC_NAMESPACE_DECL
+
+#endif // LIBC_TYPES_HAS_FLOAT16
+
+#endif // LLVM_LIBC_SHARED_MATH_EXP10M1F16_H
diff --git a/libc/src/__support/math/CMakeLists.txt b/libc/src/__support/math/CMakeLists.txt
index 84c1b15..98f9bb42 100644
--- a/libc/src/__support/math/CMakeLists.txt
+++ b/libc/src/__support/math/CMakeLists.txt
@@ -499,6 +499,23 @@ add_header_library(
)
add_header_library(
+ exp10m1f16
+ HDRS
+ exp10m1f16.h
+ DEPENDS
+ .exp10f16_utils
+ libc.src.__support.FPUtil.cast
+ libc.src.__support.FPUtil.except_value_utils
+ libc.src.__support.FPUtil.fenv_impl
+ libc.src.__support.FPUtil.fp_bits
+ libc.src.__support.FPUtil.multiply_add
+ libc.src.__support.FPUtil.polyeval
+ libc.src.__support.FPUtil.rounding_mode
+ libc.src.__support.macros.optimization
+ libc.src.__support.macros.properties.cpu_features
+)
+
+add_header_library(
erff
HDRS
erff.h
diff --git a/libc/src/__support/math/exp10m1f16.h b/libc/src/__support/math/exp10m1f16.h
new file mode 100644
index 0000000..6367a85
--- /dev/null
+++ b/libc/src/__support/math/exp10m1f16.h
@@ -0,0 +1,185 @@
+//===-- Implementation header for exp10m1f16 --------------------*- 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___SUPPORT_MATH_EXP10M1F16_H
+#define LLVM_LIBC_SRC___SUPPORT_MATH_EXP10M1F16_H
+
+#include "include/llvm-libc-macros/float16-macros.h"
+
+#ifdef LIBC_TYPES_HAS_FLOAT16
+
+#include "exp10f16_utils.h"
+#include "src/__support/FPUtil/FEnvImpl.h"
+#include "src/__support/FPUtil/FPBits.h"
+#include "src/__support/FPUtil/PolyEval.h"
+#include "src/__support/FPUtil/cast.h"
+#include "src/__support/FPUtil/except_value_utils.h"
+#include "src/__support/FPUtil/multiply_add.h"
+#include "src/__support/FPUtil/rounding_mode.h"
+#include "src/__support/common.h"
+#include "src/__support/macros/config.h"
+#include "src/__support/macros/optimization.h"
+#include "src/__support/macros/properties/cpu_features.h"
+
+namespace LIBC_NAMESPACE_DECL {
+
+namespace math {
+
+LIBC_INLINE static constexpr float16 exp10m1f16(float16 x) {
+
+#ifndef LIBC_MATH_HAS_SKIP_ACCURATE_PASS
+ constexpr fputil::ExceptValues<float16, 3> EXP10M1F16_EXCEPTS_LO = {{
+ // (input, RZ output, RU offset, RD offset, RN offset)
+ // x = 0x1.5c4p-4, exp10m1f16(x) = 0x1.bacp-3 (RZ)
+ {0x2d71U, 0x32ebU, 1U, 0U, 0U},
+ // x = -0x1.5ep-13, exp10m1f16(x) = -0x1.92cp-12 (RZ)
+ {0x8978U, 0x8e4bU, 0U, 1U, 0U},
+ // x = -0x1.e2p-10, exp10m1f16(x) = -0x1.14cp-8 (RZ)
+ {0x9788U, 0x9c53U, 0U, 1U, 0U},
+ }};
+
+#ifdef LIBC_TARGET_CPU_HAS_FMA_FLOAT
+ constexpr size_t N_EXP10M1F16_EXCEPTS_HI = 3;
+#else
+ constexpr size_t N_EXP10M1F16_EXCEPTS_HI = 6;
+#endif
+
+ constexpr fputil::ExceptValues<float16, N_EXP10M1F16_EXCEPTS_HI>
+ EXP10M1F16_EXCEPTS_HI = {{
+ // (input, RZ output, RU offset, RD offset, RN offset)
+ // x = 0x1.8f4p-2, exp10m1f16(x) = 0x1.744p+0 (RZ)
+ {0x363dU, 0x3dd1U, 1U, 0U, 0U},
+ // x = 0x1.95cp-2, exp10m1f16(x) = 0x1.7d8p+0 (RZ)
+ {0x3657U, 0x3df6U, 1U, 0U, 0U},
+ // x = 0x1.d04p-2, exp10m1f16(x) = 0x1.d7p+0 (RZ)
+ {0x3741U, 0x3f5cU, 1U, 0U, 1U},
+#ifndef LIBC_TARGET_CPU_HAS_FMA_FLOAT
+ // x = 0x1.0cp+1, exp10m1f16(x) = 0x1.ec4p+6 (RZ)
+ {0x4030U, 0x57b1U, 1U, 0U, 1U},
+ // x = 0x1.1b8p+1, exp10m1f16(x) = 0x1.45cp+7 (RZ)
+ {0x406eU, 0x5917U, 1U, 0U, 1U},
+ // x = 0x1.2f4p+2, exp10m1f16(x) = 0x1.ab8p+15 (RZ)
+ {0x44bdU, 0x7aaeU, 1U, 0U, 1U},
+#endif
+ }};
+#endif // !LIBC_MATH_HAS_SKIP_ACCURATE_PASS
+
+ using FPBits = fputil::FPBits<float16>;
+ FPBits x_bits(x);
+
+ uint16_t x_u = x_bits.uintval();
+ uint16_t x_abs = x_u & 0x7fffU;
+
+ // When |x| <= 2^(-3), or |x| >= 11 * log10(2), or x is NaN.
+ if (LIBC_UNLIKELY(x_abs <= 0x3000U || x_abs >= 0x429fU)) {
+ // exp10m1(NaN) = NaN
+ if (x_bits.is_nan()) {
+ if (x_bits.is_signaling_nan()) {
+ fputil::raise_except_if_required(FE_INVALID);
+ return FPBits::quiet_nan().get_val();
+ }
+
+ return x;
+ }
+
+ // When x >= 16 * log10(2).
+ if (x_u >= 0x44d1U && x_bits.is_pos()) {
+ // exp10m1(+inf) = +inf
+ if (x_bits.is_inf())
+ return FPBits::inf().get_val();
+
+ switch (fputil::quick_get_round()) {
+ case FE_TONEAREST:
+ case FE_UPWARD:
+ fputil::set_errno_if_required(ERANGE);
+ fputil::raise_except_if_required(FE_OVERFLOW | FE_INEXACT);
+ return FPBits::inf().get_val();
+ default:
+ return FPBits::max_normal().get_val();
+ }
+ }
+
+ // When x < -11 * log10(2).
+ if (x_u > 0xc29fU) {
+ // exp10m1(-inf) = -1
+ if (x_bits.is_inf())
+ return FPBits::one(Sign::NEG).get_val();
+
+ // When x >= -0x1.ce4p+1, round(10^x - 1, HP, RN) = -0x1.ffcp-1.
+ if (x_u <= 0xc339U) {
+ return fputil::round_result_slightly_down(
+ fputil::cast<float16>(-0x1.ffcp-1));
+ }
+
+ // When x < -0x1.ce4p+1, round(10^x - 1, HP, RN) = -1.
+ switch (fputil::quick_get_round()) {
+ case FE_TONEAREST:
+ case FE_DOWNWARD:
+ return FPBits::one(Sign::NEG).get_val();
+ default:
+ return fputil::cast<float16>(-0x1.ffcp-1);
+ }
+ }
+
+ // When |x| <= 2^(-3).
+ if (x_abs <= 0x3000U) {
+ if (LIBC_UNLIKELY(x_abs == 0))
+ return x;
+
+#ifndef LIBC_MATH_HAS_SKIP_ACCURATE_PASS
+ if (auto r = EXP10M1F16_EXCEPTS_LO.lookup(x_u);
+ LIBC_UNLIKELY(r.has_value()))
+ return r.value();
+#endif // !LIBC_MATH_HAS_SKIP_ACCURATE_PASS
+
+ float xf = x;
+ // Degree-5 minimax polynomial generated by Sollya with the following
+ // commands:
+ // > display = hexadecimal;
+ // > P = fpminimax((10^x - 1)/x, 4, [|SG...|], [-2^-3, 2^-3]);
+ // > x * P;
+ return fputil::cast<float16>(
+ xf * fputil::polyeval(xf, 0x1.26bb1cp+1f, 0x1.5351c8p+1f,
+ 0x1.04704p+1f, 0x1.2ce084p+0f, 0x1.14a6bep-1f));
+ }
+ }
+
+ // When x is 1, 2, or 3. These are hard-to-round cases with exact results.
+ // 10^4 - 1 = 9'999 is not exactly representable as a float16, but luckily the
+ // polynomial approximation gives the correct result for x = 4 in all
+ // rounding modes.
+ if (LIBC_UNLIKELY((x_u & ~(0x3c00U | 0x4000U | 0x4200U | 0x4400U)) == 0)) {
+ switch (x_u) {
+ case 0x3c00U: // x = 1.0f16
+ return fputil::cast<float16>(9.0);
+ case 0x4000U: // x = 2.0f16
+ return fputil::cast<float16>(99.0);
+ case 0x4200U: // x = 3.0f16
+ return fputil::cast<float16>(999.0);
+ }
+ }
+
+#ifndef LIBC_MATH_HAS_SKIP_ACCURATE_PASS
+ if (auto r = EXP10M1F16_EXCEPTS_HI.lookup(x_u); LIBC_UNLIKELY(r.has_value()))
+ return r.value();
+#endif // !LIBC_MATH_HAS_SKIP_ACCURATE_PASS
+
+ // exp10(x) = exp2((hi + mid) * log2(10)) * exp10(lo)
+ auto [exp2_hi_mid, exp10_lo] = exp10_range_reduction(x);
+ // exp10m1(x) = exp2((hi + mid) * log2(lo)) * exp10(lo) - 1
+ return fputil::cast<float16>(
+ fputil::multiply_add(exp2_hi_mid, exp10_lo, -1.0f));
+}
+
+} // namespace math
+
+} // namespace LIBC_NAMESPACE_DECL
+
+#endif // LIBC_TYPES_HAS_FLOAT16
+
+#endif // LLVM_LIBC_SRC___SUPPORT_MATH_EXP10M1F16_H
diff --git a/libc/src/math/generic/CMakeLists.txt b/libc/src/math/generic/CMakeLists.txt
index 8074a39..99c1b08 100644
--- a/libc/src/math/generic/CMakeLists.txt
+++ b/libc/src/math/generic/CMakeLists.txt
@@ -1603,18 +1603,7 @@ add_entrypoint_object(
HDRS
../exp10m1f16.h
DEPENDS
- libc.hdr.errno_macros
- libc.hdr.fenv_macros
- libc.src.__support.FPUtil.cast
- libc.src.__support.FPUtil.except_value_utils
- libc.src.__support.FPUtil.fenv_impl
- libc.src.__support.FPUtil.fp_bits
- libc.src.__support.FPUtil.multiply_add
- libc.src.__support.FPUtil.polyeval
- libc.src.__support.FPUtil.rounding_mode
- libc.src.__support.macros.optimization
- libc.src.__support.macros.properties.cpu_features
- libc.src.__support.math.exp10f16_utils
+ libc.src.__support.math.exp10m1f16
)
add_entrypoint_object(
diff --git a/libc/src/math/generic/exp10m1f16.cpp b/libc/src/math/generic/exp10m1f16.cpp
index 6c2fdbe..8a3c4ab 100644
--- a/libc/src/math/generic/exp10m1f16.cpp
+++ b/libc/src/math/generic/exp10m1f16.cpp
@@ -7,166 +7,12 @@
//===----------------------------------------------------------------------===//
#include "src/math/exp10m1f16.h"
-#include "hdr/errno_macros.h"
-#include "hdr/fenv_macros.h"
-#include "src/__support/FPUtil/FEnvImpl.h"
-#include "src/__support/FPUtil/FPBits.h"
-#include "src/__support/FPUtil/PolyEval.h"
-#include "src/__support/FPUtil/cast.h"
-#include "src/__support/FPUtil/except_value_utils.h"
-#include "src/__support/FPUtil/multiply_add.h"
-#include "src/__support/FPUtil/rounding_mode.h"
-#include "src/__support/common.h"
-#include "src/__support/macros/config.h"
-#include "src/__support/macros/optimization.h"
-#include "src/__support/macros/properties/cpu_features.h"
-#include "src/__support/math/exp10f16_utils.h"
+#include "src/__support/math/exp10m1f16.h"
namespace LIBC_NAMESPACE_DECL {
-#ifndef LIBC_MATH_HAS_SKIP_ACCURATE_PASS
-static constexpr fputil::ExceptValues<float16, 3> EXP10M1F16_EXCEPTS_LO = {{
- // (input, RZ output, RU offset, RD offset, RN offset)
- // x = 0x1.5c4p-4, exp10m1f16(x) = 0x1.bacp-3 (RZ)
- {0x2d71U, 0x32ebU, 1U, 0U, 0U},
- // x = -0x1.5ep-13, exp10m1f16(x) = -0x1.92cp-12 (RZ)
- {0x8978U, 0x8e4bU, 0U, 1U, 0U},
- // x = -0x1.e2p-10, exp10m1f16(x) = -0x1.14cp-8 (RZ)
- {0x9788U, 0x9c53U, 0U, 1U, 0U},
-}};
-
-#ifdef LIBC_TARGET_CPU_HAS_FMA_FLOAT
-static constexpr size_t N_EXP10M1F16_EXCEPTS_HI = 3;
-#else
-static constexpr size_t N_EXP10M1F16_EXCEPTS_HI = 6;
-#endif
-
-static constexpr fputil::ExceptValues<float16, N_EXP10M1F16_EXCEPTS_HI>
- EXP10M1F16_EXCEPTS_HI = {{
- // (input, RZ output, RU offset, RD offset, RN offset)
- // x = 0x1.8f4p-2, exp10m1f16(x) = 0x1.744p+0 (RZ)
- {0x363dU, 0x3dd1U, 1U, 0U, 0U},
- // x = 0x1.95cp-2, exp10m1f16(x) = 0x1.7d8p+0 (RZ)
- {0x3657U, 0x3df6U, 1U, 0U, 0U},
- // x = 0x1.d04p-2, exp10m1f16(x) = 0x1.d7p+0 (RZ)
- {0x3741U, 0x3f5cU, 1U, 0U, 1U},
-#ifndef LIBC_TARGET_CPU_HAS_FMA_FLOAT
- // x = 0x1.0cp+1, exp10m1f16(x) = 0x1.ec4p+6 (RZ)
- {0x4030U, 0x57b1U, 1U, 0U, 1U},
- // x = 0x1.1b8p+1, exp10m1f16(x) = 0x1.45cp+7 (RZ)
- {0x406eU, 0x5917U, 1U, 0U, 1U},
- // x = 0x1.2f4p+2, exp10m1f16(x) = 0x1.ab8p+15 (RZ)
- {0x44bdU, 0x7aaeU, 1U, 0U, 1U},
-#endif
- }};
-#endif // !LIBC_MATH_HAS_SKIP_ACCURATE_PASS
-
LLVM_LIBC_FUNCTION(float16, exp10m1f16, (float16 x)) {
- using FPBits = fputil::FPBits<float16>;
- FPBits x_bits(x);
-
- uint16_t x_u = x_bits.uintval();
- uint16_t x_abs = x_u & 0x7fffU;
-
- // When |x| <= 2^(-3), or |x| >= 11 * log10(2), or x is NaN.
- if (LIBC_UNLIKELY(x_abs <= 0x3000U || x_abs >= 0x429fU)) {
- // exp10m1(NaN) = NaN
- if (x_bits.is_nan()) {
- if (x_bits.is_signaling_nan()) {
- fputil::raise_except_if_required(FE_INVALID);
- return FPBits::quiet_nan().get_val();
- }
-
- return x;
- }
-
- // When x >= 16 * log10(2).
- if (x_u >= 0x44d1U && x_bits.is_pos()) {
- // exp10m1(+inf) = +inf
- if (x_bits.is_inf())
- return FPBits::inf().get_val();
-
- switch (fputil::quick_get_round()) {
- case FE_TONEAREST:
- case FE_UPWARD:
- fputil::set_errno_if_required(ERANGE);
- fputil::raise_except_if_required(FE_OVERFLOW | FE_INEXACT);
- return FPBits::inf().get_val();
- default:
- return FPBits::max_normal().get_val();
- }
- }
-
- // When x < -11 * log10(2).
- if (x_u > 0xc29fU) {
- // exp10m1(-inf) = -1
- if (x_bits.is_inf())
- return FPBits::one(Sign::NEG).get_val();
-
- // When x >= -0x1.ce4p+1, round(10^x - 1, HP, RN) = -0x1.ffcp-1.
- if (x_u <= 0xc339U) {
- return fputil::round_result_slightly_down(
- fputil::cast<float16>(-0x1.ffcp-1));
- }
-
- // When x < -0x1.ce4p+1, round(10^x - 1, HP, RN) = -1.
- switch (fputil::quick_get_round()) {
- case FE_TONEAREST:
- case FE_DOWNWARD:
- return FPBits::one(Sign::NEG).get_val();
- default:
- return fputil::cast<float16>(-0x1.ffcp-1);
- }
- }
-
- // When |x| <= 2^(-3).
- if (x_abs <= 0x3000U) {
- if (LIBC_UNLIKELY(x_abs == 0))
- return x;
-
-#ifndef LIBC_MATH_HAS_SKIP_ACCURATE_PASS
- if (auto r = EXP10M1F16_EXCEPTS_LO.lookup(x_u);
- LIBC_UNLIKELY(r.has_value()))
- return r.value();
-#endif // !LIBC_MATH_HAS_SKIP_ACCURATE_PASS
-
- float xf = x;
- // Degree-5 minimax polynomial generated by Sollya with the following
- // commands:
- // > display = hexadecimal;
- // > P = fpminimax((10^x - 1)/x, 4, [|SG...|], [-2^-3, 2^-3]);
- // > x * P;
- return fputil::cast<float16>(
- xf * fputil::polyeval(xf, 0x1.26bb1cp+1f, 0x1.5351c8p+1f,
- 0x1.04704p+1f, 0x1.2ce084p+0f, 0x1.14a6bep-1f));
- }
- }
-
- // When x is 1, 2, or 3. These are hard-to-round cases with exact results.
- // 10^4 - 1 = 9'999 is not exactly representable as a float16, but luckily the
- // polynomial approximation gives the correct result for x = 4 in all
- // rounding modes.
- if (LIBC_UNLIKELY((x_u & ~(0x3c00U | 0x4000U | 0x4200U | 0x4400U)) == 0)) {
- switch (x_u) {
- case 0x3c00U: // x = 1.0f16
- return fputil::cast<float16>(9.0);
- case 0x4000U: // x = 2.0f16
- return fputil::cast<float16>(99.0);
- case 0x4200U: // x = 3.0f16
- return fputil::cast<float16>(999.0);
- }
- }
-
-#ifndef LIBC_MATH_HAS_SKIP_ACCURATE_PASS
- if (auto r = EXP10M1F16_EXCEPTS_HI.lookup(x_u); LIBC_UNLIKELY(r.has_value()))
- return r.value();
-#endif // !LIBC_MATH_HAS_SKIP_ACCURATE_PASS
-
- // exp10(x) = exp2((hi + mid) * log2(10)) * exp10(lo)
- auto [exp2_hi_mid, exp10_lo] = exp10_range_reduction(x);
- // exp10m1(x) = exp2((hi + mid) * log2(lo)) * exp10(lo) - 1
- return fputil::cast<float16>(
- fputil::multiply_add(exp2_hi_mid, exp10_lo, -1.0f));
+ return math::exp10m1f16(x);
}
} // namespace LIBC_NAMESPACE_DECL
diff --git a/libc/test/shared/CMakeLists.txt b/libc/test/shared/CMakeLists.txt
index 13a0aae..ea4634c 100644
--- a/libc/test/shared/CMakeLists.txt
+++ b/libc/test/shared/CMakeLists.txt
@@ -37,6 +37,7 @@ add_fp_unittest(
libc.src.__support.math.cospif16
libc.src.__support.math.dsqrtl
libc.src.__support.math.exp10m1f
+ libc.src.__support.math.exp10m1f16
libc.src.__support.math.erff
libc.src.__support.math.exp
libc.src.__support.math.exp10
diff --git a/libc/test/shared/shared_math_test.cpp b/libc/test/shared/shared_math_test.cpp
index 25bf5ad..1722193 100644
--- a/libc/test/shared/shared_math_test.cpp
+++ b/libc/test/shared/shared_math_test.cpp
@@ -27,6 +27,7 @@ TEST(LlvmLibcSharedMathTest, AllFloat16) {
EXPECT_FP_EQ(0x1p+0f16, LIBC_NAMESPACE::shared::coshf16(0.0f16));
EXPECT_FP_EQ(0x1p+0f16, LIBC_NAMESPACE::shared::cospif16(0.0f16));
EXPECT_FP_EQ(0x1p+0f16, LIBC_NAMESPACE::shared::exp10f16(0.0f16));
+ EXPECT_FP_EQ(0x0p+0f16, LIBC_NAMESPACE::shared::exp10m1f16(0.0f16));
EXPECT_FP_EQ(0x1p+0f16, LIBC_NAMESPACE::shared::expf16(0.0f16));