aboutsummaryrefslogtreecommitdiff
path: root/libc/src/math
diff options
context:
space:
mode:
Diffstat (limited to 'libc/src/math')
-rw-r--r--libc/src/math/CMakeLists.txt1
-rw-r--r--libc/src/math/fabsbf16.h21
-rw-r--r--libc/src/math/generic/CMakeLists.txt69
-rw-r--r--libc/src/math/generic/asinf16.cpp121
-rw-r--r--libc/src/math/generic/asinhf.cpp106
-rw-r--r--libc/src/math/generic/asinhf16.cpp96
-rw-r--r--libc/src/math/generic/atan.cpp167
-rw-r--r--libc/src/math/generic/atan2.cpp3
-rw-r--r--libc/src/math/generic/atan2f128.cpp3
-rw-r--r--libc/src/math/generic/atan_utils.h241
-rw-r--r--libc/src/math/generic/fabsbf16.cpp19
11 files changed, 72 insertions, 775 deletions
diff --git a/libc/src/math/CMakeLists.txt b/libc/src/math/CMakeLists.txt
index 455ad34..0522e0e 100644
--- a/libc/src/math/CMakeLists.txt
+++ b/libc/src/math/CMakeLists.txt
@@ -189,6 +189,7 @@ add_math_entrypoint_object(fabsf)
add_math_entrypoint_object(fabsl)
add_math_entrypoint_object(fabsf16)
add_math_entrypoint_object(fabsf128)
+add_math_entrypoint_object(fabsbf16)
add_math_entrypoint_object(fadd)
add_math_entrypoint_object(faddl)
diff --git a/libc/src/math/fabsbf16.h b/libc/src/math/fabsbf16.h
new file mode 100644
index 0000000..4993668
--- /dev/null
+++ b/libc/src/math/fabsbf16.h
@@ -0,0 +1,21 @@
+//===-- Implementation header for fabsbf16 ----------------------*- 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_FABSBF16_H
+#define LLVM_LIBC_SRC_MATH_FABSBF16_H
+
+#include "src/__support/macros/config.h"
+#include "src/__support/macros/properties/types.h"
+
+namespace LIBC_NAMESPACE_DECL {
+
+bfloat16 fabsbf16(bfloat16 x);
+
+} // namespace LIBC_NAMESPACE_DECL
+
+#endif // LLVM_LIBC_SRC_MATH_FABSBF16_H
diff --git a/libc/src/math/generic/CMakeLists.txt b/libc/src/math/generic/CMakeLists.txt
index ecf0967..970908c 100644
--- a/libc/src/math/generic/CMakeLists.txt
+++ b/libc/src/math/generic/CMakeLists.txt
@@ -698,6 +698,19 @@ add_entrypoint_object(
)
add_entrypoint_object(
+ fabsbf16
+ SRCS
+ fabsbf16.cpp
+ HDRS
+ ../fabsbf16.h
+ DEPENDS
+ libc.src.__support.FPUtil.basic_operations
+ libc.src.__support.FPUtil.bfloat16
+ libc.src.__support.macros.config
+ libc.src.__support.macros.properties.types
+)
+
+add_entrypoint_object(
fadd
SRCS
fadd.cpp
@@ -3889,12 +3902,7 @@ add_entrypoint_object(
HDRS
../asinhf.h
DEPENDS
- .explogxf
- libc.src.__support.FPUtil.fp_bits
- libc.src.__support.FPUtil.multiply_add
- libc.src.__support.FPUtil.polyeval
- libc.src.__support.FPUtil.sqrt
- libc.src.__support.macros.optimization
+ libc.src.__support.math.asinhf
)
add_entrypoint_object(
@@ -3904,18 +3912,7 @@ add_entrypoint_object(
HDRS
../asinhf16.h
DEPENDS
- .explogxf
- 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.FPUtil.sqrt
- libc.src.__support.macros.optimization
- libc.src.__support.macros.properties.types
+ libc.src.__support.math.asinhf16
)
add_entrypoint_object(
@@ -3968,16 +3965,7 @@ add_entrypoint_object(
HDRS
../asinf16.h
DEPENDS
- libc.hdr.errno_macros
- libc.hdr.fenv_macros
- libc.src.__support.FPUtil.cast
- 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.sqrt
- libc.src.__support.macros.optimization
- libc.src.__support.macros.properties.types
+ libc.src.__support.math.asinf16
)
add_entrypoint_object(
@@ -4032,19 +4020,6 @@ add_entrypoint_object(
libc.src.errno.errno
)
-add_header_library(
- atan_utils
- HDRS
- atan_utils.h
- DEPENDS
- libc.src.__support.integer_literals
- libc.src.__support.FPUtil.double_double
- libc.src.__support.FPUtil.dyadic_float
- libc.src.__support.FPUtil.multiply_add
- libc.src.__support.FPUtil.polyeval
- libc.src.__support.macros.optimization
-)
-
add_entrypoint_object(
atanf
SRCS
@@ -4091,13 +4066,7 @@ add_entrypoint_object(
COMPILE_OPTIONS
-O3
DEPENDS
- .atan_utils
- libc.src.__support.FPUtil.double_double
- 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.macros.optimization
+ libc.src.__support.math.atan
)
add_entrypoint_object(
@@ -4127,7 +4096,7 @@ add_entrypoint_object(
HDRS
../atan2.h
DEPENDS
- .atan_utils
+ libc.src.__support.math.atan_utils
libc.src.__support.FPUtil.double_double
libc.src.__support.FPUtil.fenv_impl
libc.src.__support.FPUtil.fp_bits
@@ -4153,7 +4122,7 @@ add_entrypoint_object(
HDRS
../atan2f128.h
DEPENDS
- .atan_utils
+ libc.src.__support.math.atan_utils
libc.src.__support.integer_literals
libc.src.__support.uint128
libc.src.__support.FPUtil.dyadic_float
diff --git a/libc/src/math/generic/asinf16.cpp b/libc/src/math/generic/asinf16.cpp
index 518c384..af8dbfe 100644
--- a/libc/src/math/generic/asinf16.cpp
+++ b/libc/src/math/generic/asinf16.cpp
@@ -7,127 +7,10 @@
//===----------------------------------------------------------------------===//
#include "src/math/asinf16.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/multiply_add.h"
-#include "src/__support/FPUtil/sqrt.h"
-#include "src/__support/macros/optimization.h"
+#include "src/__support/math/asinf16.h"
namespace LIBC_NAMESPACE_DECL {
-// Generated by Sollya using the following command:
-// > round(pi/2, D, RN);
-static constexpr float PI_2 = 0x1.921fb54442d18p0f;
-
-LLVM_LIBC_FUNCTION(float16, asinf16, (float16 x)) {
- using FPBits = fputil::FPBits<float16>;
- FPBits xbits(x);
-
- uint16_t x_u = xbits.uintval();
- uint16_t x_abs = x_u & 0x7fff;
- float xf = x;
-
- // |x| > 0x1p0, |x| > 1, or x is NaN.
- if (LIBC_UNLIKELY(x_abs > 0x3c00)) {
- // asinf16(NaN) = NaN
- if (xbits.is_nan()) {
- if (xbits.is_signaling_nan()) {
- fputil::raise_except_if_required(FE_INVALID);
- return FPBits::quiet_nan().get_val();
- }
-
- return x;
- }
-
- // 1 < |x| <= +/-inf
- fputil::raise_except_if_required(FE_INVALID);
- fputil::set_errno_if_required(EDOM);
-
- return FPBits::quiet_nan().get_val();
- }
-
- float xsq = xf * xf;
-
- // |x| <= 0x1p-1, |x| <= 0.5
- if (x_abs <= 0x3800) {
- // asinf16(+/-0) = +/-0
- if (LIBC_UNLIKELY(x_abs == 0))
- return x;
-
- // Exhaustive tests show that,
- // for |x| <= 0x1.878p-9, when:
- // x > 0, and rounding upward, or
- // x < 0, and rounding downward, then,
- // asin(x) = x * 2^-11 + x
- // else, in other rounding modes,
- // asin(x) = x
- if (LIBC_UNLIKELY(x_abs <= 0x1a1e)) {
- int rounding = fputil::quick_get_round();
-
- if ((xbits.is_pos() && rounding == FE_UPWARD) ||
- (xbits.is_neg() && rounding == FE_DOWNWARD))
- return fputil::cast<float16>(fputil::multiply_add(xf, 0x1.0p-11f, xf));
- return x;
- }
-
- // Degree-6 minimax odd polynomial of asin(x) generated by Sollya with:
- // > P = fpminimax(asin(x)/x, [|0, 2, 4, 6, 8|], [|SG...|], [0, 0.5]);
- float result =
- fputil::polyeval(xsq, 0x1.000002p0f, 0x1.554c2ap-3f, 0x1.3541ccp-4f,
- 0x1.43b2d6p-5f, 0x1.a0d73ep-5f);
- return fputil::cast<float16>(xf * result);
- }
-
- // When |x| > 0.5, assume that 0.5 < |x| <= 1,
- //
- // Step-by-step range-reduction proof:
- // 1: Let y = asin(x), such that, x = sin(y)
- // 2: From complimentary angle identity:
- // x = sin(y) = cos(pi/2 - y)
- // 3: Let z = pi/2 - y, such that x = cos(z)
- // 4: From double angle formula; cos(2A) = 1 - sin^2(A):
- // z = 2A, z/2 = A
- // cos(z) = 1 - 2 * sin^2(z/2)
- // 5: Make sin(z/2) subject of the formula:
- // sin(z/2) = sqrt((1 - cos(z))/2)
- // 6: Recall [3]; x = cos(z). Therefore:
- // sin(z/2) = sqrt((1 - x)/2)
- // 7: Let u = (1 - x)/2
- // 8: Therefore:
- // asin(sqrt(u)) = z/2
- // 2 * asin(sqrt(u)) = z
- // 9: Recall [3], z = pi/2 - y. Therefore:
- // y = pi/2 - z
- // y = pi/2 - 2 * asin(sqrt(u))
- // 10: Recall [1], y = asin(x). Therefore:
- // asin(x) = pi/2 - 2 * asin(sqrt(u))
- //
- // WHY?
- // 11: Recall [7], u = (1 - x)/2
- // 12: Since 0.5 < x <= 1, therefore:
- // 0 <= u <= 0.25 and 0 <= sqrt(u) <= 0.5
- //
- // Hence, we can reuse the same [0, 0.5] domain polynomial approximation for
- // Step [10] as `sqrt(u)` is in range.
-
- // 0x1p-1 < |x| <= 0x1p0, 0.5 < |x| <= 1.0
- float xf_abs = (xf < 0 ? -xf : xf);
- float sign = (xbits.uintval() >> 15 == 1 ? -1.0 : 1.0);
- float u = fputil::multiply_add(-0.5f, xf_abs, 0.5f);
- float u_sqrt = fputil::sqrt<float>(u);
-
- // Degree-6 minimax odd polynomial of asin(x) generated by Sollya with:
- // > P = fpminimax(asin(x)/x, [|0, 2, 4, 6, 8|], [|SG...|], [0, 0.5]);
- float asin_sqrt_u =
- u_sqrt * fputil::polyeval(u, 0x1.000002p0f, 0x1.554c2ap-3f,
- 0x1.3541ccp-4f, 0x1.43b2d6p-5f, 0x1.a0d73ep-5f);
-
- return fputil::cast<float16>(sign *
- fputil::multiply_add(-2.0f, asin_sqrt_u, PI_2));
-}
+LLVM_LIBC_FUNCTION(float16, asinf16, (float16 x)) { return math::asinf16(x); }
} // namespace LIBC_NAMESPACE_DECL
diff --git a/libc/src/math/generic/asinhf.cpp b/libc/src/math/generic/asinhf.cpp
index 3aed3bc..45023c8 100644
--- a/libc/src/math/generic/asinhf.cpp
+++ b/libc/src/math/generic/asinhf.cpp
@@ -7,112 +7,10 @@
//===----------------------------------------------------------------------===//
#include "src/math/asinhf.h"
-#include "src/__support/FPUtil/FPBits.h"
-#include "src/__support/FPUtil/PolyEval.h"
-#include "src/__support/FPUtil/multiply_add.h"
-#include "src/__support/FPUtil/sqrt.h"
-#include "src/__support/macros/config.h"
-#include "src/__support/macros/optimization.h" // LIBC_UNLIKELY
-#include "src/math/generic/common_constants.h"
-#include "src/math/generic/explogxf.h"
+#include "src/__support/math/asinhf.h"
namespace LIBC_NAMESPACE_DECL {
-LLVM_LIBC_FUNCTION(float, asinhf, (float x)) {
- using namespace acoshf_internal;
- using FPBits_t = typename fputil::FPBits<float>;
- FPBits_t xbits(x);
- uint32_t x_u = xbits.uintval();
- uint32_t x_abs = xbits.abs().uintval();
-
- // |x| <= 2^-3
- if (LIBC_UNLIKELY(x_abs <= 0x3e80'0000U)) {
- // |x| <= 2^-26
- if (LIBC_UNLIKELY(x_abs <= 0x3280'0000U)) {
- return static_cast<float>(LIBC_UNLIKELY(x_abs == 0)
- ? x
- : (x - 0x1.5555555555555p-3 * x * x * x));
- }
-
- double x_d = x;
- double x_sq = x_d * x_d;
- // Generated by Sollya with:
- // > P = fpminimax(asinh(x)/x, [|0, 2, 4, 6, 8, 10, 12, 14, 16|], [|D...|],
- // [0, 2^-2]);
- double p = fputil::polyeval(
- x_sq, 0.0, -0x1.555555555551ep-3, 0x1.3333333325495p-4,
- -0x1.6db6db5a7622bp-5, 0x1.f1c70f82928c6p-6, -0x1.6e893934266b7p-6,
- 0x1.1c0b41d3fbe78p-6, -0x1.c0f47810b3c4fp-7, 0x1.2c8602690143dp-7);
- return static_cast<float>(fputil::multiply_add(x_d, p, x_d));
- }
-
- const double SIGN[2] = {1.0, -1.0};
- double x_sign = SIGN[x_u >> 31];
- double x_d = x;
-
-#ifndef LIBC_MATH_HAS_SKIP_ACCURATE_PASS
- // Helper functions to set results for exceptional cases.
- auto round_result_slightly_down = [x_sign](float r) -> float {
- return fputil::multiply_add(static_cast<float>(x_sign), r,
- static_cast<float>(x_sign) * (-0x1.0p-24f));
- };
- auto round_result_slightly_up = [x_sign](float r) -> float {
- return fputil::multiply_add(static_cast<float>(x_sign), r,
- static_cast<float>(x_sign) * 0x1.0p-24f);
- };
-
- if (LIBC_UNLIKELY(x_abs >= 0x4bdd'65a5U)) {
- if (LIBC_UNLIKELY(xbits.is_inf_or_nan())) {
- if (xbits.is_signaling_nan()) {
- fputil::raise_except_if_required(FE_INVALID);
- return FPBits_t::quiet_nan().get_val();
- }
-
- return x;
- }
-
- // Exceptional cases when x > 2^24.
- switch (x_abs) {
- case 0x4bdd65a5: // |x| = 0x1.bacb4ap24f
- return round_result_slightly_down(0x1.1e0696p4f);
- case 0x4c803f2c: // |x| = 0x1.007e58p26f
- return round_result_slightly_down(0x1.2b786cp4f);
- case 0x4f8ffb03: // |x| = 0x1.1ff606p32f
- return round_result_slightly_up(0x1.6fdd34p4f);
- case 0x5c569e88: // |x| = 0x1.ad3d1p57f
- return round_result_slightly_up(0x1.45c146p5f);
- case 0x5e68984e: // |x| = 0x1.d1309cp61f
- return round_result_slightly_up(0x1.5c9442p5f);
- case 0x655890d3: // |x| = 0x1.b121a6p75f
- return round_result_slightly_down(0x1.a9a3f2p5f);
- case 0x65de7ca6: // |x| = 0x1.bcf94cp76f
- return round_result_slightly_up(0x1.af66cp5f);
- case 0x6eb1a8ec: // |x| = 0x1.6351d8p94f
- return round_result_slightly_down(0x1.08b512p6f);
- case 0x7997f30a: // |x| = 0x1.2fe614p116f
- return round_result_slightly_up(0x1.451436p6f);
- }
- } else {
- // Exceptional cases when x < 2^24.
- if (LIBC_UNLIKELY(x_abs == 0x45abaf26)) {
- // |x| = 0x1.575e4cp12f
- return round_result_slightly_down(0x1.29becap3f);
- }
- if (LIBC_UNLIKELY(x_abs == 0x49d29048)) {
- // |x| = 0x1.a5209p20f
- return round_result_slightly_down(0x1.e1b92p3f);
- }
- }
-#else
- if (LIBC_UNLIKELY(xbits.is_inf_or_nan()))
- return x;
-#endif // !LIBC_MATH_HAS_SKIP_ACCURATE_PASS
-
- // asinh(x) = log(x + sqrt(x^2 + 1))
- return static_cast<float>(
- x_sign * log_eval(fputil::multiply_add(
- x_d, x_sign,
- fputil::sqrt<double>(fputil::multiply_add(x_d, x_d, 1.0)))));
-}
+LLVM_LIBC_FUNCTION(float, asinhf, (float x)) { return math::asinhf(x); }
} // namespace LIBC_NAMESPACE_DECL
diff --git a/libc/src/math/generic/asinhf16.cpp b/libc/src/math/generic/asinhf16.cpp
index 0a0b471..d517e63 100644
--- a/libc/src/math/generic/asinhf16.cpp
+++ b/libc/src/math/generic/asinhf16.cpp
@@ -7,102 +7,10 @@
//===----------------------------------------------------------------------===//
#include "src/math/asinhf16.h"
-#include "explogxf.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/FPUtil/sqrt.h"
-#include "src/__support/common.h"
-#include "src/__support/macros/config.h"
-#include "src/__support/macros/optimization.h"
+#include "src/__support/math/asinhf16.h"
namespace LIBC_NAMESPACE_DECL {
-#ifndef LIBC_MATH_HAS_SKIP_ACCURATE_PASS
-static constexpr size_t N_EXCEPTS = 8;
-
-static constexpr fputil::ExceptValues<float16, N_EXCEPTS> ASINHF16_EXCEPTS{{
- // (input, RZ output, RU offset, RD offset, RN offset)
-
- // x = 0x1.da4p-2, asinhf16(x) = 0x1.ca8p-2 (RZ)
- {0x3769, 0x372a, 1, 0, 1},
- // x = 0x1.d6cp-1, asinhf16(x) = 0x1.a58p-1 (RZ)
- {0x3b5b, 0x3a96, 1, 0, 0},
- // x = 0x1.c7cp+3, asinhf16(x) = 0x1.accp+1 (RZ)
- {0x4b1f, 0x42b3, 1, 0, 0},
- // x = 0x1.26cp+4, asinhf16(x) = 0x1.cd8p+1 (RZ)
- {0x4c9b, 0x4336, 1, 0, 1},
- // x = -0x1.da4p-2, asinhf16(x) = -0x1.ca8p-2 (RZ)
- {0xb769, 0xb72a, 0, 1, 1},
- // x = -0x1.d6cp-1, asinhf16(x) = -0x1.a58p-1 (RZ)
- {0xbb5b, 0xba96, 0, 1, 0},
- // x = -0x1.c7cp+3, asinhf16(x) = -0x1.accp+1 (RZ)
- {0xcb1f, 0xc2b3, 0, 1, 0},
- // x = -0x1.26cp+4, asinhf16(x) = -0x1.cd8p+1 (RZ)
- {0xcc9b, 0xc336, 0, 1, 1},
-}};
-#endif // !LIBC_MATH_HAS_SKIP_ACCURATE_PASS
-
-LLVM_LIBC_FUNCTION(float16, asinhf16, (float16 x)) {
- using namespace acoshf_internal;
- using FPBits = fputil::FPBits<float16>;
- FPBits xbits(x);
-
- uint16_t x_u = xbits.uintval();
- uint16_t x_abs = x_u & 0x7fff;
-
- if (LIBC_UNLIKELY(xbits.is_inf_or_nan())) {
- if (xbits.is_signaling_nan()) {
- fputil::raise_except_if_required(FE_INVALID);
- return FPBits::quiet_nan().get_val();
- }
-
- return x;
- }
-
-#ifndef LIBC_MATH_HAS_SKIP_ACCURATE_PASS
- // Handle exceptional values
- if (auto r = ASINHF16_EXCEPTS.lookup(x_u); LIBC_UNLIKELY(r.has_value()))
- return r.value();
-#endif // !LIBC_MATH_HAS_SKIP_ACCURATE_PASS
-
- float xf = x;
- const float SIGN[2] = {1.0f, -1.0f};
- float x_sign = SIGN[x_u >> 15];
-
- // |x| <= 0.25
- if (LIBC_UNLIKELY(x_abs <= 0x3400)) {
- // when |x| < 0x1.718p-5, asinhf16(x) = x. Adjust by 1 ULP for certain
- // rounding types.
- if (LIBC_UNLIKELY(x_abs < 0x29c6)) {
- int rounding = fputil::quick_get_round();
- if ((rounding == FE_UPWARD || rounding == FE_TOWARDZERO) && xf < 0)
- return fputil::cast<float16>(xf + 0x1p-24f);
- if ((rounding == FE_DOWNWARD || rounding == FE_TOWARDZERO) && xf > 0)
- return fputil::cast<float16>(xf - 0x1p-24f);
- return fputil::cast<float16>(xf);
- }
-
- float x_sq = xf * xf;
- // Generated by Sollya with:
- // > P = fpminimax(asinh(x)/x, [|0, 2, 4, 6, 8|], [|SG...|], [0, 2^-2]);
- // The last coefficient 0x1.bd114ep-6f has been changed to 0x1.bd114ep-5f
- // for better accuracy.
- float p = fputil::polyeval(x_sq, 1.0f, -0x1.555552p-3f, 0x1.332f6ap-4f,
- -0x1.6c53dep-5f, 0x1.bd114ep-5f);
-
- return fputil::cast<float16>(xf * p);
- }
-
- // General case: asinh(x) = ln(x + sqrt(x^2 + 1))
- float sqrt_term = fputil::sqrt<float>(fputil::multiply_add(xf, xf, 1.0f));
- return fputil::cast<float16>(
- x_sign * log_eval(fputil::multiply_add(xf, x_sign, sqrt_term)));
-}
+LLVM_LIBC_FUNCTION(float16, asinhf16, (float16 x)) { return math::asinhf16(x); }
} // namespace LIBC_NAMESPACE_DECL
diff --git a/libc/src/math/generic/atan.cpp b/libc/src/math/generic/atan.cpp
index cbca605..93bf2e1 100644
--- a/libc/src/math/generic/atan.cpp
+++ b/libc/src/math/generic/atan.cpp
@@ -7,173 +7,10 @@
//===----------------------------------------------------------------------===//
#include "src/math/atan.h"
-#include "atan_utils.h"
-#include "src/__support/FPUtil/FEnvImpl.h"
-#include "src/__support/FPUtil/FPBits.h"
-#include "src/__support/FPUtil/double_double.h"
-#include "src/__support/FPUtil/multiply_add.h"
-#include "src/__support/FPUtil/nearest_integer.h"
-#include "src/__support/macros/config.h"
-#include "src/__support/macros/optimization.h" // LIBC_UNLIKELY
+#include "src/__support/math/atan.h"
namespace LIBC_NAMESPACE_DECL {
-// To compute atan(x), we divided it into the following cases:
-// * |x| < 2^-26:
-// Since |x| > atan(|x|) > |x| - |x|^3/3, and |x|^3/3 < ulp(x)/2, we simply
-// return atan(x) = x - sign(x) * epsilon.
-// * 2^-26 <= |x| < 1:
-// We perform range reduction mod 2^-6 = 1/64 as follow:
-// Let k = 2^(-6) * round(|x| * 2^6), then
-// atan(x) = sign(x) * atan(|x|)
-// = sign(x) * (atan(k) + atan((|x| - k) / (1 + |x|*k)).
-// We store atan(k) in a look up table, and perform intermediate steps in
-// double-double.
-// * 1 < |x| < 2^53:
-// First we perform the transformation y = 1/|x|:
-// atan(x) = sign(x) * (pi/2 - atan(1/|x|))
-// = sign(x) * (pi/2 - atan(y)).
-// Then we compute atan(y) using range reduction mod 2^-6 = 1/64 as the
-// previous case:
-// Let k = 2^(-6) * round(y * 2^6), then
-// atan(y) = atan(k) + atan((y - k) / (1 + y*k))
-// = atan(k) + atan((1/|x| - k) / (1 + k/|x|)
-// = atan(k) + atan((1 - k*|x|) / (|x| + k)).
-// * |x| >= 2^53:
-// Using the reciprocal transformation:
-// atan(x) = sign(x) * (pi/2 - atan(1/|x|)).
-// We have that:
-// atan(1/|x|) <= 1/|x| <= 2^-53,
-// which is smaller than ulp(pi/2) / 2.
-// So we can return:
-// atan(x) = sign(x) * (pi/2 - epsilon)
-
-LLVM_LIBC_FUNCTION(double, atan, (double x)) {
- using FPBits = fputil::FPBits<double>;
-
- constexpr double IS_NEG[2] = {1.0, -1.0};
- constexpr DoubleDouble PI_OVER_2 = {0x1.1a62633145c07p-54,
- 0x1.921fb54442d18p0};
- constexpr DoubleDouble MPI_OVER_2 = {-0x1.1a62633145c07p-54,
- -0x1.921fb54442d18p0};
-
- FPBits xbits(x);
- bool x_sign = xbits.is_neg();
- xbits = xbits.abs();
- uint64_t x_abs = xbits.uintval();
- int x_exp =
- static_cast<int>(x_abs >> FPBits::FRACTION_LEN) - FPBits::EXP_BIAS;
-
- // |x| < 1.
- if (x_exp < 0) {
- if (LIBC_UNLIKELY(x_exp < -26)) {
-#ifdef LIBC_MATH_HAS_SKIP_ACCURATE_PASS
- return x;
-#else
- if (x == 0.0)
- return x;
- // |x| < 2^-26
- return fputil::multiply_add(-0x1.0p-54, x, x);
-#endif // LIBC_MATH_HAS_SKIP_ACCURATE_PASS
- }
-
- double x_d = xbits.get_val();
- // k = 2^-6 * round(2^6 * |x|)
- double k = fputil::nearest_integer(0x1.0p6 * x_d);
- unsigned idx = static_cast<unsigned>(k);
- k *= 0x1.0p-6;
-
- // numerator = |x| - k
- DoubleDouble num, den;
- num.lo = 0.0;
- num.hi = x_d - k;
-
- // denominator = 1 - k * |x|
- den.hi = fputil::multiply_add(x_d, k, 1.0);
- DoubleDouble prod = fputil::exact_mult(x_d, k);
- // Using Dekker's 2SUM algorithm to compute the lower part.
- den.lo = ((1.0 - den.hi) + prod.hi) + prod.lo;
-
- // x_r = (|x| - k) / (1 + k * |x|)
- DoubleDouble x_r = fputil::div(num, den);
-
- // Approximating atan(x_r) using Taylor polynomial.
- DoubleDouble p = atan_eval(x_r);
-
- // atan(x) = sign(x) * (atan(k) + atan(x_r))
- // = sign(x) * (atan(k) + atan( (|x| - k) / (1 + k * |x|) ))
-#ifdef LIBC_MATH_HAS_SKIP_ACCURATE_PASS
- return IS_NEG[x_sign] * (ATAN_I[idx].hi + (p.hi + (p.lo + ATAN_I[idx].lo)));
-#else
-
- DoubleDouble c0 = fputil::exact_add(ATAN_I[idx].hi, p.hi);
- double c1 = c0.lo + (ATAN_I[idx].lo + p.lo);
- double r = IS_NEG[x_sign] * (c0.hi + c1);
-
- return r;
-#endif // LIBC_MATH_HAS_SKIP_ACCURATE_PASS
- }
-
- // |x| >= 2^53 or x is NaN.
- if (LIBC_UNLIKELY(x_exp >= 53)) {
- // x is nan
- if (xbits.is_nan()) {
- if (xbits.is_signaling_nan()) {
- fputil::raise_except_if_required(FE_INVALID);
- return FPBits::quiet_nan().get_val();
- }
- return x;
- }
- // |x| >= 2^53
- // atan(x) ~ sign(x) * pi/2.
- if (x_exp >= 53)
-#ifdef LIBC_MATH_HAS_SKIP_ACCURATE_PASS
- return IS_NEG[x_sign] * PI_OVER_2.hi;
-#else
- return fputil::multiply_add(IS_NEG[x_sign], PI_OVER_2.hi,
- IS_NEG[x_sign] * PI_OVER_2.lo);
-#endif // LIBC_MATH_HAS_SKIP_ACCURATE_PASS
- }
-
- double x_d = xbits.get_val();
- double y = 1.0 / x_d;
-
- // k = 2^-6 * round(2^6 / |x|)
- double k = fputil::nearest_integer(0x1.0p6 * y);
- unsigned idx = static_cast<unsigned>(k);
- k *= 0x1.0p-6;
-
- // denominator = |x| + k
- DoubleDouble den = fputil::exact_add(x_d, k);
- // numerator = 1 - k * |x|
- DoubleDouble num;
- num.hi = fputil::multiply_add(-x_d, k, 1.0);
- DoubleDouble prod = fputil::exact_mult(x_d, k);
- // Using Dekker's 2SUM algorithm to compute the lower part.
- num.lo = ((1.0 - num.hi) - prod.hi) - prod.lo;
-
- // x_r = (1/|x| - k) / (1 - k/|x|)
- // = (1 - k * |x|) / (|x| - k)
- DoubleDouble x_r = fputil::div(num, den);
-
- // Approximating atan(x_r) using Taylor polynomial.
- DoubleDouble p = atan_eval(x_r);
-
- // atan(x) = sign(x) * (pi/2 - atan(1/|x|))
- // = sign(x) * (pi/2 - atan(k) - atan(x_r))
- // = (-sign(x)) * (-pi/2 + atan(k) + atan((1 - k*|x|)/(|x| - k)))
-#ifdef LIBC_MATH_HAS_SKIP_ACCURATE_PASS
- double lo_part = p.lo + ATAN_I[idx].lo + MPI_OVER_2.lo;
- return IS_NEG[!x_sign] * (MPI_OVER_2.hi + ATAN_I[idx].hi + (p.hi + lo_part));
-#else
- DoubleDouble c0 = fputil::exact_add(MPI_OVER_2.hi, ATAN_I[idx].hi);
- DoubleDouble c1 = fputil::exact_add(c0.hi, p.hi);
- double c2 = c1.lo + (c0.lo + p.lo) + (ATAN_I[idx].lo + MPI_OVER_2.lo);
-
- double r = IS_NEG[!x_sign] * (c1.hi + c2);
-
- return r;
-#endif
-}
+LLVM_LIBC_FUNCTION(double, atan, (double x)) { return math::atan(x); }
} // namespace LIBC_NAMESPACE_DECL
diff --git a/libc/src/math/generic/atan2.cpp b/libc/src/math/generic/atan2.cpp
index aa770de..58042d3 100644
--- a/libc/src/math/generic/atan2.cpp
+++ b/libc/src/math/generic/atan2.cpp
@@ -7,7 +7,6 @@
//===----------------------------------------------------------------------===//
#include "src/math/atan2.h"
-#include "atan_utils.h"
#include "src/__support/FPUtil/FEnvImpl.h"
#include "src/__support/FPUtil/FPBits.h"
#include "src/__support/FPUtil/double_double.h"
@@ -15,6 +14,7 @@
#include "src/__support/FPUtil/nearest_integer.h"
#include "src/__support/macros/config.h"
#include "src/__support/macros/optimization.h" // LIBC_UNLIKELY
+#include "src/__support/math/atan_utils.h"
namespace LIBC_NAMESPACE_DECL {
@@ -72,6 +72,7 @@ namespace LIBC_NAMESPACE_DECL {
// |(atan(u) - P(u)) / P(u)| < u^10 / 11 < 2^-73.
LLVM_LIBC_FUNCTION(double, atan2, (double y, double x)) {
+ using namespace atan_internal;
using FPBits = fputil::FPBits<double>;
constexpr double IS_NEG[2] = {1.0, -1.0};
diff --git a/libc/src/math/generic/atan2f128.cpp b/libc/src/math/generic/atan2f128.cpp
index a3aba0b..8838d94 100644
--- a/libc/src/math/generic/atan2f128.cpp
+++ b/libc/src/math/generic/atan2f128.cpp
@@ -7,7 +7,6 @@
//===----------------------------------------------------------------------===//
#include "src/math/atan2f128.h"
-#include "atan_utils.h"
#include "src/__support/FPUtil/FPBits.h"
#include "src/__support/FPUtil/dyadic_float.h"
#include "src/__support/FPUtil/multiply_add.h"
@@ -16,6 +15,7 @@
#include "src/__support/macros/config.h"
#include "src/__support/macros/optimization.h" // LIBC_UNLIKELY
#include "src/__support/macros/properties/types.h"
+#include "src/__support/math/atan_utils.h"
#include "src/__support/uint128.h"
namespace LIBC_NAMESPACE_DECL {
@@ -103,6 +103,7 @@ static constexpr Float128 CONST_ADJ[2][2][2] = {
// |(atan(u) - P(u)) / P(u)| < 2^-114.
LLVM_LIBC_FUNCTION(float128, atan2f128, (float128 y, float128 x)) {
+ using namespace atan_internal;
using FPBits = fputil::FPBits<float128>;
using Float128 = fputil::DyadicFloat<128>;
diff --git a/libc/src/math/generic/atan_utils.h b/libc/src/math/generic/atan_utils.h
deleted file mode 100644
index 24c7271..0000000
--- a/libc/src/math/generic/atan_utils.h
+++ /dev/null
@@ -1,241 +0,0 @@
-//===-- Collection of utils for atan/atan2 ----------------------*- 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_GENERIC_ATAN_UTILS_H
-#define LLVM_LIBC_SRC_MATH_GENERIC_ATAN_UTILS_H
-
-#include "src/__support/FPUtil/PolyEval.h"
-#include "src/__support/FPUtil/double_double.h"
-#include "src/__support/FPUtil/dyadic_float.h"
-#include "src/__support/FPUtil/multiply_add.h"
-#include "src/__support/integer_literals.h"
-#include "src/__support/macros/config.h"
-
-namespace LIBC_NAMESPACE_DECL {
-
-namespace {
-
-using DoubleDouble = fputil::DoubleDouble;
-using Float128 = fputil::DyadicFloat<128>;
-
-// atan(i/64) with i = 0..64, generated by Sollya with:
-// > for i from 0 to 64 do {
-// a = round(atan(i/64), D, RN);
-// b = round(atan(i/64) - a, D, RN);
-// print("{", b, ",", a, "},");
-// };
-constexpr DoubleDouble ATAN_I[65] = {
- {0.0, 0.0},
- {-0x1.220c39d4dff5p-61, 0x1.fff555bbb729bp-7},
- {-0x1.5ec431444912cp-60, 0x1.ffd55bba97625p-6},
- {-0x1.86ef8f794f105p-63, 0x1.7fb818430da2ap-5},
- {-0x1.c934d86d23f1dp-60, 0x1.ff55bb72cfdeap-5},
- {0x1.ac4ce285df847p-58, 0x1.3f59f0e7c559dp-4},
- {-0x1.cfb654c0c3d98p-58, 0x1.7ee182602f10fp-4},
- {0x1.f7b8f29a05987p-58, 0x1.be39ebe6f07c3p-4},
- {-0x1.cd37686760c17p-59, 0x1.fd5ba9aac2f6ep-4},
- {-0x1.b485914dacf8cp-59, 0x1.1e1fafb043727p-3},
- {0x1.61a3b0ce9281bp-57, 0x1.3d6eee8c6626cp-3},
- {-0x1.054ab2c010f3dp-58, 0x1.5c9811e3ec26ap-3},
- {0x1.347b0b4f881cap-58, 0x1.7b97b4bce5b02p-3},
- {0x1.cf601e7b4348ep-59, 0x1.9a6a8e96c8626p-3},
- {0x1.17b10d2e0e5abp-61, 0x1.b90d7529260a2p-3},
- {0x1.c648d1534597ep-57, 0x1.d77d5df205736p-3},
- {0x1.8ab6e3cf7afbdp-57, 0x1.f5b75f92c80ddp-3},
- {0x1.62e47390cb865p-56, 0x1.09dc597d86362p-2},
- {0x1.30ca4748b1bf9p-57, 0x1.18bf5a30bf178p-2},
- {-0x1.077cdd36dfc81p-56, 0x1.278372057ef46p-2},
- {-0x1.963a544b672d8p-57, 0x1.362773707ebccp-2},
- {-0x1.5d5e43c55b3bap-56, 0x1.44aa436c2af0ap-2},
- {-0x1.2566480884082p-57, 0x1.530ad9951cd4ap-2},
- {-0x1.a725715711fp-56, 0x1.614840309cfe2p-2},
- {-0x1.c63aae6f6e918p-56, 0x1.6f61941e4def1p-2},
- {0x1.69c885c2b249ap-56, 0x1.7d5604b63b3f7p-2},
- {0x1.b6d0ba3748fa8p-56, 0x1.8b24d394a1b25p-2},
- {0x1.9e6c988fd0a77p-56, 0x1.98cd5454d6b18p-2},
- {-0x1.24dec1b50b7ffp-56, 0x1.a64eec3cc23fdp-2},
- {0x1.ae187b1ca504p-56, 0x1.b3a911da65c6cp-2},
- {-0x1.cc1ce70934c34p-56, 0x1.c0db4c94ec9fp-2},
- {-0x1.a2cfa4418f1adp-56, 0x1.cde53432c1351p-2},
- {0x1.a2b7f222f65e2p-56, 0x1.dac670561bb4fp-2},
- {0x1.0e53dc1bf3435p-56, 0x1.e77eb7f175a34p-2},
- {-0x1.a3992dc382a23p-57, 0x1.f40dd0b541418p-2},
- {-0x1.b32c949c9d593p-55, 0x1.0039c73c1a40cp-1},
- {-0x1.d5b495f6349e6p-56, 0x1.0657e94db30dp-1},
- {0x1.974fa13b5404fp-58, 0x1.0c6145b5b43dap-1},
- {-0x1.2bdaee1c0ee35p-58, 0x1.1255d9bfbd2a9p-1},
- {0x1.c621cec00c301p-55, 0x1.1835a88be7c13p-1},
- {-0x1.928df287a668fp-58, 0x1.1e00babdefeb4p-1},
- {0x1.c421c9f38224ep-57, 0x1.23b71e2cc9e6ap-1},
- {-0x1.09e73b0c6c087p-56, 0x1.2958e59308e31p-1},
- {0x1.c5d5e9ff0cf8dp-55, 0x1.2ee628406cbcap-1},
- {0x1.1021137c71102p-55, 0x1.345f01cce37bbp-1},
- {-0x1.2304331d8bf46p-55, 0x1.39c391cd4171ap-1},
- {0x1.ecf8b492644fp-56, 0x1.3f13fb89e96f4p-1},
- {-0x1.f76d0163f79c8p-56, 0x1.445065b795b56p-1},
- {0x1.2419a87f2a458p-56, 0x1.4978fa3269ee1p-1},
- {0x1.4a33dbeb3796cp-55, 0x1.4e8de5bb6ec04p-1},
- {-0x1.1bb74abda520cp-55, 0x1.538f57b89061fp-1},
- {-0x1.5e5c9d8c5a95p-56, 0x1.587d81f732fbbp-1},
- {0x1.0028e4bc5e7cap-57, 0x1.5d58987169b18p-1},
- {-0x1.2b785350ee8c1p-57, 0x1.6220d115d7b8ep-1},
- {-0x1.6ea6febe8bbbap-56, 0x1.66d663923e087p-1},
- {-0x1.a80386188c50ep-55, 0x1.6b798920b3d99p-1},
- {-0x1.8c34d25aadef6p-56, 0x1.700a7c5784634p-1},
- {0x1.7b2a6165884a1p-59, 0x1.748978fba8e0fp-1},
- {0x1.406a08980374p-55, 0x1.78f6bbd5d315ep-1},
- {0x1.560821e2f3aa9p-55, 0x1.7d528289fa093p-1},
- {-0x1.bf76229d3b917p-56, 0x1.819d0b7158a4dp-1},
- {0x1.6b66e7fc8b8c3p-57, 0x1.85d69576cc2c5p-1},
- {-0x1.55b9a5e177a1bp-55, 0x1.89ff5ff57f1f8p-1},
- {-0x1.ec182ab042f61p-56, 0x1.8e17aa99cc05ep-1},
- {0x1.1a62633145c07p-55, 0x1.921fb54442d18p-1},
-};
-
-// Approximate atan(x) for |x| <= 2^-7.
-// Using degree-9 Taylor polynomial:
-// P = x - x^3/3 + x^5/5 -x^7/7 + x^9/9;
-// Then the absolute error is bounded by:
-// |atan(x) - P(x)| < |x|^11/11 < 2^(-7*11) / 11 < 2^-80.
-// And the relative error is bounded by:
-// |(atan(x) - P(x))/atan(x)| < |x|^10 / 10 < 2^-73.
-// For x = x_hi + x_lo, fully expand the polynomial and drop any terms less than
-// ulp(x_hi^3 / 3) gives us:
-// P(x) ~ x_hi - x_hi^3/3 + x_hi^5/5 - x_hi^7/7 + x_hi^9/9 +
-// + x_lo * (1 - x_hi^2 + x_hi^4)
-// Since p.lo is ~ x^3/3, the relative error from rounding is bounded by:
-// |(atan(x) - P(x))/atan(x)| < ulp(x^2) <= 2^(-14-52) = 2^-66.
-[[maybe_unused]] DoubleDouble atan_eval(const DoubleDouble &x) {
- DoubleDouble p;
- p.hi = x.hi;
- double x_hi_sq = x.hi * x.hi;
- // c0 ~ x_hi^2 * 1/5 - 1/3
- double c0 = fputil::multiply_add(x_hi_sq, 0x1.999999999999ap-3,
- -0x1.5555555555555p-2);
- // c1 ~ x_hi^2 * 1/9 - 1/7
- double c1 = fputil::multiply_add(x_hi_sq, 0x1.c71c71c71c71cp-4,
- -0x1.2492492492492p-3);
- // x_hi^3
- double x_hi_3 = x_hi_sq * x.hi;
- // x_hi^4
- double x_hi_4 = x_hi_sq * x_hi_sq;
- // d0 ~ 1/3 - x_hi^2 / 5 + x_hi^4 / 7 - x_hi^6 / 9
- double d0 = fputil::multiply_add(x_hi_4, c1, c0);
- // x_lo - x_lo * x_hi^2 + x_lo * x_hi^4
- double d1 = fputil::multiply_add(x_hi_4 - x_hi_sq, x.lo, x.lo);
- // p.lo ~ -x_hi^3/3 + x_hi^5/5 - x_hi^7/7 + x_hi^9/9 +
- // + x_lo * (1 - x_hi^2 + x_hi^4)
- p.lo = fputil::multiply_add(x_hi_3, d0, d1);
- return p;
-}
-
-// Float128 versions.
-// atan(i/64) with i = 0..64, generated by Sollya with:
-// > for i from 1 to 64 do {
-// a = round(atan(i/64), 128, RN);
-// ll = ceil(log2(a));
-// b = 2^ll + a;
-// print("{Sign::POS, ", 2^(ll - 128), ",", b, "},");
-// };
-constexpr Float128 ATAN_I_F128[65] = {
- {Sign::POS, 0, 0_u128},
- {Sign::POS, -134, 0xfffaaadd'db94d5bb'e78c5640'15f76048_u128},
- {Sign::POS, -133, 0xffeaaddd'4bb12542'779d776d'da8c6214_u128},
- {Sign::POS, -132, 0xbfdc0c21'86d14fcf'220e10d6'1df56ec7_u128},
- {Sign::POS, -132, 0xffaaddb9'67ef4e36'cb2792dc'0e2e0d51_u128},
- {Sign::POS, -131, 0x9facf873'e2aceb58'99c50bbf'08e6cdf6_u128},
- {Sign::POS, -131, 0xbf70c130'17887460'93567e78'4cf83676_u128},
- {Sign::POS, -131, 0xdf1cf5f3'783e1bef'71e5340b'30e5d9ef_u128},
- {Sign::POS, -131, 0xfeadd4d5'617b6e32'c897989f'3e888ef8_u128},
- {Sign::POS, -130, 0x8f0fd7d8'21b93725'bd375929'83a0af9a_u128},
- {Sign::POS, -130, 0x9eb77746'331362c3'47619d25'0360fe85_u128},
- {Sign::POS, -130, 0xae4c08f1'f6134efa'b54d3fef'0c2de994_u128},
- {Sign::POS, -130, 0xbdcbda5e'72d81134'7b0b4f88'1c9c7488_u128},
- {Sign::POS, -130, 0xcd35474b'643130e7'b00f3da1'a46eeb3b_u128},
- {Sign::POS, -130, 0xdc86ba94'93051022'f621a5c1'cb552f03_u128},
- {Sign::POS, -130, 0xebbeaef9'02b9b38c'91a2a68b'2fbd78e8_u128},
- {Sign::POS, -130, 0xfadbafc9'6406eb15'6dc79ef5'f7a217e6_u128},
- {Sign::POS, -129, 0x84ee2cbe'c31b12c5'c8e72197'0cabd3a3_u128},
- {Sign::POS, -129, 0x8c5fad18'5f8bc130'ca4748b1'bf88298d_u128},
- {Sign::POS, -129, 0x93c1b902'bf7a2df1'06459240'6fe1447a_u128},
- {Sign::POS, -129, 0x9b13b9b8'3f5e5e69'c5abb498'd27af328_u128},
- {Sign::POS, -129, 0xa25521b6'15784d45'43787549'88b8d9e3_u128},
- {Sign::POS, -129, 0xa9856cca'8e6a4eda'99b7f77b'f7d9e8c1_u128},
- {Sign::POS, -129, 0xb0a42018'4e7f0cb1'b51d51dc'200a0fc3_u128},
- {Sign::POS, -129, 0xb7b0ca0f'26f78473'8aa32122'dcfe4483_u128},
- {Sign::POS, -129, 0xbeab025b'1d9fbad3'910b8564'93411026_u128},
- {Sign::POS, -129, 0xc59269ca'50d92b6d'a1746e91'f50a28de_u128},
- {Sign::POS, -129, 0xcc66aa2a'6b58c33c'd9311fa1'4ed9b7c4_u128},
- {Sign::POS, -129, 0xd327761e'611fe5b6'427c95e9'001e7136_u128},
- {Sign::POS, -129, 0xd9d488ed'32e3635c'30f6394a'0806345d_u128},
- {Sign::POS, -129, 0xe06da64a'764f7c67'c631ed96'798cb804_u128},
- {Sign::POS, -129, 0xe6f29a19'609a84ba'60b77ce1'ca6dc2c8_u128},
- {Sign::POS, -129, 0xed63382b'0dda7b45'6fe445ec'bc3a8d03_u128},
- {Sign::POS, -129, 0xf3bf5bf8'bad1a21c'a7b837e6'86adf3fa_u128},
- {Sign::POS, -129, 0xfa06e85a'a0a0be5c'66d23c7d'5dc8ecc2_u128},
- {Sign::POS, -128, 0x801ce39e'0d205c99'a6d6c6c5'4d938596_u128},
- {Sign::POS, -128, 0x832bf4a6'd9867e2a'4b6a09cb'61a515c1_u128},
- {Sign::POS, -128, 0x8630a2da'da1ed065'd3e84ed5'013ca37e_u128},
- {Sign::POS, -128, 0x892aecdf'de9547b5'094478fc'472b4afc_u128},
- {Sign::POS, -128, 0x8c1ad445'f3e09b8c'439d8018'60205921_u128},
- {Sign::POS, -128, 0x8f005d5e'f7f59f9b'5c835e16'65c43748_u128},
- {Sign::POS, -128, 0x91db8f16'64f350e2'10e4f9c1'126e0220_u128},
- {Sign::POS, -128, 0x94ac72c9'847186f6'18c4f393'f78a32f9_u128},
- {Sign::POS, -128, 0x97731420'365e538b'abd3fe19'f1aeb6b3_u128},
- {Sign::POS, -128, 0x9a2f80e6'71bdda20'4226f8e2'204ff3bd_u128},
- {Sign::POS, -128, 0x9ce1c8e6'a0b8cdb9'f799c4e8'174cf11c_u128},
- {Sign::POS, -128, 0x9f89fdc4'f4b7a1ec'f8b49264'4f0701e0_u128},
- {Sign::POS, -128, 0xa22832db'cadaae08'92fe9c08'637af0e6_u128},
- {Sign::POS, -128, 0xa4bc7d19'34f70924'19a87f2a'457dac9f_u128},
- {Sign::POS, -128, 0xa746f2dd'b7602294'67b7d66f'2d74e019_u128},
- {Sign::POS, -128, 0xa9c7abdc'4830f5c8'916a84b5'be7933f6_u128},
- {Sign::POS, -128, 0xac3ec0fb'997dd6a1'a36273a5'6afa8ef4_u128},
- {Sign::POS, -128, 0xaeac4c38'b4d8c080'14725e2f'3e52070a_u128},
- {Sign::POS, -128, 0xb110688a'ebdc6f6a'43d65788'b9f6a7b5_u128},
- {Sign::POS, -128, 0xb36b31c9'1f043691'59014174'4462f93a_u128},
- {Sign::POS, -128, 0xb5bcc490'59ecc4af'f8f3cee7'5e3907d5_u128},
- {Sign::POS, -128, 0xb8053e2b'c2319e73'cb2da552'10a4443d_u128},
- {Sign::POS, -128, 0xba44bc7d'd470782f'654c2cb1'0942e386_u128},
- {Sign::POS, -128, 0xbc7b5dea'e98af280'd4113006'e80fb290_u128},
- {Sign::POS, -128, 0xbea94144'fd049aac'1043c5e7'55282e7d_u128},
- {Sign::POS, -128, 0xc0ce85b8'ac526640'89dd62c4'6e92fa25_u128},
- {Sign::POS, -128, 0xc2eb4abb'661628b5'b373fe45'c61bb9fb_u128},
- {Sign::POS, -128, 0xc4ffaffa'bf8fbd54'8cb43d10'bc9e0221_u128},
- {Sign::POS, -128, 0xc70bd54c'e602ee13'e7d54fbd'09f2be38_u128},
- {Sign::POS, -128, 0xc90fdaa2'2168c234'c4c6628b'80dc1cd1_u128},
-};
-
-// Degree-13 minimax polynomial generated by Sollya with:
-// > P = fpminimax(atan(x), [|1, 3, 5, 7, 9, 11, 13|], [|1, 128...|],
-// [0, 2^-7]);
-// > dirtyinfnorm(atan(x) - P, [0, 2^-7]);
-// 0x1.26016ad97f323875760f869684c0898d7b7bb8bep-122
-constexpr Float128 ATAN_POLY_F128[] = {
- {Sign::NEG, -129, 0xaaaaaaaa'aaaaaaaa'aaaaaaa6'003c5d1d_u128},
- {Sign::POS, -130, 0xcccccccc'cccccccc'cca00232'8776b063_u128},
- {Sign::NEG, -130, 0x92492492'49249201'27f5268a'cb24aec0_u128},
- {Sign::POS, -131, 0xe38e38e3'8dce3d96'626a1643'f8eb68f3_u128},
- {Sign::NEG, -131, 0xba2e8b7a'ea4ad00f'005a35c7'6ef609b1_u128},
- {Sign::POS, -131, 0x9d82765e'd22a7d92'ac09c405'c0a69214_u128},
-};
-
-// Approximate atan for |x| <= 2^-7.
-[[maybe_unused]] Float128 atan_eval(const Float128 &x) {
- Float128 x_sq = fputil::quick_mul(x, x);
- Float128 x3 = fputil::quick_mul(x, x_sq);
- Float128 p = fputil::polyeval(x_sq, ATAN_POLY_F128[0], ATAN_POLY_F128[1],
- ATAN_POLY_F128[2], ATAN_POLY_F128[3],
- ATAN_POLY_F128[4], ATAN_POLY_F128[5]);
- return fputil::multiply_add(x3, p, x);
-}
-
-} // anonymous namespace
-
-} // namespace LIBC_NAMESPACE_DECL
-
-#endif // LLVM_LIBC_SRC_MATH_GENERIC_ATAN_UTILS_H
diff --git a/libc/src/math/generic/fabsbf16.cpp b/libc/src/math/generic/fabsbf16.cpp
new file mode 100644
index 0000000..ea39719
--- /dev/null
+++ b/libc/src/math/generic/fabsbf16.cpp
@@ -0,0 +1,19 @@
+//===-- Implementation of fabsbf16 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/fabsbf16.h"
+
+#include "src/__support/FPUtil/BasicOperations.h"
+#include "src/__support/FPUtil/bfloat16.h"
+#include "src/__support/macros/config.h"
+
+namespace LIBC_NAMESPACE_DECL {
+
+LLVM_LIBC_FUNCTION(bfloat16, fabsbf16, (bfloat16 x)) { return fputil::abs(x); }
+
+} // namespace LIBC_NAMESPACE_DECL