[libc-commits] [libc] [libc][math][c23] Add rsqrtf16() function (PR #137545)
via libc-commits
libc-commits at lists.llvm.org
Sun Apr 27 13:08:24 PDT 2025
https://github.com/amemov updated https://github.com/llvm/llvm-project/pull/137545
>From 43fefedc6b8c644fb71aa5a12b17e9c76aac6cef Mon Sep 17 00:00:00 2001
From: amemov <shepelev777 at gmail.com>
Date: Sun, 27 Apr 2025 19:54:02 +0000
Subject: [PATCH 1/2] Initial implementation of rsqrtf16 - Defined and declared
rsqrtf16 - Added entrypoints to support rsqrt in float16 format - Added tests
- both unit and exhaustive - Added MPFR support for Rsqrt to verify the test
results
TODO:
- Write computation as a polynomial in the end of rsqrtf16.cpp
- Check if errors are generated correctly (Reference from C23: The rsqrt functions compute the reciprocal of the nonnegative square root of the argument. A domain error occurs if the argument is less than zero. A pole error may occur if the argument equals zero.)
---
libc/config/linux/x86_64/entrypoints.txt | 1 +
libc/docs/headers/math/index.rst | 2 +-
libc/include/math.yaml | 7 +++
libc/src/math/CMakeLists.txt | 2 +
libc/src/math/generic/CMakeLists.txt | 19 ++++++
libc/src/math/generic/rsqrtf16.cpp | 67 ++++++++++++++++++++++
libc/src/math/rsqrtf16.h | 21 +++++++
libc/test/src/math/CMakeLists.txt | 11 ++++
libc/test/src/math/rsqrtf16_test.cpp | 42 ++++++++++++++
libc/test/src/math/smoke/CMakeLists.txt | 11 ++++
libc/test/src/math/smoke/rsqrtf16_test.cpp | 38 ++++++++++++
libc/utils/MPFRWrapper/MPCommon.cpp | 6 ++
libc/utils/MPFRWrapper/MPCommon.h | 1 +
libc/utils/MPFRWrapper/MPFRUtils.cpp | 2 +
libc/utils/MPFRWrapper/MPFRUtils.h | 1 +
15 files changed, 230 insertions(+), 1 deletion(-)
create mode 100644 libc/src/math/generic/rsqrtf16.cpp
create mode 100644 libc/src/math/rsqrtf16.h
create mode 100644 libc/test/src/math/rsqrtf16_test.cpp
create mode 100644 libc/test/src/math/smoke/rsqrtf16_test.cpp
diff --git a/libc/config/linux/x86_64/entrypoints.txt b/libc/config/linux/x86_64/entrypoints.txt
index 0e4fbdd2f2f57..dc3120ca913a3 100644
--- a/libc/config/linux/x86_64/entrypoints.txt
+++ b/libc/config/linux/x86_64/entrypoints.txt
@@ -740,6 +740,7 @@ if(LIBC_TYPES_HAS_FLOAT16)
libc.src.math.rintf16
libc.src.math.roundevenf16
libc.src.math.roundf16
+ libc.src.math.rsqrtf16
libc.src.math.scalblnf16
libc.src.math.scalbnf16
libc.src.math.setpayloadf16
diff --git a/libc/docs/headers/math/index.rst b/libc/docs/headers/math/index.rst
index cf88e6237d1e3..a29e70a7b625b 100644
--- a/libc/docs/headers/math/index.rst
+++ b/libc/docs/headers/math/index.rst
@@ -333,7 +333,7 @@ Higher Math Functions
+-----------+------------------+-----------------+------------------------+----------------------+------------------------+------------------------+----------------------------+
| rootn | | | | | | 7.12.7.8 | F.10.4.8 |
+-----------+------------------+-----------------+------------------------+----------------------+------------------------+------------------------+----------------------------+
-| rsqrt | | | | | | 7.12.7.9 | F.10.4.9 |
+| rsqrt | | | | |check| | | 7.12.7.9 | F.10.4.9 |
+-----------+------------------+-----------------+------------------------+----------------------+------------------------+------------------------+----------------------------+
| sin | |check| | |check| | | |check| | | 7.12.4.6 | F.10.1.6 |
+-----------+------------------+-----------------+------------------------+----------------------+------------------------+------------------------+----------------------------+
diff --git a/libc/include/math.yaml b/libc/include/math.yaml
index b725c33c0bb06..35aee5f54d8cd 100644
--- a/libc/include/math.yaml
+++ b/libc/include/math.yaml
@@ -2237,6 +2237,13 @@ functions:
return_type: long double
arguments:
- type: long double
+ - name: rsqrtf16
+ standards:
+ - stdc
+ return_type: _Float16
+ arguments:
+ - type: _Float16
+ guard: LIBC_TYPES_HAS_FLOAT16
- name: scalbln
standards:
- stdc
diff --git a/libc/src/math/CMakeLists.txt b/libc/src/math/CMakeLists.txt
index 5161b2b61aa94..689327b2f1249 100644
--- a/libc/src/math/CMakeLists.txt
+++ b/libc/src/math/CMakeLists.txt
@@ -467,6 +467,8 @@ add_math_entrypoint_object(roundevenl)
add_math_entrypoint_object(roundevenf16)
add_math_entrypoint_object(roundevenf128)
+add_math_entrypoint_object(rsqrtf16)
+
add_math_entrypoint_object(scalbln)
add_math_entrypoint_object(scalblnf)
add_math_entrypoint_object(scalblnl)
diff --git a/libc/src/math/generic/CMakeLists.txt b/libc/src/math/generic/CMakeLists.txt
index 25f583035fbc2..ae0977ce255c2 100644
--- a/libc/src/math/generic/CMakeLists.txt
+++ b/libc/src/math/generic/CMakeLists.txt
@@ -955,6 +955,25 @@ add_entrypoint_object(
libc.src.__support.FPUtil.nearest_integer_operations
)
+add_entrypoint_object(
+ rsqrtf16
+ SRCS
+ rsqrtf16.cpp
+ HDRS
+ ../rsqrtf16.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
+)
+
add_entrypoint_object(
lround
SRCS
diff --git a/libc/src/math/generic/rsqrtf16.cpp b/libc/src/math/generic/rsqrtf16.cpp
new file mode 100644
index 0000000000000..930c8f6fc7b4e
--- /dev/null
+++ b/libc/src/math/generic/rsqrtf16.cpp
@@ -0,0 +1,67 @@
+//===-- Half-precision rsqrt 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/rsqrtf16.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"
+
+namespace LIBC_NAMESPACE_DECL {
+
+LLVM_LIBC_FUNCTION(float16, rsqrtf16, (float16 x)) {
+ using FPBits = fputil::FPBits<float16>;
+ FPBits xbits(x);
+
+ uint16_t x_u = xbits.uintval();
+ uint16_t x_abs = x_u & 0x7fff;
+ uint16_t x_sign = x_u >> 15;
+
+ // x is NaN
+ if (LIBC_UNLIKELY(xbits.is_nan())) {
+ if (xbits.is_signaling_nan()) {
+ fputil::raise_except_if_required(FE_INVALID);
+ return FPBits::quiet_nan().get_val();
+ }
+ return x;
+ }
+
+ // |x| = 0
+ if (LIBC_UNLIKELY(x_abs == 0x0)) {
+ fputil::raise_except_if_required(FE_DIVBYZERO);
+ fputil::set_errno_if_required(ERANGE);
+ return FPBits::quiet_nan().get_val();
+ }
+
+ // -inf <= x < 0
+ if (LIBC_UNLIKELY(x_sign == 1)) {
+ fputil::raise_except_if_required(FE_INVALID);
+ fputil::set_errno_if_required(EDOM);
+ return FPBits::quiet_nan().get_val();
+ }
+
+ // x = +inf => rsqrt(x) = 0
+ if (LIBC_UNLIKELY(xbits.is_inf())) {
+ return fputil::cast<float16>(0.0f);
+ }
+
+ // x = 1 => rsqrt(x) = 1
+ if (LIBC_UNLIKELY(x_u == 0x1)) {
+ return fputil::cast<float16>(1.0f);
+ }
+
+ // x is valid, estimate the result - below is temporary solution for just testing
+ float xf = x;
+ return fputil::cast<float16>(1.0f / xf);
+}
+} // namespace LIBC_NAMESPACE_DECL
diff --git a/libc/src/math/rsqrtf16.h b/libc/src/math/rsqrtf16.h
new file mode 100644
index 0000000000000..c88ab5256ce88
--- /dev/null
+++ b/libc/src/math/rsqrtf16.h
@@ -0,0 +1,21 @@
+//===-- Implementation header for rsqrtf16 ----------------------*- 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_RSQRTF16_H
+#define LLVM_LIBC_SRC_MATH_RSQRTF16_H
+
+#include "src/__support/macros/config.h"
+#include "src/__support/macros/properties/types.h"
+
+namespace LIBC_NAMESPACE_DECL {
+
+float16 rsqrtf16(float16 x);
+
+} // namespace LIBC_NAMESPACE_DECL
+
+#endif // LLVM_LIBC_SRC_MATH_RSQRTF16_H
diff --git a/libc/test/src/math/CMakeLists.txt b/libc/test/src/math/CMakeLists.txt
index ebf9f1c86cf15..8daf59032622c 100644
--- a/libc/test/src/math/CMakeLists.txt
+++ b/libc/test/src/math/CMakeLists.txt
@@ -1560,6 +1560,17 @@ add_fp_unittest(
libc.src.math.sqrtl
)
+add_fp_unittest(
+ rsqrtf16_test
+ NEED_MPFR
+ SUITE
+ libc-math-unittests
+ SRCS
+ rsqrtf16_test.cpp
+ DEPENDS
+ libc.src.math.rsqrtf16
+)
+
add_fp_unittest(
sqrtf16_test
NEED_MPFR
diff --git a/libc/test/src/math/rsqrtf16_test.cpp b/libc/test/src/math/rsqrtf16_test.cpp
new file mode 100644
index 0000000000000..ec22df5bd31c7
--- /dev/null
+++ b/libc/test/src/math/rsqrtf16_test.cpp
@@ -0,0 +1,42 @@
+//===-- Exhaustive test for rsqrtf16 --------------------------------------===//
+//
+// 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/rsqrtf16.h"
+#include "test/UnitTest/FPMatcher.h"
+#include "test/UnitTest/Test.h"
+#include "utils/MPFRWrapper/MPFRUtils.h"
+
+using LlvmLibcRsqrtf16Test = LIBC_NAMESPACE::testing::FPTest<float16>;
+
+namespace mpfr = LIBC_NAMESPACE::testing::mpfr;
+
+// Range: [0, Inf]
+static constexpr uint16_t POS_START = 0x0000U;
+static constexpr uint16_t POS_STOP = 0x7c00U;
+
+// Range: [-Inf, 0]
+static constexpr uint16_t NEG_START = 0x8000U;
+static constexpr uint16_t NEG_STOP = 0xfc00U;
+
+TEST_F(LlvmLibcRsqrtf16Test, PositiveRange) {
+ for (uint16_t v = POS_START; v <= POS_STOP; ++v) {
+ float16 x = FPBits(v).get_val();
+
+ EXPECT_MPFR_MATCH_ALL_ROUNDING(mpfr::Operation::Rsqrt, x,
+ LIBC_NAMESPACE::rsqrtf16(x), 0.5);
+ }
+}
+
+TEST_F(LlvmLibcRsqrtf16Test, NegativeRange) {
+ for (uint16_t v = NEG_START; v <= NEG_STOP; ++v) {
+ float16 x = FPBits(v).get_val();
+
+ EXPECT_MPFR_MATCH_ALL_ROUNDING(mpfr::Operation::Rsqrt, x,
+ LIBC_NAMESPACE::rsqrtf16(x), 0.5);
+ }
+}
\ No newline at end of file
diff --git a/libc/test/src/math/smoke/CMakeLists.txt b/libc/test/src/math/smoke/CMakeLists.txt
index 79b189159e9d8..fe5d8a6b3b864 100644
--- a/libc/test/src/math/smoke/CMakeLists.txt
+++ b/libc/test/src/math/smoke/CMakeLists.txt
@@ -2958,6 +2958,17 @@ add_fp_unittest(
libc.src.math.sqrtl
)
+add_fp_unittest(
+ rsqrtf16_test
+ SUITE
+ libc-math-smoke-tests
+ SRCS
+ rsqrtf16_test.cpp
+ DEPENDS
+ libc.src.errno.errno
+ libc.src.math.rsqrtf16
+)
+
add_fp_unittest(
sqrtf16_test
SUITE
diff --git a/libc/test/src/math/smoke/rsqrtf16_test.cpp b/libc/test/src/math/smoke/rsqrtf16_test.cpp
new file mode 100644
index 0000000000000..eb429d3dbf6c8
--- /dev/null
+++ b/libc/test/src/math/smoke/rsqrtf16_test.cpp
@@ -0,0 +1,38 @@
+//===-- Unittests for rsqrtf16 --------------------------------------------===//
+//
+// 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/errno/libc_errno.h"
+#include "src/math/rsqrtf16.h"
+#include "test/UnitTest/FPMatcher.h"
+#include "test/UnitTest/Test.h"
+
+using LlvmLibcRsqrtf16Test = LIBC_NAMESPACE::testing::FPTest<float16>;
+TEST_F(LlvmLibcRsqrtf16Test, SpecialNumbers) {
+ LIBC_NAMESPACE::libc_errno = 0;
+ EXPECT_FP_EQ(aNaN, LIBC_NAMESPACE::rsqrtf16(aNaN));
+ EXPECT_MATH_ERRNO(0);
+
+ EXPECT_FP_EQ_WITH_EXCEPTION(aNaN, LIBC_NAMESPACE::rsqrtf16(sNaN),
+ FE_INVALID);
+ EXPECT_MATH_ERRNO(0);
+
+ EXPECT_FP_EQ(aNaN, LIBC_NAMESPACE::rsqrtf16(0.0f));
+ EXPECT_MATH_ERRNO(ERANGE);
+
+ EXPECT_FP_EQ(1.0f, LIBC_NAMESPACE::rsqrtf16(1.0f));
+ EXPECT_MATH_ERRNO(0);
+
+ EXPECT_FP_EQ(0.0f, LIBC_NAMESPACE::rsqrtf16(inf));
+ EXPECT_MATH_ERRNO(0);
+
+ EXPECT_FP_EQ(aNaN, LIBC_NAMESPACE::rsqrtf16(neg_inf));
+ EXPECT_MATH_ERRNO(EDOM);
+
+ EXPECT_FP_EQ(aNaN, LIBC_NAMESPACE::rsqrtf16(-2.0f));
+ EXPECT_MATH_ERRNO(EDOM);
+}
\ No newline at end of file
diff --git a/libc/utils/MPFRWrapper/MPCommon.cpp b/libc/utils/MPFRWrapper/MPCommon.cpp
index ccd4d2d01a4e2..1a78ca5c24ba2 100644
--- a/libc/utils/MPFRWrapper/MPCommon.cpp
+++ b/libc/utils/MPFRWrapper/MPCommon.cpp
@@ -366,6 +366,12 @@ MPFRNumber MPFRNumber::rint(mpfr_rnd_t rnd) const {
return result;
}
+MPFRNumber MPFRNumber::rsqrt() const {
+ MPFRNumber result(*this);
+ mpfr_rec_sqrt(result.value, value, mpfr_rounding);
+ return result;
+}
+
MPFRNumber MPFRNumber::mod_2pi() const {
MPFRNumber result(0.0, 1280);
MPFRNumber _2pi(0.0, 1280);
diff --git a/libc/utils/MPFRWrapper/MPCommon.h b/libc/utils/MPFRWrapper/MPCommon.h
index 99cb7ec66a2ca..43218ee7662db 100644
--- a/libc/utils/MPFRWrapper/MPCommon.h
+++ b/libc/utils/MPFRWrapper/MPCommon.h
@@ -216,6 +216,7 @@ class MPFRNumber {
bool round_to_long(long &result) const;
bool round_to_long(mpfr_rnd_t rnd, long &result) const;
MPFRNumber rint(mpfr_rnd_t rnd) const;
+ MPFRNumber rsqrt() const;
MPFRNumber mod_2pi() const;
MPFRNumber mod_pi_over_2() const;
MPFRNumber mod_pi_over_4() const;
diff --git a/libc/utils/MPFRWrapper/MPFRUtils.cpp b/libc/utils/MPFRWrapper/MPFRUtils.cpp
index 8853f96ef8f92..a68d21650f004 100644
--- a/libc/utils/MPFRWrapper/MPFRUtils.cpp
+++ b/libc/utils/MPFRWrapper/MPFRUtils.cpp
@@ -86,6 +86,8 @@ unary_operation(Operation op, InputType input, unsigned int precision,
return mpfrInput.round();
case Operation::RoundEven:
return mpfrInput.roundeven();
+ case Operation::Rsqrt:
+ return mpfrInput.rsqrt();
case Operation::Sin:
return mpfrInput.sin();
case Operation::Sinpi:
diff --git a/libc/utils/MPFRWrapper/MPFRUtils.h b/libc/utils/MPFRWrapper/MPFRUtils.h
index c77a6aa3adeae..532b2c8d9d819 100644
--- a/libc/utils/MPFRWrapper/MPFRUtils.h
+++ b/libc/utils/MPFRWrapper/MPFRUtils.h
@@ -55,6 +55,7 @@ enum class Operation : int {
ModPIOver4,
Round,
RoundEven,
+ Rsqrt,
Sin,
Sinpi,
Sinh,
>From 10f24106d6947ce0c8408c527e81d9201589be10 Mon Sep 17 00:00:00 2001
From: amemov <shepelev777 at gmail.com>
Date: Sun, 27 Apr 2025 20:07:27 +0000
Subject: [PATCH 2/2] Clang-formated the files
---
libc/src/math/generic/rsqrtf16.cpp | 9 +++++----
libc/test/src/math/smoke/rsqrtf16_test.cpp | 3 +--
2 files changed, 6 insertions(+), 6 deletions(-)
diff --git a/libc/src/math/generic/rsqrtf16.cpp b/libc/src/math/generic/rsqrtf16.cpp
index 930c8f6fc7b4e..8e97ec896aa0c 100644
--- a/libc/src/math/generic/rsqrtf16.cpp
+++ b/libc/src/math/generic/rsqrtf16.cpp
@@ -35,7 +35,7 @@ LLVM_LIBC_FUNCTION(float16, rsqrtf16, (float16 x)) {
}
return x;
}
-
+
// |x| = 0
if (LIBC_UNLIKELY(x_abs == 0x0)) {
fputil::raise_except_if_required(FE_DIVBYZERO);
@@ -44,7 +44,7 @@ LLVM_LIBC_FUNCTION(float16, rsqrtf16, (float16 x)) {
}
// -inf <= x < 0
- if (LIBC_UNLIKELY(x_sign == 1)) {
+ if (LIBC_UNLIKELY(x_sign == 1)) {
fputil::raise_except_if_required(FE_INVALID);
fputil::set_errno_if_required(EDOM);
return FPBits::quiet_nan().get_val();
@@ -54,13 +54,14 @@ LLVM_LIBC_FUNCTION(float16, rsqrtf16, (float16 x)) {
if (LIBC_UNLIKELY(xbits.is_inf())) {
return fputil::cast<float16>(0.0f);
}
-
+
// x = 1 => rsqrt(x) = 1
if (LIBC_UNLIKELY(x_u == 0x1)) {
return fputil::cast<float16>(1.0f);
}
- // x is valid, estimate the result - below is temporary solution for just testing
+ // x is valid, estimate the result - below is temporary solution for just
+ // testing
float xf = x;
return fputil::cast<float16>(1.0f / xf);
}
diff --git a/libc/test/src/math/smoke/rsqrtf16_test.cpp b/libc/test/src/math/smoke/rsqrtf16_test.cpp
index eb429d3dbf6c8..d50d8ee0ad2ff 100644
--- a/libc/test/src/math/smoke/rsqrtf16_test.cpp
+++ b/libc/test/src/math/smoke/rsqrtf16_test.cpp
@@ -17,8 +17,7 @@ TEST_F(LlvmLibcRsqrtf16Test, SpecialNumbers) {
EXPECT_FP_EQ(aNaN, LIBC_NAMESPACE::rsqrtf16(aNaN));
EXPECT_MATH_ERRNO(0);
- EXPECT_FP_EQ_WITH_EXCEPTION(aNaN, LIBC_NAMESPACE::rsqrtf16(sNaN),
- FE_INVALID);
+ EXPECT_FP_EQ_WITH_EXCEPTION(aNaN, LIBC_NAMESPACE::rsqrtf16(sNaN), FE_INVALID);
EXPECT_MATH_ERRNO(0);
EXPECT_FP_EQ(aNaN, LIBC_NAMESPACE::rsqrtf16(0.0f));
More information about the libc-commits
mailing list