[libc-commits] [libc] [libc][math] Implement double-precision sinh (PR #199956)
Aayush Shrivastava via libc-commits
libc-commits at lists.llvm.org
Wed May 27 03:58:47 PDT 2026
https://github.com/iamaayushrivastava created https://github.com/llvm/llvm-project/pull/199956
This PR adds an implementation of the double-precision hyperbolic sine function `sinh` to the `math` library.
>From 287e81e0848ffed80f29d49f53b7612386fa1ac4 Mon Sep 17 00:00:00 2001
From: Aayush Shrivastava <iamaayushrivastava at gmail.com>
Date: Wed, 27 May 2026 01:40:32 +0530
Subject: [PATCH] [libc][math] Implement double-precision sinh
---
libc/config/linux/x86_64/entrypoints.txt | 1 +
libc/src/__support/math/CMakeLists.txt | 14 ++++
libc/src/__support/math/sinh.h | 88 +++++++++++++++++++
libc/src/math/generic/CMakeLists.txt | 10 +++
libc/src/math/generic/sinh.cpp | 16 ++++
libc/test/src/math/CMakeLists.txt | 13 +++
libc/test/src/math/sinh_test.cpp | 102 +++++++++++++++++++++++
libc/test/src/math/smoke/CMakeLists.txt | 12 +++
libc/test/src/math/smoke/sinh_test.cpp | 75 +++++++++++++++++
9 files changed, 331 insertions(+)
create mode 100644 libc/src/__support/math/sinh.h
create mode 100644 libc/src/math/generic/sinh.cpp
create mode 100644 libc/test/src/math/sinh_test.cpp
create mode 100644 libc/test/src/math/smoke/sinh_test.cpp
diff --git a/libc/config/linux/x86_64/entrypoints.txt b/libc/config/linux/x86_64/entrypoints.txt
index 5545790fecd85..6c8557e1e57b7 100644
--- a/libc/config/linux/x86_64/entrypoints.txt
+++ b/libc/config/linux/x86_64/entrypoints.txt
@@ -746,6 +746,7 @@ set(TARGET_LIBM_ENTRYPOINTS
libc.src.math.sincos
libc.src.math.sincosf
libc.src.math.sinf
+ libc.src.math.sinh
libc.src.math.sinhf
libc.src.math.sinpif
libc.src.math.sqrt
diff --git a/libc/src/__support/math/CMakeLists.txt b/libc/src/__support/math/CMakeLists.txt
index 9e3ec26cdc881..c62800ff79a19 100644
--- a/libc/src/__support/math/CMakeLists.txt
+++ b/libc/src/__support/math/CMakeLists.txt
@@ -5221,6 +5221,20 @@ add_header_library(
libc.src.__support.macros.properties.types
)
+add_header_library(
+ sinh
+ HDRS
+ sinh.h
+ DEPENDS
+ .exp
+ .expm1
+ libc.src.__support.FPUtil.fenv_impl
+ libc.src.__support.FPUtil.fp_bits
+ libc.src.__support.FPUtil.rounding_mode
+ libc.src.__support.macros.config
+ libc.src.__support.macros.optimization
+)
+
add_header_library(
sinhf
HDRS
diff --git a/libc/src/__support/math/sinh.h b/libc/src/__support/math/sinh.h
new file mode 100644
index 0000000000000..daaf7e3f7c9ac
--- /dev/null
+++ b/libc/src/__support/math/sinh.h
@@ -0,0 +1,88 @@
+//===-- Implementation header for sinh --------------------------*- 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_SINH_H
+#define LLVM_LIBC_SRC___SUPPORT_MATH_SINH_H
+
+#include "exp.h"
+#include "expm1.h"
+#include "src/__support/FPUtil/FEnvImpl.h"
+#include "src/__support/FPUtil/FPBits.h"
+#include "src/__support/FPUtil/rounding_mode.h"
+#include "src/__support/macros/config.h"
+#include "src/__support/macros/optimization.h" // LIBC_UNLIKELY
+
+namespace LIBC_NAMESPACE_DECL {
+
+namespace math {
+
+LIBC_INLINE double sinh(double x) {
+ using FPBits = fputil::FPBits<double>;
+ FPBits xbits(x);
+ bool is_neg = xbits.is_neg();
+ // Work with |x|.
+ xbits.set_sign(Sign::POS);
+ double x_abs = xbits.get_val();
+ uint64_t x_abs_u = xbits.uintval();
+
+ // Handle NaN: sinh(NaN) = 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;
+ }
+
+ // sinh(+/-inf) = +/-inf.
+ if (LIBC_UNLIKELY(xbits.is_inf()))
+ return x;
+
+ // For very small |x| (|x| <= 2^-26), sinh(x) ~ x.
+ // Use FP arithmetic to ensure FTZ/DAZ mode behavior.
+ if (LIBC_UNLIKELY(x_abs_u <= 0x3e50000000000000ULL)) {
+ // sinh(+/-0) = +/-0, preserve sign of zero exactly.
+ if (LIBC_UNLIKELY(x_abs_u == 0))
+ return x;
+ double x2 = x_abs * x_abs;
+ double result = x_abs + x2 * x_abs * (1.0 / 6.0);
+ return is_neg ? -result : result;
+ }
+
+ // For |x| >= 710, overflow.
+ if (LIBC_UNLIKELY(x_abs_u >= 0x4086340000000000ULL)) {
+ int rounding = fputil::quick_get_round();
+ if ((rounding == FE_DOWNWARD && !is_neg) ||
+ (rounding == FE_UPWARD && is_neg) || rounding == FE_TOWARDZERO)
+ return is_neg ? -FPBits::max_normal().get_val()
+ : FPBits::max_normal().get_val();
+ fputil::set_errno_if_required(ERANGE);
+ fputil::raise_except_if_required(FE_OVERFLOW);
+ return x + (is_neg ? -FPBits::inf().get_val() : FPBits::inf().get_val());
+ }
+
+ double result;
+ if (x_abs_u < 0x3ff0000000000000ULL) {
+ // |x| < 1: use expm1 to avoid catastrophic cancellation.
+ // sinh(x) = expm1(x)/2 + expm1(x)/(2*(1+expm1(x)))
+ double em1 = math::expm1(x_abs);
+ result = em1 / 2.0 + em1 / (2.0 * (1.0 + em1));
+ } else {
+ // |x| >= 1: sinh(x) = (exp(x) - exp(-x)) / 2.
+ double ex = math::exp(x_abs);
+ result = (ex - 1.0 / ex) * 0.5;
+ }
+
+ return is_neg ? -result : result;
+}
+
+} // namespace math
+
+} // namespace LIBC_NAMESPACE_DECL
+
+#endif // LLVM_LIBC_SRC___SUPPORT_MATH_SINH_H
diff --git a/libc/src/math/generic/CMakeLists.txt b/libc/src/math/generic/CMakeLists.txt
index 7ccbddba07b8d..99a0a6939360d 100644
--- a/libc/src/math/generic/CMakeLists.txt
+++ b/libc/src/math/generic/CMakeLists.txt
@@ -3651,6 +3651,16 @@ add_entrypoint_object(
libc.src.__support.math.coshf16
)
+add_entrypoint_object(
+ sinh
+ SRCS
+ sinh.cpp
+ HDRS
+ ../sinh.h
+ DEPENDS
+ libc.src.__support.math.sinh
+)
+
add_entrypoint_object(
sinhf
SRCS
diff --git a/libc/src/math/generic/sinh.cpp b/libc/src/math/generic/sinh.cpp
new file mode 100644
index 0000000000000..ff1b9f5234ba8
--- /dev/null
+++ b/libc/src/math/generic/sinh.cpp
@@ -0,0 +1,16 @@
+//===-- Double-precision sinh implementation ------------------------------===//
+//
+// 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/sinh.h"
+#include "src/__support/math/sinh.h"
+
+namespace LIBC_NAMESPACE_DECL {
+
+LLVM_LIBC_FUNCTION(double, sinh, (double x)) { return math::sinh(x); }
+
+} // namespace LIBC_NAMESPACE_DECL
diff --git a/libc/test/src/math/CMakeLists.txt b/libc/test/src/math/CMakeLists.txt
index 4213c11eca515..bd57e7fabb9ec 100644
--- a/libc/test/src/math/CMakeLists.txt
+++ b/libc/test/src/math/CMakeLists.txt
@@ -2326,6 +2326,19 @@ add_fp_unittest(
libc.src.math.coshf16
)
+add_fp_unittest(
+ sinh_test
+ NEED_MPFR
+ SUITE
+ libc-math-unittests
+ SRCS
+ sinh_test.cpp
+ DEPENDS
+ libc.hdr.errno_macros
+ libc.src.math.sinh
+ libc.src.__support.FPUtil.fp_bits
+)
+
add_fp_unittest(
sinhf_test
NEED_MPFR
diff --git a/libc/test/src/math/sinh_test.cpp b/libc/test/src/math/sinh_test.cpp
new file mode 100644
index 0000000000000..aacf238609d0e
--- /dev/null
+++ b/libc/test/src/math/sinh_test.cpp
@@ -0,0 +1,102 @@
+//===-- Unittests for sinh ------------------------------------------------===//
+//
+// 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/sinh.h"
+#include "test/UnitTest/FPMatcher.h"
+#include "test/UnitTest/Test.h"
+#include "utils/MPFRWrapper/MPFRUtils.h"
+
+using LlvmLibcSinhTest = LIBC_NAMESPACE::testing::FPTest<double>;
+
+namespace mpfr = LIBC_NAMESPACE::testing::mpfr;
+
+using LIBC_NAMESPACE::testing::tlog;
+
+TEST_F(LlvmLibcSinhTest, InDoubleRange) {
+ // sinh is odd; test [2^-26, 710) and negate for negatives.
+ // Values >= 710 overflow to +/-inf.
+ constexpr uint64_t COUNT = 123'451;
+ uint64_t START = FPBits(0x1.0p-26).uintval();
+ uint64_t STOP = FPBits(710.0).uintval();
+ uint64_t STEP = (STOP - START) / COUNT;
+
+ auto test = [&](mpfr::RoundingMode rounding_mode) {
+ mpfr::ForceRoundingMode __r(rounding_mode);
+ if (!__r.success)
+ return;
+
+ uint64_t fails = 0;
+ uint64_t count = 0;
+ uint64_t cc = 0;
+ double mx = 0.0, mr = 0.0;
+ double tol = 0.5;
+
+ for (uint64_t i = 0, v = START; i <= COUNT; ++i, v += STEP) {
+ double x = FPBits(v).get_val();
+ if (FPBits(v).is_inf_or_nan())
+ continue;
+
+ // Test positive x.
+ double result = LIBC_NAMESPACE::sinh(x);
+ ++cc;
+ if (!FPBits(result).is_inf_or_nan())
+ ++count;
+ if (!TEST_MPFR_MATCH_ROUNDING_SILENTLY(mpfr::Operation::Sinh, x, result,
+ 0.5, rounding_mode)) {
+ ++fails;
+ while (!TEST_MPFR_MATCH_ROUNDING_SILENTLY(mpfr::Operation::Sinh, x,
+ result, tol, rounding_mode)) {
+ mx = x;
+ mr = result;
+ if (tol > 1000.0)
+ break;
+ tol *= 2.0;
+ }
+ }
+
+ // Test negative x (sinh is odd).
+ double neg_x = -x;
+ double neg_result = LIBC_NAMESPACE::sinh(neg_x);
+ ++cc;
+ if (!FPBits(neg_result).is_inf_or_nan())
+ ++count;
+ if (!TEST_MPFR_MATCH_ROUNDING_SILENTLY(mpfr::Operation::Sinh, neg_x,
+ neg_result, 0.5, rounding_mode)) {
+ ++fails;
+ while (!TEST_MPFR_MATCH_ROUNDING_SILENTLY(mpfr::Operation::Sinh, neg_x,
+ neg_result, tol,
+ rounding_mode)) {
+ mx = neg_x;
+ mr = neg_result;
+ if (tol > 1000.0)
+ break;
+ tol *= 2.0;
+ }
+ }
+ }
+ if (fails) {
+ tlog << " Sinh failed: " << fails << "/" << count << "/" << cc
+ << " tests.\n";
+ tlog << " Max ULPs is at most: " << static_cast<uint64_t>(tol) << ".\n";
+ EXPECT_MPFR_MATCH(mpfr::Operation::Sinh, mx, mr, 0.5, rounding_mode);
+ }
+ };
+
+ tlog << " Test Rounding To Nearest...\n";
+ test(mpfr::RoundingMode::Nearest);
+
+ tlog << " Test Rounding Downward...\n";
+ test(mpfr::RoundingMode::Downward);
+
+ tlog << " Test Rounding Upward...\n";
+ test(mpfr::RoundingMode::Upward);
+
+ tlog << " Test Rounding Toward Zero...\n";
+ test(mpfr::RoundingMode::TowardZero);
+}
diff --git a/libc/test/src/math/smoke/CMakeLists.txt b/libc/test/src/math/smoke/CMakeLists.txt
index 28b85b1a25bbd..0b691027b00c5 100644
--- a/libc/test/src/math/smoke/CMakeLists.txt
+++ b/libc/test/src/math/smoke/CMakeLists.txt
@@ -4648,6 +4648,18 @@ add_fp_unittest(
libc.src.__support.FPUtil.cast
)
+add_fp_unittest(
+ sinh_test
+ SUITE
+ libc-math-smoke-tests
+ SRCS
+ sinh_test.cpp
+ DEPENDS
+ libc.hdr.errno_macros
+ libc.src.math.sinh
+ libc.src.__support.FPUtil.fp_bits
+)
+
add_fp_unittest(
sinhf_test
SUITE
diff --git a/libc/test/src/math/smoke/sinh_test.cpp b/libc/test/src/math/smoke/sinh_test.cpp
new file mode 100644
index 0000000000000..802ca6801323c
--- /dev/null
+++ b/libc/test/src/math/smoke/sinh_test.cpp
@@ -0,0 +1,75 @@
+//===-- Unittests for sinh ------------------------------------------------===//
+//
+// 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 "hdr/errno_macros.h"
+#include "hdr/math_macros.h"
+#include "src/__support/FPUtil/FPBits.h"
+#include "src/math/sinh.h"
+#include "test/UnitTest/FPMatcher.h"
+#include "test/UnitTest/Test.h"
+
+using LlvmLibcSinhTest = LIBC_NAMESPACE::testing::FPTest<double>;
+
+TEST_F(LlvmLibcSinhTest, SpecialNumbers) {
+ EXPECT_FP_EQ_WITH_EXCEPTION(aNaN, LIBC_NAMESPACE::sinh(sNaN), FE_INVALID);
+ EXPECT_MATH_ERRNO(0);
+
+ EXPECT_FP_EQ_ALL_ROUNDING(aNaN, LIBC_NAMESPACE::sinh(aNaN));
+ EXPECT_MATH_ERRNO(0);
+
+ EXPECT_FP_EQ_ALL_ROUNDING(inf, LIBC_NAMESPACE::sinh(inf));
+ EXPECT_MATH_ERRNO(0);
+
+ EXPECT_FP_EQ_ALL_ROUNDING(neg_inf, LIBC_NAMESPACE::sinh(neg_inf));
+ EXPECT_MATH_ERRNO(0);
+
+ EXPECT_FP_EQ_ALL_ROUNDING(0.0, LIBC_NAMESPACE::sinh(0.0));
+ EXPECT_MATH_ERRNO(0);
+
+ EXPECT_FP_EQ_ALL_ROUNDING(neg_zero, LIBC_NAMESPACE::sinh(neg_zero));
+ EXPECT_MATH_ERRNO(0);
+}
+
+TEST_F(LlvmLibcSinhTest, Overflow) {
+ EXPECT_FP_EQ_WITH_EXCEPTION(inf, LIBC_NAMESPACE::sinh(0x1.0p10), FE_OVERFLOW);
+ EXPECT_MATH_ERRNO(ERANGE);
+
+ EXPECT_FP_EQ_WITH_EXCEPTION(neg_inf, LIBC_NAMESPACE::sinh(-0x1.0p10),
+ FE_OVERFLOW);
+ EXPECT_MATH_ERRNO(ERANGE);
+}
+
+#ifdef LIBC_TEST_FTZ_DAZ
+
+using namespace LIBC_NAMESPACE::testing;
+
+TEST_F(LlvmLibcSinhTest, FTZMode) {
+ ModifyMXCSR mxcsr(FTZ);
+
+ const double neg_min_denormal = FPBits::min_subnormal(Sign::NEG).get_val();
+ EXPECT_FP_EQ(0.0, LIBC_NAMESPACE::sinh(min_denormal));
+ EXPECT_FP_EQ(neg_zero, LIBC_NAMESPACE::sinh(neg_min_denormal));
+}
+
+TEST_F(LlvmLibcSinhTest, DAZMode) {
+ ModifyMXCSR mxcsr(DAZ);
+
+ const double neg_min_denormal = FPBits::min_subnormal(Sign::NEG).get_val();
+ EXPECT_FP_EQ(0.0, LIBC_NAMESPACE::sinh(min_denormal));
+ EXPECT_FP_EQ(neg_zero, LIBC_NAMESPACE::sinh(neg_min_denormal));
+}
+
+TEST_F(LlvmLibcSinhTest, FTZDAZMode) {
+ ModifyMXCSR mxcsr(FTZ | DAZ);
+
+ const double neg_min_denormal = FPBits::min_subnormal(Sign::NEG).get_val();
+ EXPECT_FP_EQ(0.0, LIBC_NAMESPACE::sinh(min_denormal));
+ EXPECT_FP_EQ(neg_zero, LIBC_NAMESPACE::sinh(neg_min_denormal));
+}
+
+#endif
More information about the libc-commits
mailing list