[libc][math] Refactor asinf16 implementation to header-only in src/__support/math folder. #150800
Conversation
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
…support/math folder.
bassiounix
requested review from
aaronmondal,
keith and
rupprecht
as code owners
|
@llvm/pr-subscribers-libc Author: Muhammad Bassiouni (bassiounix) ChangesPart of #147386 in preparation for: https://discourse.llvm.org/t/rfc-make-clang-builtin-math-functions-constexpr-with-llvm-libc-to-support-c-23-constexpr-math-functions/86450 Full diff: https://github.com/llvm/llvm-project/pull/150800.diff 9 Files Affected:
diff --git a/libc/shared/math.h b/libc/shared/math.h
index 042daf65e45bd..e01c8db81eb80 100644
--- a/libc/shared/math.h
+++ b/libc/shared/math.h
@@ -19,6 +19,7 @@
#include "math/acospif16.h"
#include "math/asin.h"
#include "math/asinf.h"
+#include "math/asinf16.h"
#include "math/erff.h"
#include "math/exp.h"
#include "math/exp10.h"
diff --git a/libc/shared/math/asinf16.h b/libc/shared/math/asinf16.h
new file mode 100644
index 0000000000000..d545e269a6402
--- /dev/null
+++ b/libc/shared/math/asinf16.h
@@ -0,0 +1,28 @@
+//===-- Shared asinf16 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_ASINF16_H
+#define LLVM_LIBC_SHARED_MATH_ASINF16_H
+
+#include "shared/libc_common.h"
+
+#ifdef LIBC_TYPES_HAS_FLOAT16
+
+#include "src/__support/math/asinf16.h"
+
+namespace LIBC_NAMESPACE_DECL {
+namespace shared {
+
+using math::asinf16;
+
+} // namespace shared
+} // namespace LIBC_NAMESPACE_DECL
+
+#endif // LIBC_TYPES_HAS_FLOAT16
+
+#endif // LLVM_LIBC_SHARED_MATH_ASINF_H
diff --git a/libc/src/__support/math/CMakeLists.txt b/libc/src/__support/math/CMakeLists.txt
index b096c616b0927..105093887980d 100644
--- a/libc/src/__support/math/CMakeLists.txt
+++ b/libc/src/__support/math/CMakeLists.txt
@@ -156,6 +156,20 @@ add_header_library(
libc.src.__support.macros.properties.cpu_features
)
+add_header_library(
+ asinf16
+ HDRS
+ asinf16.h
+ DEPENDS
+ libc.src.__support.FPUtil.fenv_impl
+ libc.src.__support.FPUtil.fp_bits
+ libc.src.__support.FPUtil.polyeval
+ libc.src.__support.FPUtil.cast
+ libc.src.__support.FPUtil.multiply_add
+ libc.src.__support.FPUtil.sqrt
+ libc.src.__support.macros.optimization
+)
+
add_header_library(
erff
HDRS
diff --git a/libc/src/__support/math/asinf16.h b/libc/src/__support/math/asinf16.h
new file mode 100644
index 0000000000000..3d032a41a2bb3
--- /dev/null
+++ b/libc/src/__support/math/asinf16.h
@@ -0,0 +1,146 @@
+//===-- Implementation header for asinf16 -----------------------*- 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_ASINF16_H
+#define LLVM_LIBC_SRC___SUPPORT_MATH_ASINF16_H
+
+#include "include/llvm-libc-macros/float16-macros.h"
+
+#ifdef LIBC_TYPES_HAS_FLOAT16
+
+#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"
+
+namespace LIBC_NAMESPACE_DECL {
+
+namespace math {
+
+LIBC_INLINE static constexpr float16 asinf16(float16 x) {
+
+ // Generated by Sollya using the following command:
+ // > round(pi/2, D, RN);
+ constexpr float PI_2 = 0x1.921fb54442d18p0f;
+
+ 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));
+}
+
+} // namespace math
+
+} // namespace LIBC_NAMESPACE_DECL
+
+#endif // LIBC_TYPES_HAS_FLOAT16
+
+#endif // LLVM_LIBC_SRC___SUPPORT_MATH_ASINF16_H
diff --git a/libc/src/math/generic/CMakeLists.txt b/libc/src/math/generic/CMakeLists.txt
index ecf09674b4d4c..d4d268cff97e8 100644
--- a/libc/src/math/generic/CMakeLists.txt
+++ b/libc/src/math/generic/CMakeLists.txt
@@ -3968,16 +3968,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(
diff --git a/libc/src/math/generic/asinf16.cpp b/libc/src/math/generic/asinf16.cpp
index 518c384a61530..af8dbfebef843 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/test/shared/CMakeLists.txt b/libc/test/shared/CMakeLists.txt
index e5dfac98b5e69..e93fbb9ffbbdf 100644
--- a/libc/test/shared/CMakeLists.txt
+++ b/libc/test/shared/CMakeLists.txt
@@ -15,6 +15,7 @@ add_fp_unittest(
libc.src.__support.math.acospif16
libc.src.__support.math.asin
libc.src.__support.math.asinf
+ libc.src.__support.math.asinf16
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 7881d6841087b..51d38ecc22300 100644
--- a/libc/test/shared/shared_math_test.cpp
+++ b/libc/test/shared/shared_math_test.cpp
@@ -17,6 +17,7 @@ TEST(LlvmLibcSharedMathTest, AllFloat16) {
EXPECT_FP_EQ(0x0p+0f16, LIBC_NAMESPACE::shared::acoshf16(1.0f16));
EXPECT_FP_EQ(0x0p+0f16, LIBC_NAMESPACE::shared::acospif16(1.0f16));
+ EXPECT_FP_EQ(0x0p+0f16, LIBC_NAMESPACE::shared::asinf16(0.0f16));
EXPECT_FP_EQ(0x1p+0f16, LIBC_NAMESPACE::shared::exp10f16(0.0f16));
diff --git a/utils/bazel/llvm-project-overlay/libc/BUILD.bazel b/utils/bazel/llvm-project-overlay/libc/BUILD.bazel
index 25353e18845e1..d424026c2723f 100644
--- a/utils/bazel/llvm-project-overlay/libc/BUILD.bazel
+++ b/utils/bazel/llvm-project-overlay/libc/BUILD.bazel
@@ -2247,6 +2247,19 @@ libc_support_library(
],
)
+libc_support_library(
+ name = "__support_math_asinf16",
+ hdrs = ["src/__support/math/asinf16.h"],
+ deps = [
+ ":__support_fputil_fenv_impl",
+ ":__support_fputil_fp_bits",
+ ":__support_fputil_polyeval",
+ ":__support_fputil_multiply_add",
+ ":__support_fputil_sqrt",
+ ":__support_macros_optimization",
+ ],
+)
+
libc_support_library(
name = "__support_math_erff",
hdrs = ["src/__support/math/erff.h"],
@@ -2809,14 +2822,7 @@ libc_math_function(
libc_math_function(
name = "asinf16",
additional_deps = [
- ":__support_fputil_cast",
- ":__support_fputil_fma",
- ":__support_fputil_multiply_add",
- ":__support_fputil_nearest_integer",
- ":__support_fputil_polyeval",
- ":__support_fputil_sqrt",
- ":__support_macros_optimization",
- ":__support_math_inv_trigf_utils",
+ ":__support_math_asinf16",
],
)
|
bassiounix
requested review from
lntue
and removed request for
aaronmondal,
keith and
rupprecht
lntue
approved these changes
Jul 27, 2025
|
LLVM Buildbot has detected a new failure on builder Full details are available at: https://lab.llvm.org/buildbot/#/builders/11/builds/20482 Here is the relevant piece of the build log for the reference |
mahesh-attarde
pushed a commit
to mahesh-attarde/llvm-project
that referenced
this pull request
Jul 28, 2025
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
4 participants
Add this suggestion to a batch that can be applied as a single commit.
This suggestion is invalid because no changes were made to the code.
Suggestions cannot be applied while the pull request is closed.
Suggestions cannot be applied while viewing a subset of changes.
Only one suggestion per line can be applied in a batch.
Add this suggestion to a batch that can be applied as a single commit.
Applying suggestions on deleted lines is not supported.
You must change the existing code in this line in order to create a valid suggestion.
Outdated suggestions cannot be applied.
This suggestion has been applied or marked resolved.
Suggestions cannot be applied from pending reviews.
Suggestions cannot be applied on multi-line comments.
Suggestions cannot be applied while the pull request is queued to merge.
Suggestion cannot be applied right now. Please check back later.
Part of #147386
in preparation for: https://discourse.llvm.org/t/rfc-make-clang-builtin-math-functions-constexpr-with-llvm-libc-to-support-c-23-constexpr-math-functions/86450