[libc-commits] [libc] [llvm] [libc][math] Add float-only implementation of tanf (PR #233046)

Prabhu Rajasekaran via libc-commits libc-commits at lists.llvm.org
Sun Oct 18 12:34:15 PDT 2026


https://github.com/Prabhuk updated https://github.com/llvm/llvm-project/pull/233046

>From e865a2d46f41e4ae74956994e38410fefe56d89d Mon Sep 17 00:00:00 2001
From: prabhukr <prabhukr at google.com>
Date: Sat, 17 Oct 2026 07:44:44 +0000
Subject: [PATCH 1/3] [libc][math] Add float-only implementation of tanf

Reorganize tanf in the same way as sinf and cosf (#226529):
src/__support/math/tanf.h now selects the implementation, and
tanf_<type>_eval.h implements tanf with <type> as the intermediate
computational type. The existing implementation moves unchanged to
tanf_double_eval.h.

The new float-only implementation in tanf_float_eval.h is selected
under the same conditions as sinf/cosf float_eval:
LIBC_MATH_SKIP_ACCURATE_PASS, LIBC_MATH_SMALL_TABLES or
LIBC_MATH_INTERMEDIATE_COMP_IN_FLOAT, and FMA. It reuses the float range
reduction and the SIN_K_PI_OVER_8 table from sincosf_float_eval.h and
computes:

  tan(x) = (sin_k + cos_k * tan(y)) / (cos_k - sin_k * tan(y)),

with x = k * pi/8 + y, |y| <= pi/16, and tan(y) approximated by a degree-9
odd Taylor polynomial (relative truncation error < 2^-30). Cancellation
can only happen when k = 0 or 4 (mod 8), where sin_k or cos_k is exactly
zero.

With this change, tanf no longer pulls in the double-precision
large_range_reduction and the 512-byte SIN_K_PI_OVER_32 table when the
float implementations are selected, so it shares the reduction tables
with sinf/cosf. For x86-64 with FMA and
LIBC_MATH_SKIP_ACCURATE_PASS | LIBC_MATH_SMALL_TABLES |
LIBC_MATH_INTERMEDIATE_COMP_IN_FLOAT, the tanf object goes from 1391 to
1029 bytes standalone (288 bytes of which are tables shared with
sinf/cosf).

Accuracy: the maximum error of the float implementation, checked
exhaustively in round-to-nearest mode against MPFR, is 3.37 ULP, which is
within the 3.5 ULP bound used for sinf/cosf float_eval. The unit and
exhaustive tests now exercise both implementations independently of the
configured selector.

Assisted-by: Gemini
---
 libc/src/__support/math/CMakeLists.txt        |  30 ++-
 libc/src/__support/math/tanf.h                | 167 +++--------------
 libc/src/__support/math/tanf_double_eval.h    | 172 ++++++++++++++++++
 libc/src/__support/math/tanf_float_eval.h     | 106 +++++++++++
 libc/test/src/math/CMakeLists.txt             |   3 +
 libc/test/src/math/exhaustive/CMakeLists.txt  |   1 +
 libc/test/src/math/exhaustive/tanf_test.cpp   |  19 ++
 libc/test/src/math/tanf_test.cpp              |  60 +++++-
 .../llvm-project-overlay/libc/BUILD.bazel     |   7 +-
 9 files changed, 407 insertions(+), 158 deletions(-)
 create mode 100644 libc/src/__support/math/tanf_double_eval.h
 create mode 100644 libc/src/__support/math/tanf_float_eval.h

diff --git a/libc/src/__support/math/CMakeLists.txt b/libc/src/__support/math/CMakeLists.txt
index baac96e76d5435..e339d3422c9ae0 100644
--- a/libc/src/__support/math/CMakeLists.txt
+++ b/libc/src/__support/math/CMakeLists.txt
@@ -5821,9 +5821,9 @@ add_header_library(
 )
 
 add_header_library(
-  tanf
+  tanf_double_eval
   HDRS
-    tanf.h
+    tanf_double_eval.h
   DEPENDS
     .range_reduction
     .sincosf_utils
@@ -5836,6 +5836,32 @@ add_header_library(
     libc.src.__support.macros.optimization
 )
 
+add_header_library(
+  tanf_float_eval
+  HDRS
+    tanf_float_eval.h
+  DEPENDS
+    .sincosf_utils
+    libc.src.__support.FPUtil.fenv_impl
+    libc.src.__support.FPUtil.fp_bits
+    libc.src.__support.FPUtil.multiply_add
+    libc.src.__support.common
+    libc.src.__support.macros.config
+    libc.src.__support.macros.optimization
+)
+
+add_header_library(
+  tanf
+  HDRS
+    tanf.h
+  DEPENDS
+    .tanf_double_eval
+    .tanf_float_eval
+    libc.src.__support.macros.config
+    libc.src.__support.macros.optimization
+    libc.src.__support.macros.properties.cpu_features
+)
+
 add_header_library(
   tanf16
   HDRS
diff --git a/libc/src/__support/math/tanf.h b/libc/src/__support/math/tanf.h
index 3a1d8e8bd53b62..9574156801189a 100644
--- a/libc/src/__support/math/tanf.h
+++ b/libc/src/__support/math/tanf.h
@@ -1,165 +1,42 @@
-//===-- Single-precision tan 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
 //
 //===----------------------------------------------------------------------===//
+///
+/// \file
+/// Implementation selector for tanf.
+///
+//===----------------------------------------------------------------------===//
 
 #ifndef LLVM_LIBC_SRC___SUPPORT_MATH_TANF_H
 #define LLVM_LIBC_SRC___SUPPORT_MATH_TANF_H
 
-#include "sincosf_utils.h"
-#include "src/__support/FPUtil/FEnvImpl.h"
-#include "src/__support/FPUtil/FPBits.h"
-#include "src/__support/FPUtil/PolyEval.h"
-#include "src/__support/FPUtil/except_value_utils.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/macros/properties/cpu_features.h" // LIBC_TARGET_CPU_HAS_FMA
+#include "src/__support/macros/optimization.h"
+#include "src/__support/macros/properties/cpu_features.h"
+
+#if defined(LIBC_MATH_HAS_SKIP_ACCURATE_PASS) &&                               \
+    (defined(LIBC_MATH_HAS_SMALL_TABLES) ||                                    \
+     defined(LIBC_MATH_HAS_INTERMEDIATE_COMP_IN_FLOAT)) &&                     \
+    defined(LIBC_TARGET_CPU_HAS_FMA_FLOAT)
+#include "tanf_float_eval.h"
+#define LIBC_MATH_TANF_IMPL float_eval
+#else
+#include "tanf_double_eval.h"
+#define LIBC_MATH_TANF_IMPL double_eval
+#endif
 
 namespace LIBC_NAMESPACE_DECL {
-
 namespace math {
 
-namespace tanf_internal {
-
-#ifndef LIBC_MATH_HAS_SKIP_ACCURATE_PASS
-// Exceptional cases for tanf.
-LIBC_INLINE_VAR constexpr size_t N_EXCEPTS = 6;
-
-LIBC_INLINE_VAR constexpr fputil::ExceptValues<float, N_EXCEPTS> TANF_EXCEPTS{{
-    // (inputs, RZ output, RU offset, RD offset, RN offset)
-    // x = 0x1.ada6aap27, tan(x) = 0x1.e80304p-3 (RZ)
-    {0x4d56d355, 0x3e740182, 1, 0, 0},
-    // x = 0x1.862064p33, tan(x) = -0x1.8dee56p-3 (RZ)
-    {0x50431032, 0xbe46f72b, 0, 1, 1},
-    // x = 0x1.af61dap48, tan(x) = 0x1.60d1c6p-2 (RZ)
-    {0x57d7b0ed, 0x3eb068e3, 1, 0, 1},
-    // x = 0x1.0088bcp52, tan(x) = 0x1.ca1edp0 (RZ)
-    {0x5980445e, 0x3fe50f68, 1, 0, 0},
-    // x = 0x1.f90dfcp72, tan(x) = 0x1.597f9cp-1 (RZ)
-    {0x63fc86fe, 0x3f2cbfce, 1, 0, 0},
-    // x = 0x1.a6ce12p86, tan(x) = -0x1.c5612ep-1 (RZ)
-    {0x6ad36709, 0xbf62b097, 0, 1, 0},
-}};
-#endif // !LIBC_MATH_HAS_SKIP_ACCURATE_PASS
-
-} // namespace tanf_internal
-
-LIBC_INLINE float tanf(float x) {
-  using namespace sincosf_utils_internal;
-  using namespace tanf_internal;
-  using FPBits = typename fputil::FPBits<float>;
-  FPBits xbits(x);
-  uint32_t x_abs = xbits.uintval() & 0x7fff'ffffU;
-
-  // |x| < pi/32
-  if (LIBC_UNLIKELY(x_abs <= 0x3dc9'0fdbU)) {
-    double xd = static_cast<double>(x);
-
-    // |x| < 0x1.0p-12f
-    if (LIBC_UNLIKELY(x_abs < 0x3980'0000U)) {
-      if (LIBC_UNLIKELY(x_abs == 0U)) {
-        // For signed zeros.
-        return x;
-      }
-      // When |x| < 2^-12, the relative error of the approximation tan(x) ~ x
-      // is:
-      //   |tan(x) - x| / |tan(x)| < |x^3| / (3|x|)
-      //                           = x^2 / 3
-      //                           < 2^-25
-      //                           < epsilon(1)/2.
-      // So the correctly rounded values of tan(x) are:
-      //   = x + sign(x)*eps(x) if rounding mode = FE_UPWARD and x is positive,
-      //                        or (rounding mode = FE_DOWNWARD and x is
-      //                        negative),
-      //   = x otherwise.
-      // To simplify the rounding decision and make it more efficient, we use
-      //   fma(x, 2^-25, x) instead.
-      // Note: to use the formula x + 2^-25*x to decide the correct rounding, we
-      // do need fma(x, 2^-25, x) to prevent underflow caused by 2^-25*x when
-      // |x| < 2^-125. For targets without FMA instructions, we simply use
-      // double for intermediate results as it is more efficient than using an
-      // emulated version of FMA.
-#if defined(LIBC_TARGET_CPU_HAS_FMA_FLOAT)
-      return fputil::multiply_add(x, 0x1.0p-25f, x);
-#else
-      return static_cast<float>(fputil::multiply_add(xd, 0x1.0p-25, xd));
-#endif // LIBC_TARGET_CPU_HAS_FMA_FLOAT
-    }
-
-    // |x| < pi/32
-    double xsq = xd * xd;
-
-    // Degree-9 minimax odd polynomial of tan(x) generated by Sollya with:
-    // > P = fpminimax(tan(x)/x, [|0, 2, 4, 6, 8|], [|1, D...|], [0, pi/32]);
-    double result =
-        fputil::polyeval(xsq, 1.0, 0x1.555555553d022p-2, 0x1.111111ce442c1p-3,
-                         0x1.ba180a6bbdecdp-5, 0x1.69c0a88a0b71fp-6);
-    return static_cast<float>(xd * result);
-  }
-
-#ifndef LIBC_MATH_HAS_SKIP_ACCURATE_PASS
-  bool x_sign = xbits.uintval() >> 31;
-  // Check for exceptional values
-  if (LIBC_UNLIKELY(x_abs == 0x3f8a1f62U)) {
-    // |x| = 0x1.143ec4p0
-    float sign = x_sign ? -1.0f : 1.0f;
-
-    // volatile is used to prevent compiler (gcc) from optimizing the
-    // computation, making the results incorrect in different rounding modes.
-    volatile float tmp = 0x1.ddf9f4p0f;
-    tmp = fputil::multiply_add(sign, tmp, sign * 0x1.1p-24f);
-
-    return tmp;
-  }
-#endif // !LIBC_MATH_HAS_SKIP_ACCURATE_PASS
-
-  // |x| > 0x1.ada6a8p+27f
-  if (LIBC_UNLIKELY(x_abs > 0x4d56'd354U)) {
-    // Inf or NaN
-    if (LIBC_UNLIKELY(x_abs >= 0x7f80'0000U)) {
-      if (xbits.is_signaling_nan()) {
-        fputil::raise_except_if_required(FE_INVALID);
-        return FPBits::quiet_nan().get_val();
-      }
-
-      if (x_abs == 0x7f80'0000U) {
-        fputil::set_errno_if_required(EDOM);
-        fputil::raise_except_if_required(FE_INVALID);
-      }
-      return x + FPBits::quiet_nan().get_val();
-    }
-#ifndef LIBC_MATH_HAS_SKIP_ACCURATE_PASS
-    // Other large exceptional values
-    if (auto r = TANF_EXCEPTS.lookup_odd(x_abs, x_sign);
-        LIBC_UNLIKELY(r.has_value()))
-      return r.value();
-#endif // !LIBC_MATH_HAS_SKIP_ACCURATE_PASS
-  }
-
-  // For |x| >= pi/32, we use the definition of tan(x) function:
-  //   tan(x) = sin(x) / cos(x)
-  // The we follow the same computations of sin(x) and cos(x) as sinf, cosf,
-  // and sincosf.
-
-  double xd = static_cast<double>(x);
-  double sin_k, cos_k, sin_y, cosm1_y;
-
-  sincosf_eval(xd, x_abs, sin_k, cos_k, sin_y, cosm1_y);
-  // tan(x) = sin(x) / cos(x)
-  //        = (sin_y * cos_k + cos_y * sin_k) / (cos_y * cos_k - sin_y * sin_k)
-  using fputil::multiply_add;
-  return static_cast<float>(
-      multiply_add(sin_y, cos_k, multiply_add(cosm1_y, sin_k, sin_k)) /
-      multiply_add(sin_y, -sin_k, multiply_add(cosm1_y, cos_k, cos_k)));
-}
+using LIBC_MATH_TANF_IMPL::tanf;
 
 } // namespace math
-
 } // namespace LIBC_NAMESPACE_DECL
 
+#undef LIBC_MATH_TANF_IMPL
+
 #endif // LLVM_LIBC_SRC___SUPPORT_MATH_TANF_H
diff --git a/libc/src/__support/math/tanf_double_eval.h b/libc/src/__support/math/tanf_double_eval.h
new file mode 100644
index 00000000000000..42a8a7ee9a62a1
--- /dev/null
+++ b/libc/src/__support/math/tanf_double_eval.h
@@ -0,0 +1,172 @@
+//===----------------------------------------------------------------------===//
+//
+// 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
+//
+//===----------------------------------------------------------------------===//
+///
+/// \file
+/// Double-precision implementation of tanf.
+///
+//===----------------------------------------------------------------------===//
+
+#ifndef LLVM_LIBC_SRC___SUPPORT_MATH_TANF_DOUBLE_EVAL_H
+#define LLVM_LIBC_SRC___SUPPORT_MATH_TANF_DOUBLE_EVAL_H
+
+#include "sincosf_utils.h"
+#include "src/__support/FPUtil/FEnvImpl.h"
+#include "src/__support/FPUtil/FPBits.h"
+#include "src/__support/FPUtil/PolyEval.h"
+#include "src/__support/FPUtil/except_value_utils.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/macros/properties/cpu_features.h" // LIBC_TARGET_CPU_HAS_FMA
+
+namespace LIBC_NAMESPACE_DECL {
+
+namespace math {
+namespace double_eval {
+
+namespace tanf_internal {
+
+#ifndef LIBC_MATH_HAS_SKIP_ACCURATE_PASS
+// Exceptional cases for tanf.
+LIBC_INLINE_VAR constexpr size_t N_EXCEPTS = 6;
+
+LIBC_INLINE_VAR constexpr fputil::ExceptValues<float, N_EXCEPTS> TANF_EXCEPTS{{
+    // (inputs, RZ output, RU offset, RD offset, RN offset)
+    // x = 0x1.ada6aap27, tan(x) = 0x1.e80304p-3 (RZ)
+    {0x4d56d355, 0x3e740182, 1, 0, 0},
+    // x = 0x1.862064p33, tan(x) = -0x1.8dee56p-3 (RZ)
+    {0x50431032, 0xbe46f72b, 0, 1, 1},
+    // x = 0x1.af61dap48, tan(x) = 0x1.60d1c6p-2 (RZ)
+    {0x57d7b0ed, 0x3eb068e3, 1, 0, 1},
+    // x = 0x1.0088bcp52, tan(x) = 0x1.ca1edp0 (RZ)
+    {0x5980445e, 0x3fe50f68, 1, 0, 0},
+    // x = 0x1.f90dfcp72, tan(x) = 0x1.597f9cp-1 (RZ)
+    {0x63fc86fe, 0x3f2cbfce, 1, 0, 0},
+    // x = 0x1.a6ce12p86, tan(x) = -0x1.c5612ep-1 (RZ)
+    {0x6ad36709, 0xbf62b097, 0, 1, 0},
+}};
+#endif // !LIBC_MATH_HAS_SKIP_ACCURATE_PASS
+
+} // namespace tanf_internal
+
+LIBC_INLINE float tanf(float x) {
+  using namespace sincosf_utils_internal;
+  using namespace tanf_internal;
+  using FPBits = typename fputil::FPBits<float>;
+  FPBits xbits(x);
+  uint32_t x_abs = xbits.uintval() & 0x7fff'ffffU;
+
+  // |x| < pi/32
+  if (LIBC_UNLIKELY(x_abs <= 0x3dc9'0fdbU)) {
+    double xd = static_cast<double>(x);
+
+    // |x| < 0x1.0p-12f
+    if (LIBC_UNLIKELY(x_abs < 0x3980'0000U)) {
+      if (LIBC_UNLIKELY(x_abs == 0U)) {
+        // For signed zeros.
+        return x;
+      }
+      // When |x| < 2^-12, the relative error of the approximation tan(x) ~ x
+      // is:
+      //   |tan(x) - x| / |tan(x)| < |x^3| / (3|x|)
+      //                           = x^2 / 3
+      //                           < 2^-25
+      //                           < epsilon(1)/2.
+      // So the correctly rounded values of tan(x) are:
+      //   = x + sign(x)*eps(x) if rounding mode = FE_UPWARD and x is positive,
+      //                        or (rounding mode = FE_DOWNWARD and x is
+      //                        negative),
+      //   = x otherwise.
+      // To simplify the rounding decision and make it more efficient, we use
+      //   fma(x, 2^-25, x) instead.
+      // Note: to use the formula x + 2^-25*x to decide the correct rounding, we
+      // do need fma(x, 2^-25, x) to prevent underflow caused by 2^-25*x when
+      // |x| < 2^-125. For targets without FMA instructions, we simply use
+      // double for intermediate results as it is more efficient than using an
+      // emulated version of FMA.
+#if defined(LIBC_TARGET_CPU_HAS_FMA_FLOAT)
+      return fputil::multiply_add(x, 0x1.0p-25f, x);
+#else
+      return static_cast<float>(fputil::multiply_add(xd, 0x1.0p-25, xd));
+#endif // LIBC_TARGET_CPU_HAS_FMA_FLOAT
+    }
+
+    // |x| < pi/32
+    double xsq = xd * xd;
+
+    // Degree-9 minimax odd polynomial of tan(x) generated by Sollya with:
+    // > P = fpminimax(tan(x)/x, [|0, 2, 4, 6, 8|], [|1, D...|], [0, pi/32]);
+    double result =
+        fputil::polyeval(xsq, 1.0, 0x1.555555553d022p-2, 0x1.111111ce442c1p-3,
+                         0x1.ba180a6bbdecdp-5, 0x1.69c0a88a0b71fp-6);
+    return static_cast<float>(xd * result);
+  }
+
+#ifndef LIBC_MATH_HAS_SKIP_ACCURATE_PASS
+  bool x_sign = xbits.uintval() >> 31;
+  // Check for exceptional values
+  if (LIBC_UNLIKELY(x_abs == 0x3f8a1f62U)) {
+    // |x| = 0x1.143ec4p0
+    float sign = x_sign ? -1.0f : 1.0f;
+
+    // volatile is used to prevent compiler (gcc) from optimizing the
+    // computation, making the results incorrect in different rounding modes.
+    volatile float tmp = 0x1.ddf9f4p0f;
+    tmp = fputil::multiply_add(sign, tmp, sign * 0x1.1p-24f);
+
+    return tmp;
+  }
+#endif // !LIBC_MATH_HAS_SKIP_ACCURATE_PASS
+
+  // |x| > 0x1.ada6a8p+27f
+  if (LIBC_UNLIKELY(x_abs > 0x4d56'd354U)) {
+    // Inf or NaN
+    if (LIBC_UNLIKELY(x_abs >= 0x7f80'0000U)) {
+      if (xbits.is_signaling_nan()) {
+        fputil::raise_except_if_required(FE_INVALID);
+        return FPBits::quiet_nan().get_val();
+      }
+
+      if (x_abs == 0x7f80'0000U) {
+        fputil::set_errno_if_required(EDOM);
+        fputil::raise_except_if_required(FE_INVALID);
+      }
+      return x + FPBits::quiet_nan().get_val();
+    }
+#ifndef LIBC_MATH_HAS_SKIP_ACCURATE_PASS
+    // Other large exceptional values
+    if (auto r = TANF_EXCEPTS.lookup_odd(x_abs, x_sign);
+        LIBC_UNLIKELY(r.has_value()))
+      return r.value();
+#endif // !LIBC_MATH_HAS_SKIP_ACCURATE_PASS
+  }
+
+  // For |x| >= pi/32, we use the definition of tan(x) function:
+  //   tan(x) = sin(x) / cos(x)
+  // The we follow the same computations of sin(x) and cos(x) as sinf, cosf,
+  // and sincosf.
+
+  double xd = static_cast<double>(x);
+  double sin_k, cos_k, sin_y, cosm1_y;
+
+  sincosf_eval(xd, x_abs, sin_k, cos_k, sin_y, cosm1_y);
+  // tan(x) = sin(x) / cos(x)
+  //        = (sin_y * cos_k + cos_y * sin_k) / (cos_y * cos_k - sin_y * sin_k)
+  using fputil::multiply_add;
+  return static_cast<float>(
+      multiply_add(sin_y, cos_k, multiply_add(cosm1_y, sin_k, sin_k)) /
+      multiply_add(sin_y, -sin_k, multiply_add(cosm1_y, cos_k, cos_k)));
+}
+
+} // namespace double_eval
+} // namespace math
+
+} // namespace LIBC_NAMESPACE_DECL
+
+#endif // LLVM_LIBC_SRC___SUPPORT_MATH_TANF_DOUBLE_EVAL_H
diff --git a/libc/src/__support/math/tanf_float_eval.h b/libc/src/__support/math/tanf_float_eval.h
new file mode 100644
index 00000000000000..dadb30ace1b82f
--- /dev/null
+++ b/libc/src/__support/math/tanf_float_eval.h
@@ -0,0 +1,106 @@
+//===----------------------------------------------------------------------===//
+//
+// 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
+//
+//===----------------------------------------------------------------------===//
+///
+/// \file
+/// Float-only implementation of tanf.
+///
+//===----------------------------------------------------------------------===//
+
+#ifndef LLVM_LIBC_SRC___SUPPORT_MATH_TANF_FLOAT_EVAL_H
+#define LLVM_LIBC_SRC___SUPPORT_MATH_TANF_FLOAT_EVAL_H
+
+#include "sincosf_float_eval.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/macros/optimization.h" // LIBC_UNLIKELY
+
+namespace LIBC_NAMESPACE_DECL {
+namespace math {
+namespace float_eval {
+
+// Float-only tanf, requiring FMA. The maximum error, checked exhaustively in
+// round-to-nearest mode against MPFR, is 3.37 ULP.
+LIBC_INLINE float tanf(float x) {
+  using namespace sincosf_float_eval;
+  using FPBits = fputil::FPBits<float>;
+  FPBits xbits(x);
+  uint32_t x_abs = xbits.uintval() & 0x7fff'ffffU;
+
+  // |x| < 2^-12
+  if (LIBC_UNLIKELY(x_abs < 0x3980'0000U)) {
+    // For signed zeros.
+    if (LIBC_UNLIKELY(x_abs == 0U))
+      return x;
+    // When |x| < 2^-12, the relative error of the approximation tan(x) ~ x is:
+    //   |tan(x) - x| / |tan(x)| < x^2 / 3 < 2^-25 < epsilon(1)/2.
+    // fma(x, 2^-25, x) gives the correctly rounded result in all rounding modes
+    // without underflowing for subnormal x.
+    return fputil::multiply_add(x, 0x1.0p-25f, x);
+  }
+
+  float y;
+  unsigned k;
+  if (x_abs < 0x4880'0000U) {
+    k = sincosf_range_reduction_small(x, y);
+  } else {
+    // Inf or NaN
+    if (LIBC_UNLIKELY(x_abs >= 0x7f80'0000U)) {
+      if (xbits.is_signaling_nan()) {
+        fputil::raise_except_if_required(FE_INVALID);
+        return FPBits::quiet_nan().get_val();
+      }
+
+      if (x_abs == 0x7f80'0000U) {
+        fputil::set_errno_if_required(EDOM);
+        fputil::raise_except_if_required(FE_INVALID);
+      }
+      return x + FPBits::quiet_nan().get_val();
+    }
+
+    k = sincosf_range_reduction_large(x, y);
+  }
+
+  // x = k * pi/8 + y with |y| <= pi/16.
+  float sin_k = SIN_K_PI_OVER_8[k & 15];
+  // cos(k * pi/8) = sin((k + 4) * pi/8).
+  float cos_k = SIN_K_PI_OVER_8[(k + 4) & 15];
+
+  // Taylor polynomial of tan(y) for |y| <= pi/16:
+  //   tan(y) ~ y + y^3/3 + 2y^5/15 + 17y^7/315 + 62y^9/2835.
+  // The coefficients of tan's Taylor series decrease by a factor of about
+  // (2/pi)^2 per term, so the relative truncation error is bounded by:
+  //   1.02 * 1382/155925 * (pi/16)^10 < 2^-30.
+  float y_sq = y * y;
+  float y3 = y_sq * y;
+  float p = fputil::multiply_add(y_sq, 0x1.664f48p-6f, 0x1.ba1ba2p-5f);
+  p = fputil::multiply_add(y_sq, p, 0x1.111112p-3f);
+  p = fputil::multiply_add(y_sq, p, 0x1.555556p-2f);
+  // t ~ tan(y)
+  float t = fputil::multiply_add(y3, p, y);
+
+  // tan(x) = sin(x) / cos(x)
+  //        = (sin_k * cos(y) + cos_k * sin(y)) /
+  //          (cos_k * cos(y) - sin_k * sin(y))
+  //        = (sin_k + cos_k * tan(y)) / (cos_k - sin_k * tan(y)).
+  // Cancellation can only happen in the numerator when k = 0 (mod 8) and in
+  // the denominator when k = 4 (mod 8). In those cases, sin_k or cos_k is
+  // exactly 0 respectively, so the remaining term is computed with a small
+  // relative error.
+  float num = fputil::multiply_add(cos_k, t, sin_k);
+  float den = fputil::multiply_add(-sin_k, t, cos_k);
+  return num / den;
+}
+
+} // namespace float_eval
+} // namespace math
+} // namespace LIBC_NAMESPACE_DECL
+
+#endif // LLVM_LIBC_SRC___SUPPORT_MATH_TANF_FLOAT_EVAL_H
diff --git a/libc/test/src/math/CMakeLists.txt b/libc/test/src/math/CMakeLists.txt
index 921dbf220eb372..bb3e568c2d4c9f 100644
--- a/libc/test/src/math/CMakeLists.txt
+++ b/libc/test/src/math/CMakeLists.txt
@@ -206,6 +206,9 @@ add_fp_unittest(
   DEPENDS
     libc.hdr.errno_macros
     libc.src.math.tanf
+    libc.src.__support.math.tanf_double_eval
+    libc.src.__support.math.tanf_float_eval
+    libc.src.__support.macros.properties.cpu_features
     libc.src.__support.CPP.array
     libc.src.__support.FPUtil.fp_bits
 )
diff --git a/libc/test/src/math/exhaustive/CMakeLists.txt b/libc/test/src/math/exhaustive/CMakeLists.txt
index e1a91fca9b4c81..064a8a81548119 100644
--- a/libc/test/src/math/exhaustive/CMakeLists.txt
+++ b/libc/test/src/math/exhaustive/CMakeLists.txt
@@ -188,6 +188,7 @@ add_fp_unittest(
   DEPENDS
     .exhaustive_test
     libc.src.math.tanf
+    libc.src.__support.math.tanf_float_eval
     libc.src.__support.FPUtil.fp_bits
   LINK_LIBRARIES
     -lpthread
diff --git a/libc/test/src/math/exhaustive/tanf_test.cpp b/libc/test/src/math/exhaustive/tanf_test.cpp
index a3f606afb1c778..2954674ef0ede6 100644
--- a/libc/test/src/math/exhaustive/tanf_test.cpp
+++ b/libc/test/src/math/exhaustive/tanf_test.cpp
@@ -7,6 +7,7 @@
 //===----------------------------------------------------------------------===//
 
 #include "exhaustive_test.h"
+#include "src/__support/math/tanf_float_eval.h"
 #include "src/math/tanf.h"
 #include "utils/MPFRWrapper/MPFRUtils.h"
 
@@ -31,3 +32,21 @@ static constexpr uint32_t NEG_STOP = 0xff80'0000U;
 TEST_F(LlvmLibcTanfExhaustiveTest, NegativeRange) {
   test_full_range_all_roundings(NEG_START, NEG_STOP);
 }
+
+// Preserve the float implementation's 3.5 ULP MPFR bound in
+// round-to-nearest mode.
+static float tanf_float_eval(float x) {
+  return LIBC_NAMESPACE::math::float_eval::tanf(x);
+}
+
+using LlvmLibcTanfFloatExhaustiveTest =
+    LlvmLibcUnaryOpExhaustiveMathTest<float, mpfr::Operation::Tan,
+                                      tanf_float_eval, 3>;
+
+TEST_F(LlvmLibcTanfFloatExhaustiveTest, PositiveRange) {
+  test_full_range(mpfr::RoundingMode::Nearest, POS_START, POS_STOP);
+}
+
+TEST_F(LlvmLibcTanfFloatExhaustiveTest, NegativeRange) {
+  test_full_range(mpfr::RoundingMode::Nearest, 0x8000'0000U, NEG_STOP);
+}
diff --git a/libc/test/src/math/tanf_test.cpp b/libc/test/src/math/tanf_test.cpp
index fa04d870912e78..50873769b81d4b 100644
--- a/libc/test/src/math/tanf_test.cpp
+++ b/libc/test/src/math/tanf_test.cpp
@@ -10,6 +10,9 @@
 #include "hdr/math_macros.h"
 #include "src/__support/FPUtil/FPBits.h"
 #include "src/__support/macros/optimization.h"
+#include "src/__support/macros/properties/cpu_features.h"
+#include "src/__support/math/tanf_double_eval.h"
+#include "src/__support/math/tanf_float_eval.h"
 #include "src/math/tanf.h"
 #include "test/UnitTest/FPMatcher.h"
 #include "test/UnitTest/Test.h"
@@ -17,19 +20,41 @@
 #include "utils/MPFRWrapper/MPFRUtils.h"
 
 #ifdef LIBC_MATH_HAS_SKIP_ACCURATE_PASS
-#define TOLERANCE 1
+#define TOLERANCE 3
+#define FP_ASSERT ASSERT_MPFR_MATCH
+#define FP_EXPECT EXPECT_MPFR_MATCH
 #else
 #define TOLERANCE 0
+#define FP_ASSERT ASSERT_MPFR_MATCH_ALL_ROUNDING
+#define FP_EXPECT EXPECT_MPFR_MATCH_ALL_ROUNDING
 #endif // LIBC_MATH_HAS_SKIP_ACCURATE_PASS
 
 #include "hdr/stdint_proxy.h"
 
-using LlvmLibcTanfTest = LIBC_NAMESPACE::testing::FPTest<float>;
-
 using LIBC_NAMESPACE::testing::SDCOMP26094_VALUES;
 
 namespace mpfr = LIBC_NAMESPACE::testing::mpfr;
 
+class LlvmLibcTanfTest : public LIBC_NAMESPACE::testing::FPTest<float> {
+public:
+  void test_eval_in_float_range(float (*func)(float), double tolerance,
+                                bool all_rounding) {
+    constexpr uint32_t COUNT = 1'231;
+    constexpr uint32_t STEP = UINT32_MAX / COUNT;
+    for (uint32_t i = 0, v = 0; i <= COUNT; ++i, v += STEP) {
+      float x = FPBits(v).get_val();
+      if (FPBits(v).is_nan() || FPBits(v).is_inf())
+        continue;
+      if (all_rounding) {
+        ASSERT_MPFR_MATCH_ALL_ROUNDING(mpfr::Operation::Tan, x, func(x),
+                                       tolerance);
+      } else {
+        ASSERT_MPFR_MATCH(mpfr::Operation::Tan, x, func(x), tolerance);
+      }
+    }
+  }
+};
+
 TEST_F(LlvmLibcTanfTest, SpecialNumbers) {
   EXPECT_FP_EQ(aNaN, LIBC_NAMESPACE::tanf(aNaN));
   EXPECT_MATH_ERRNO(0);
@@ -54,8 +79,8 @@ TEST_F(LlvmLibcTanfTest, InFloatRange) {
     float x = FPBits(v).get_val();
     if (FPBits(v).is_nan() || FPBits(v).is_inf())
       continue;
-    ASSERT_MPFR_MATCH_ALL_ROUNDING(mpfr::Operation::Tan, x,
-                                   LIBC_NAMESPACE::tanf(x), 0.5);
+    FP_ASSERT(mpfr::Operation::Tan, x, LIBC_NAMESPACE::tanf(x),
+              TOLERANCE + 0.5);
   }
 }
 
@@ -120,10 +145,10 @@ TEST_F(LlvmLibcTanfTest, SpecificBitPatterns) {
 
   for (int i = 0; i < N; ++i) {
     float x = FPBits(INPUTS[i]).get_val();
-    EXPECT_MPFR_MATCH_ALL_ROUNDING(mpfr::Operation::Tan, x,
-                                   LIBC_NAMESPACE::tanf(x), TOLERANCE + 0.5);
-    EXPECT_MPFR_MATCH_ALL_ROUNDING(mpfr::Operation::Tan, -x,
-                                   LIBC_NAMESPACE::tanf(-x), TOLERANCE + 0.5);
+    FP_EXPECT(mpfr::Operation::Tan, x, LIBC_NAMESPACE::tanf(x),
+              TOLERANCE + 0.5);
+    FP_EXPECT(mpfr::Operation::Tan, -x, LIBC_NAMESPACE::tanf(-x),
+              TOLERANCE + 0.5);
   }
 }
 
@@ -132,6 +157,21 @@ TEST_F(LlvmLibcTanfTest, SpecificBitPatterns) {
 TEST_F(LlvmLibcTanfTest, SDCOMP_26094) {
   for (uint32_t v : SDCOMP26094_VALUES) {
     float x = FPBits(v).get_val();
-    ASSERT_MPFR_MATCH(mpfr::Operation::Tan, x, LIBC_NAMESPACE::tanf(x), 0.5);
+    ASSERT_MPFR_MATCH(mpfr::Operation::Tan, x, LIBC_NAMESPACE::tanf(x),
+                      TOLERANCE + 0.5);
   }
 }
+
+// Exercise both implementations independently of the configured selector.
+TEST_F(LlvmLibcTanfTest, DoubleEvalInFloatRange) {
+  test_eval_in_float_range(&LIBC_NAMESPACE::math::double_eval::tanf,
+                           TOLERANCE + 0.5,
+                           /*all_rounding=*/TOLERANCE == 0);
+}
+
+#ifdef LIBC_TARGET_CPU_HAS_FMA_FLOAT
+TEST_F(LlvmLibcTanfTest, FloatEvalInFloatRange) {
+  test_eval_in_float_range(&LIBC_NAMESPACE::math::float_eval::tanf, 3.5,
+                           /*all_rounding=*/false);
+}
+#endif // LIBC_TARGET_CPU_HAS_FMA_FLOAT
diff --git a/utils/bazel/llvm-project-overlay/libc/BUILD.bazel b/utils/bazel/llvm-project-overlay/libc/BUILD.bazel
index f2df1ef9466183..1e6ede0039b740 100644
--- a/utils/bazel/llvm-project-overlay/libc/BUILD.bazel
+++ b/utils/bazel/llvm-project-overlay/libc/BUILD.bazel
@@ -11475,8 +11475,13 @@ libc_support_library(
 
 libc_support_library(
     name = "__support_math_tanf",
-    hdrs = ["src/__support/math/tanf.h"],
+    hdrs = [
+        "src/__support/math/tanf.h",
+        "src/__support/math/tanf_double_eval.h",
+        "src/__support/math/tanf_float_eval.h",
+    ],
     deps = [
+        ":__support_common",
         ":__support_fputil_except_value_utils",
         ":__support_fputil_fenv_impl",
         ":__support_fputil_fma",

>From ac9626e033d026069a456d2190fc32275fc7bcbe Mon Sep 17 00:00:00 2001
From: prabhukr <prabhukr at google.com>
Date: Sat, 17 Oct 2026 16:22:52 +0000
Subject: [PATCH 2/3] [libc][math] Use include paths relative to libc/ in tanf
 headers

Address review comments.

Assisted-by: Gemini
---
 libc/src/__support/math/tanf.h             | 4 ++--
 libc/src/__support/math/tanf_double_eval.h | 2 +-
 libc/src/__support/math/tanf_float_eval.h  | 2 +-
 3 files changed, 4 insertions(+), 4 deletions(-)

diff --git a/libc/src/__support/math/tanf.h b/libc/src/__support/math/tanf.h
index 9574156801189a..90a717e5cbe6eb 100644
--- a/libc/src/__support/math/tanf.h
+++ b/libc/src/__support/math/tanf.h
@@ -22,10 +22,10 @@
     (defined(LIBC_MATH_HAS_SMALL_TABLES) ||                                    \
      defined(LIBC_MATH_HAS_INTERMEDIATE_COMP_IN_FLOAT)) &&                     \
     defined(LIBC_TARGET_CPU_HAS_FMA_FLOAT)
-#include "tanf_float_eval.h"
+#include "src/__support/math/tanf_float_eval.h"
 #define LIBC_MATH_TANF_IMPL float_eval
 #else
-#include "tanf_double_eval.h"
+#include "src/__support/math/tanf_double_eval.h"
 #define LIBC_MATH_TANF_IMPL double_eval
 #endif
 
diff --git a/libc/src/__support/math/tanf_double_eval.h b/libc/src/__support/math/tanf_double_eval.h
index 42a8a7ee9a62a1..e5ef67ed0d157c 100644
--- a/libc/src/__support/math/tanf_double_eval.h
+++ b/libc/src/__support/math/tanf_double_eval.h
@@ -14,7 +14,6 @@
 #ifndef LLVM_LIBC_SRC___SUPPORT_MATH_TANF_DOUBLE_EVAL_H
 #define LLVM_LIBC_SRC___SUPPORT_MATH_TANF_DOUBLE_EVAL_H
 
-#include "sincosf_utils.h"
 #include "src/__support/FPUtil/FEnvImpl.h"
 #include "src/__support/FPUtil/FPBits.h"
 #include "src/__support/FPUtil/PolyEval.h"
@@ -24,6 +23,7 @@
 #include "src/__support/macros/config.h"
 #include "src/__support/macros/optimization.h"            // LIBC_UNLIKELY
 #include "src/__support/macros/properties/cpu_features.h" // LIBC_TARGET_CPU_HAS_FMA
+#include "src/__support/math/sincosf_utils.h"
 
 namespace LIBC_NAMESPACE_DECL {
 
diff --git a/libc/src/__support/math/tanf_float_eval.h b/libc/src/__support/math/tanf_float_eval.h
index dadb30ace1b82f..27f49c8e35bb6d 100644
--- a/libc/src/__support/math/tanf_float_eval.h
+++ b/libc/src/__support/math/tanf_float_eval.h
@@ -14,13 +14,13 @@
 #ifndef LLVM_LIBC_SRC___SUPPORT_MATH_TANF_FLOAT_EVAL_H
 #define LLVM_LIBC_SRC___SUPPORT_MATH_TANF_FLOAT_EVAL_H
 
-#include "sincosf_float_eval.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/macros/optimization.h" // LIBC_UNLIKELY
+#include "src/__support/math/sincosf_float_eval.h"
 
 namespace LIBC_NAMESPACE_DECL {
 namespace math {

>From 7d55698a833d2f71193f5ffb99d0ce68d6ea0120 Mon Sep 17 00:00:00 2001
From: prabhukr <prabhukr at google.com>
Date: Sun, 18 Oct 2026 19:02:00 +0000
Subject: [PATCH 3/3] [libc][bazel] Add __support_math_tanf dependency to
 tanf_test

tanf_test now includes tanf_double_eval.h and tanf_float_eval.h directly
to test both implementations, so it needs a direct dependency on the
library that exports them to pass the layering check.

Assisted-by: Gemini
---
 utils/bazel/llvm-project-overlay/libc/test/src/math/BUILD.bazel | 1 +
 1 file changed, 1 insertion(+)

diff --git a/utils/bazel/llvm-project-overlay/libc/test/src/math/BUILD.bazel b/utils/bazel/llvm-project-overlay/libc/test/src/math/BUILD.bazel
index 3c6468dfd3b1e9..b8c0baec3839ac 100644
--- a/utils/bazel/llvm-project-overlay/libc/test/src/math/BUILD.bazel
+++ b/utils/bazel/llvm-project-overlay/libc/test/src/math/BUILD.bazel
@@ -887,6 +887,7 @@ math_mpfr_test(
     hdrs = ["sdcomp26094.h"],
     deps = [
         "//libc:__support_cpp_array",
+        "//libc:__support_math_tanf",
     ],
 )
 



More information about the libc-commits mailing list