[libc][math][c23] add sinpi function

[ghost]

This is the implementation of sinpi.

[github-actions]

Thank you for submitting a Pull Request (PR) to the LLVM Project!

This PR will be automatically labeled and the relevant teams will be notified.

If you wish to, you can add reviewers by using the "Reviewers" section on this page.

If this is not working for you, it is probably because you do not have write permissions for the repository. In which case you can instead tag reviewers by name in a comment by using @ followed by their GitHub username.

If you have received no comments on your PR for a week, you can request a review by "ping"ing the PR by adding a comment “Ping”. The common courtesy "ping" rate is once a week. Please remember that you are asking for valuable time from other developers.

If you have further questions, they may be answered by the LLVM GitHub User Guide.

You can also ask questions in a comment on this PR, on the LLVM Discord or on the forums.

[llvmbot]

@llvm/pr-subscribers-libc

Author: None (lara254)

Changes

This is the implementation of sinpi.

---

Full diff: https://github.com/llvm/llvm-project/pull/129379.diff

1 Files Affected:

• (added) libc/src/math/generic/sinpi.cpp (+63)

diff --git a/libc/src/math/generic/sinpi.cpp b/libc/src/math/generic/sinpi.cpp
new file mode 100644
index 0000000000000..69116ced9fbc7
--- /dev/null
+++ b/libc/src/math/generic/sinpi.cpp
@@ -0,0 +1,63 @@
+//===-- double-precision sinpi 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/sinpi.h"
+#include "sincos_eval.h"
+#include "src/__support/FPUtil/FEnvImpl.h"
+#include "src/__support/FPUtil/FPBits.h"
+#include "src/__support/FPUtil/PolyEval.h"
+#include "src/__support/FPUtil/multiply_add.h"
+#include "src/__support/common.h"
+#include "src/__support/macros/config.h"
+#include "src/__support/macros/optimization.h" // LIBC_UNLIKELY
+
+namespace LIBC_NAMESPACE_DECL {
+
+LLVM_LIBC_FUNCTION(double, sinpi, (double x)) {
+  using FPBits = typename fputil::FPBits<double>;
+  FPBits xbits(x);
+
+  uint64_t x_u = xbits.uintval();
+  uint64_t x_abs = x_u & 0x7fff'ffffU;
+  long double xd = static_cast<long double>(x);
+
+  if (LIBC_UNLIKELY(x_abs >= 0x7c00)) {
+    if (x_abs == 0x7c00) {
+      fputil::set_errno_if_required(EDOM);
+      fputil::raise_except_if_required(FE_INVALID);
+    }
+    return x + FPBits::quiet_nan().get_val();
+  } else {
+    return FPBits::zero(xbits.sign()).get_val();
+  }
+
+  if (LIBC_UNLIKELY(x_abs <= 0x3d80'0000U)) {
+    if (LIBC_UNLIKELY(x_abs < 0x33CD'01D7U)) {
+      if (LIBC_UNLIKELY(x_abs == 0U)) {
+	return x;
+    }
+      long double xdpi = xd * 0x1.921fb5444d18p1;
+       return  static_cast<double>(xdpi);
+    }
+  long double xsq = xd * xd;
+
+  long double result = fputil::polyeval(xsq,0x1.921fb54442d183f07b2385653d8p1, -0x1.4abbce625bd95cdc955aeed9abcp2, 0x1.466bc6769ddfdb085486c0ff3ep1, -0x1.32d2c4a48bfd71fa9cdf60a0e4p-1,  0x1.502cbd2c72e3168ff209bc7656cp-4);
+
+  return static_cast<double>(xd * result);
+  }
+
+  long double sin_k, cos_k, sin_y, cosm1_y;
+  sincospi_eval(xd, sin_k, cos_k, sin_y, cosm1_y);
+
+  if (LIBC_UNLIKELY(sin_y == 0 && sin_k == 0))
+    return FPBits::zero(xbits.sign()).get_val();
+
+  return static_cast<double>(fputil::multiply_add(
+						  sin_y, cos_k, fputil::multiply_add(cosm1_y, sin_k, sin_k)));
+}
+}

[github-actions]

⚠️ C/C++ code formatter, clang-format found issues in your code. ⚠️

You can test this locally with the following command:

git-clang-format --diff HEAD~1 HEAD --extensions h,cpp -- libc/src/math/generic/sinpi.cpp libc/src/math/sinpi.h libc/test/src/math/sinpi_test.cpp libc/test/src/math/smoke/sinpi_test.cpp

View the diff from clang-format here.

diff --git a/libc/src/math/generic/sinpi.cpp b/libc/src/math/generic/sinpi.cpp
index f51249da0..fdc16548c 100644
--- a/libc/src/math/generic/sinpi.cpp
+++ b/libc/src/math/generic/sinpi.cpp
@@ -3,22 +3,22 @@
 #include "src/__support/FPUtil/BasicOperations.h"
 #include "src/__support/FPUtil/FEnvImpl.h"
 #include "src/__support/FPUtil/FPBits.h"
-#include "src/__support/FPUtil/multiply_add.h"
-#include "src/__support/common.h"
-#include "src/__support/macros/config.h"
 #include "src/__support/FPUtil/double_double.h"
 #include "src/__support/FPUtil/generic/mul.h"
+#include "src/__support/FPUtil/multiply_add.h"
 #include "src/__support/FPUtil/nearest_integer.h"
+#include "src/__support/common.h"
+#include "src/__support/macros/config.h"
 #include "src/math/pow.h"
-//#include "range_reduction_double_nofma.h"
-//#include "src/__support/FPUtil/multiply_add.h"
-//#include "src/math/generic/range_reduction_double_common.h"
+// #include "range_reduction_double_nofma.h"
+// #include "src/__support/FPUtil/multiply_add.h"
+// #include "src/math/generic/range_reduction_double_common.h"
 #include <iostream>
 
 namespace LIBC_NAMESPACE_DECL {
 
 LLVM_LIBC_FUNCTION(double, sinpi, (double x)) {
-  // Given x * pi = y - (k * (pi/128)) 
+  // Given x * pi = y - (k * (pi/128))
   // find y and k such that
   // y = x * pi - (k * pi/128) = x * pi - kpi/128
   // k = round(x, 128)
@@ -29,15 +29,17 @@ LLVM_LIBC_FUNCTION(double, sinpi, (double x)) {
   FPBits xbits(x);
 
   double k = fputil::nearest_integer(x * 128);
-  int  k_int = static_cast<int>(k);
-  
+  int k_int = static_cast<int>(k);
+
   std::cout << "k" << k << std::endl;
-  
-  double yk = x - k/128;
 
-  DoubleDouble yy = fputil::exact_mult(yk, 3.141592653589793115997963468544185161590576171875);
-  yy.lo = fputil::multiply_add(yk, 1.2246467991473532071737640294583966046256921246776e-16, yy.lo);
-  
+  double yk = x - k / 128;
+
+  DoubleDouble yy = fputil::exact_mult(
+      yk, 3.141592653589793115997963468544185161590576171875);
+  yy.lo = fputil::multiply_add(
+      yk, 1.2246467991473532071737640294583966046256921246776e-16, yy.lo);
+
   uint64_t abs_u = xbits.uintval();
 
   uint64_t x_abs = abs_u & 0xFFFFFFFFFFFFFFFF;
@@ -60,36 +62,37 @@ LLVM_LIBC_FUNCTION(double, sinpi, (double x)) {
   [[maybe_unused]] double err = generic::sincos_eval(yy, sin_y, cos_y);
   DoubleDouble sin_k = SIN_K_PI_OVER_128[k_int & 255];
   DoubleDouble cos_k = SIN_K_PI_OVER_128[(k_int + 64) & 255];
-  
+
   std::cout << "sin_k: " << sin_k.hi << std::endl;
   std::cout << "sin_klo: " << sin_k.lo << std::endl;
   std::cout << "sin_y: " << sin_y.hi << std::endl;
   std::cout << "cos_y: " << cos_y.hi << std::endl;
   std::cout << "sin_y.lo: " << sin_y.lo << std::endl;
   std::cout << "cos_y.o: " << cos_y.lo << std::endl;
-  
+
   double cosm1_y = cos_y.hi - 1.0;
   DoubleDouble sin_y_cos_k = fputil::quick_mult(sin_y, cos_k);
-  
+
   std::cout << "cosm1" << cosm1_y << std::endl;
   DoubleDouble cosm1_yy;
   cosm1_yy.hi = cosm1_y;
   cosm1_yy.lo = 0.0;
-  
+
   DoubleDouble cos_y_sin_k = fputil::quick_mult(cos_y, sin_k);
   DoubleDouble rr = fputil::exact_add<false>(sin_y_cos_k.hi, cos_y_sin_k.hi);
-  
+
   std::cout << "r.hi:" << rr.hi << std::endl;
   std::cout << "r.lo" << rr.lo << std::endl;
-  
+
   rr.lo += sin_y_cos_k.lo + cos_y_sin_k.lo;
- 
-  std::cout << "rrlo2: " << rr.lo << std::endl;  
+
+  std::cout << "rrlo2: " << rr.lo << std::endl;
   std::cout << "cos_y_sin_k:" << cos_y_sin_k.hi << std::endl;
   std::cout << "siny*cosk.lo:" << sin_y_cos_k.lo << std::endl;
-  std::cout << "rrhi + rrlo + sink.hi " << rr.hi + rr.lo + sin_k.hi + sin_k.lo << std::endl;
+  std::cout << "rrhi + rrlo + sink.hi " << rr.hi + rr.lo + sin_k.hi + sin_k.lo
+            << std::endl;
   std::cout << "rrhi + rrlo " << rr.hi + rr.lo << std::endl;
-  
+
   return rr.hi + rr.lo;
- }
+}
 } // namespace LIBC_NAMESPACE_DECL
diff --git a/libc/test/src/math/sinpi_test.cpp b/libc/test/src/math/sinpi_test.cpp
index fda46173f..54dd4fa5f 100644
--- a/libc/test/src/math/sinpi_test.cpp
+++ b/libc/test/src/math/sinpi_test.cpp
@@ -48,8 +48,7 @@ TEST_F(LlvmLibcSinpiTest, InDoubleRange) {
         continue;
 
       ++count;
-     
-      
+
       if (!TEST_MPFR_MATCH_ROUNDING_SILENTLY(mpfr::Operation::Sinpi, x, result,
                                              2.0, rounding_mode)) {
         ++fails;
@@ -71,7 +70,6 @@ TEST_F(LlvmLibcSinpiTest, InDoubleRange) {
       tlog << "   Max ULPs is at most: " << static_cast<uint64_t>(tol) << ".\n";
       EXPECT_MPFR_MATCH(mpfr::Operation::Sinpi, mx, mr, 2.0, rounding_mode);
     }
-      
   };
   tlog << " Test Rounding To Nearest...\n";
   test(mpfr::RoundingMode::Nearest);
diff --git a/libc/test/src/math/smoke/sinpi_test.cpp b/libc/test/src/math/smoke/sinpi_test.cpp
index 0c6071ee6..1f1bb1466 100644
--- a/libc/test/src/math/smoke/sinpi_test.cpp
+++ b/libc/test/src/math/smoke/sinpi_test.cpp
@@ -35,7 +35,7 @@ TEST_F(LlvmLibcSinpiTest, SpecialNumbers) {
 
 TEST_F(LlvmLibcSinpiTest, Integers) {
   EXPECT_FP_EQ(neg_zero, LIBC_NAMESPACE::sinpi(-0x1.0000000000003p52));
-  
+
   ASSERT_FP_EQ(-1.0, LIBC_NAMESPACE::sinpi(4499003037990983.5));
 
   EXPECT_FP_EQ(neg_zero, LIBC_NAMESPACE::sinpi(-0x1.0000000000005p52));
@@ -53,7 +53,11 @@ TEST_F(LlvmLibcSinpiTest, Integers) {
 */
 
 TEST_F(LlvmLibcSinpiTest, Integers) {
-  //ASSERT_FP_EQ(zero, LIBC_NAMESPACE::sinpi(4503563146482784.00000000000000000000000000000000000000000000000000));
+  // ASSERT_FP_EQ(zero,
+  // LIBC_NAMESPACE::sinpi(4503563146482784.00000000000000000000000000000000000000000000000000));
 
-  ASSERT_FP_EQ(-1.0, LIBC_NAMESPACE::sinpi(4499003037990983.50000000000000000000000000000000000000000000000000));
+  ASSERT_FP_EQ(
+      -1.0,
+      LIBC_NAMESPACE::sinpi(
+          4499003037990983.50000000000000000000000000000000000000000000000000));
 }

[lntue]

0x1.0p-52 is the one that you need instead of pow, and remove dependency / include of pow.

Also are you doing range reduction mod 2^-52?