https://github.com/bassiounix created https://github.com/llvm/llvm-project/pull/151883
None >From a9b861153d10cc62444d1d893a23d0292848937e Mon Sep 17 00:00:00 2001 From: bassiounix <muhammad.m.bassio...@gmail.com> Date: Sun, 3 Aug 2025 22:24:35 +0300 Subject: [PATCH] [libc][math] Refactor cos implementation to header-only in src/__support/math folder. --- libc/shared/math.h | 1 + libc/shared/math/cos.h | 23 +++ libc/src/__support/math/CMakeLists.txt | 47 +++++ libc/src/__support/math/cos.h | 173 ++++++++++++++++++ .../math}/range_reduction_double_common.h | 18 +- .../math}/range_reduction_double_fma.h | 10 +- .../math}/range_reduction_double_nofma.h | 10 +- .../generic => __support/math}/sincos_eval.h | 8 +- libc/src/math/generic/CMakeLists.txt | 50 +---- libc/src/math/generic/cos.cpp | 155 +--------------- libc/src/math/generic/sin.cpp | 14 +- libc/src/math/generic/sincos.cpp | 14 +- libc/src/math/generic/tan.cpp | 7 +- libc/test/shared/CMakeLists.txt | 1 + libc/test/shared/shared_math_test.cpp | 3 +- .../llvm-project-overlay/libc/BUILD.bazel | 101 +++++----- 16 files changed, 368 insertions(+), 267 deletions(-) create mode 100644 libc/shared/math/cos.h create mode 100644 libc/src/__support/math/cos.h rename libc/src/{math/generic => __support/math}/range_reduction_double_common.h (98%) rename libc/src/{math/generic => __support/math}/range_reduction_double_fma.h (98%) rename libc/src/{math/generic => __support/math}/range_reduction_double_nofma.h (98%) rename libc/src/{math/generic => __support/math}/sincos_eval.h (98%) diff --git a/libc/shared/math.h b/libc/shared/math.h index ea645f0afedbc..a5581ed4272a3 100644 --- a/libc/shared/math.h +++ b/libc/shared/math.h @@ -32,6 +32,7 @@ #include "math/atanhf16.h" #include "math/cbrt.h" #include "math/cbrtf.h" +#include "math/cos.h" #include "math/erff.h" #include "math/exp.h" #include "math/exp10.h" diff --git a/libc/shared/math/cos.h b/libc/shared/math/cos.h new file mode 100644 index 0000000000000..c498550f098b4 --- /dev/null +++ b/libc/shared/math/cos.h @@ -0,0 +1,23 @@ +//===-- Shared cos 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_COS_H +#define LLVM_LIBC_SHARED_MATH_COS_H + +#include "shared/libc_common.h" +#include "src/__support/math/cos.h" + +namespace LIBC_NAMESPACE_DECL { +namespace shared { + +using math::cos; + +} // namespace shared +} // namespace LIBC_NAMESPACE_DECL + +#endif // LLVM_LIBC_SHARED_MATH_COS_H diff --git a/libc/src/__support/math/CMakeLists.txt b/libc/src/__support/math/CMakeLists.txt index fe928a8fadd5e..24844063fcd24 100644 --- a/libc/src/__support/math/CMakeLists.txt +++ b/libc/src/__support/math/CMakeLists.txt @@ -357,6 +357,24 @@ add_header_library( libc.src.__support.macros.optimization ) +add_header_library( + cos + HDRS + cos.h + DEPENDS + libc.src.__support.math.sincos_eval + libc.hdr.errno_macros + libc.src.errno.errno + libc.src.__support.FPUtil.double_double + libc.src.__support.FPUtil.dyadic_float + libc.src.__support.FPUtil.except_value_utils + libc.src.__support.FPUtil.fenv_impl + libc.src.__support.FPUtil.fp_bits + libc.src.__support.math.range_reduction_double + libc.src.__support.macros.optimization +) + + add_header_library( erff HDRS @@ -613,3 +631,32 @@ add_header_library( libc.src.__support.macros.optimization libc.src.__support.macros.properties.cpu_features ) + +add_header_library( + range_reduction_double + HDRS + range_reduction_double_common.h + range_reduction_double_fma.h + range_reduction_double_nofma.h + DEPENDS + libc.src.__support.FPUtil.double_double + libc.src.__support.FPUtil.dyadic_float + libc.src.__support.FPUtil.fp_bits + libc.src.__support.FPUtil.fma + libc.src.__support.FPUtil.multiply_add + libc.src.__support.FPUtil.nearest_integer + libc.src.__support.common + libc.src.__support.integer_literals +) + +add_header_library( + sincos_eval + HDRS + sincos_eval.h + DEPENDS + libc.src.__support.FPUtil.double_double + libc.src.__support.FPUtil.dyadic_float + libc.src.__support.FPUtil.multiply_add + libc.src.__support.FPUtil.polyeval + libc.src.__support.integer_literals +) diff --git a/libc/src/__support/math/cos.h b/libc/src/__support/math/cos.h new file mode 100644 index 0000000000000..0802f9e4f6e49 --- /dev/null +++ b/libc/src/__support/math/cos.h @@ -0,0 +1,173 @@ +//===-- Implementation header for cos ---------------------------*- 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 LIBC_SRC___SUPPORT_MATH_COS_H +#define LIBC_SRC___SUPPORT_MATH_COS_H + +#include "range_reduction_double_common.h" +#include "sincos_eval.h" +#include "src/__support/FPUtil/FEnvImpl.h" +#include "src/__support/FPUtil/FPBits.h" +#include "src/__support/FPUtil/double_double.h" +#include "src/__support/FPUtil/dyadic_float.h" +#include "src/__support/FPUtil/except_value_utils.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 + +#ifdef LIBC_TARGET_CPU_HAS_FMA_DOUBLE +#include "range_reduction_double_fma.h" +#else +#include "range_reduction_double_nofma.h" +#endif // LIBC_TARGET_CPU_HAS_FMA_DOUBLE + +namespace LIBC_NAMESPACE_DECL { + +namespace math { + +LIBC_INLINE static constexpr double cos(double x) { + using namespace range_reduction_double_internal; + using DoubleDouble = fputil::DoubleDouble; + using FPBits = typename fputil::FPBits<double>; + FPBits xbits(x); + + uint16_t x_e = xbits.get_biased_exponent(); + + DoubleDouble y; + unsigned k = 0; + LargeRangeReduction range_reduction_large; + + // |x| < 2^16. + if (LIBC_LIKELY(x_e < FPBits::EXP_BIAS + FAST_PASS_EXPONENT)) { + // |x| < 2^-7 + if (LIBC_UNLIKELY(x_e < FPBits::EXP_BIAS - 7)) { + // |x| < 2^-27 + if (LIBC_UNLIKELY(x_e < FPBits::EXP_BIAS - 27)) { + // Signed zeros. + if (LIBC_UNLIKELY(x == 0.0)) + return 1.0; + + // For |x| < 2^-27, |cos(x) - 1| < |x|^2/2 < 2^-54 = ulp(1 - 2^-53)/2. + return fputil::round_result_slightly_down(1.0); + } + // No range reduction needed. + k = 0; + y.lo = 0.0; + y.hi = x; + } else { + // Small range reduction. + k = range_reduction_small(x, y); + } + } else { + // Inf or NaN + if (LIBC_UNLIKELY(x_e > 2 * FPBits::EXP_BIAS)) { + if (xbits.is_signaling_nan()) { + fputil::raise_except_if_required(FE_INVALID); + return FPBits::quiet_nan().get_val(); + } + // cos(+-Inf) = NaN + if (xbits.get_mantissa() == 0) { + fputil::set_errno_if_required(EDOM); + fputil::raise_except_if_required(FE_INVALID); + } + return x + FPBits::quiet_nan().get_val(); + } + + // Large range reduction. + k = range_reduction_large.fast(x, y); + } + + DoubleDouble sin_y, cos_y; + + [[maybe_unused]] double err = + math::sincos_eval_internal::sincos_eval(y, sin_y, cos_y); + + // Look up sin(k * pi/128) and cos(k * pi/128) +#ifdef LIBC_MATH_HAS_SMALL_TABLES + // Memory saving versions. Use 65-entry table. + auto get_idx_dd = [](unsigned kk) -> DoubleDouble { + unsigned idx = (kk & 64) ? 64 - (kk & 63) : (kk & 63); + DoubleDouble ans = SIN_K_PI_OVER_128[idx]; + if (kk & 128) { + ans.hi = -ans.hi; + ans.lo = -ans.lo; + } + return ans; + }; + DoubleDouble msin_k = get_idx_dd(k + 128); + DoubleDouble cos_k = get_idx_dd(k + 64); +#else + // Fast look up version, but needs 256-entry table. + // -sin(k * pi/128) = sin((k + 128) * pi/128) + // cos(k * pi/128) = sin(k * pi/128 + pi/2) = sin((k + 64) * pi/128). + DoubleDouble msin_k = SIN_K_PI_OVER_128[(k + 128) & 255]; + DoubleDouble cos_k = SIN_K_PI_OVER_128[(k + 64) & 255]; +#endif // LIBC_MATH_HAS_SMALL_TABLES + + // After range reduction, k = round(x * 128 / pi) and y = x - k * (pi / 128). + // So k is an integer and -pi / 256 <= y <= pi / 256. + // Then cos(x) = cos((k * pi/128 + y) + // = cos(y) * cos(k*pi/128) - sin(y) * sin(k*pi/128) + DoubleDouble cos_k_cos_y = fputil::quick_mult(cos_y, cos_k); + DoubleDouble msin_k_sin_y = fputil::quick_mult(sin_y, msin_k); + + DoubleDouble rr = fputil::exact_add<false>(cos_k_cos_y.hi, msin_k_sin_y.hi); + rr.lo += msin_k_sin_y.lo + cos_k_cos_y.lo; + +#ifdef LIBC_MATH_HAS_SKIP_ACCURATE_PASS + return rr.hi + rr.lo; +#else + using Float128 = typename fputil::DyadicFloat<128>; + double rlp = rr.lo + err; + double rlm = rr.lo - err; + + double r_upper = rr.hi + rlp; // (rr.lo + ERR); + double r_lower = rr.hi + rlm; // (rr.lo - ERR); + + // Ziv's rounding test. + if (LIBC_LIKELY(r_upper == r_lower)) + return r_upper; + + Float128 u_f128, sin_u, cos_u; + if (LIBC_LIKELY(x_e < FPBits::EXP_BIAS + FAST_PASS_EXPONENT)) + u_f128 = range_reduction_small_f128(x); + else + u_f128 = range_reduction_large.accurate(); + + math::sincos_eval_internal::sincos_eval(u_f128, sin_u, cos_u); + + auto get_sin_k = [](unsigned kk) -> Float128 { + unsigned idx = (kk & 64) ? 64 - (kk & 63) : (kk & 63); + Float128 ans = SIN_K_PI_OVER_128_F128[idx]; + if (kk & 128) + ans.sign = Sign::NEG; + return ans; + }; + + // -sin(k * pi/128) = sin((k + 128) * pi/128) + // cos(k * pi/128) = sin(k * pi/128 + pi/2) = sin((k + 64) * pi/128). + Float128 msin_k_f128 = get_sin_k(k + 128); + Float128 cos_k_f128 = get_sin_k(k + 64); + + // cos(x) = cos((k * pi/128 + u) + // = cos(u) * cos(k*pi/128) - sin(u) * sin(k*pi/128) + Float128 r = fputil::quick_add(fputil::quick_mul(cos_k_f128, cos_u), + fputil::quick_mul(msin_k_f128, sin_u)); + + // TODO: Add assertion if Ziv's accuracy tests fail in debug mode. + // https://github.com/llvm/llvm-project/issues/96452. + + return static_cast<double>(r); +#endif // !LIBC_MATH_HAS_SKIP_ACCURATE_PASS +} + +} // namespace math + +} // namespace LIBC_NAMESPACE_DECL + +#endif // LIBC_SRC___SUPPORT_MATH_COS_H diff --git a/libc/src/math/generic/range_reduction_double_common.h b/libc/src/__support/math/range_reduction_double_common.h similarity index 98% rename from libc/src/math/generic/range_reduction_double_common.h rename to libc/src/__support/math/range_reduction_double_common.h index a93ee25201813..0cace009abf90 100644 --- a/libc/src/math/generic/range_reduction_double_common.h +++ b/libc/src/__support/math/range_reduction_double_common.h @@ -20,6 +20,10 @@ namespace LIBC_NAMESPACE_DECL { +namespace math { + +namespace range_reduction_double_internal { + #ifdef LIBC_TARGET_CPU_HAS_FMA_DOUBLE static constexpr unsigned SPLIT = fputil::DefaultSplit<double>::VALUE; #else @@ -40,7 +44,7 @@ using Float128 = LIBC_NAMESPACE::fputil::DyadicFloat<128>; // Error bound: // |(x - k * pi/128) - (u_hi + u_lo)| <= max(ulp(ulp(u_hi)), 2^-119) // <= 2^-111. -LIBC_INLINE unsigned range_reduction_small(double x, DoubleDouble &u) { +LIBC_INLINE static unsigned range_reduction_small(double x, DoubleDouble &u) { // Values of -pi/128 used for inputs with absolute value <= 2^16. // The first 3 parts are generated with (53 - 21 = 32)-bit precision, so that // the product k * MPI_OVER_128[i] is exact. @@ -267,13 +271,15 @@ struct LargeRangeReduction { } #endif // !LIBC_MATH_HAS_SKIP_ACCURATE_PASS + LIBC_INLINE LargeRangeReduction() = default; + private: // Index of x in the look-up table ONE_TWENTY_EIGHT_OVER_PI. - unsigned idx; + unsigned idx = 0; // x scaled down by 2^(-16 *(idx - 3))). - double x_reduced; + double x_reduced = 0; // Parts of (x * 128/pi) mod 1. - double y_hi, y_lo; + double y_hi = 0, y_lo = 0; DoubleDouble y_mid; }; @@ -369,6 +375,10 @@ static constexpr Float128 SIN_K_PI_OVER_128_F128[65] = { }; #endif // !LIBC_MATH_HAS_SKIP_ACCURATE_PASS +} // namespace range_reduction_double_internal + +} // namespace math + } // namespace LIBC_NAMESPACE_DECL #endif // LLVM_LIBC_SRC_MATH_GENERIC_RANGE_REDUCTION_DOUBLE_COMMON_H diff --git a/libc/src/math/generic/range_reduction_double_fma.h b/libc/src/__support/math/range_reduction_double_fma.h similarity index 98% rename from libc/src/math/generic/range_reduction_double_fma.h rename to libc/src/__support/math/range_reduction_double_fma.h index 160fb2461fe21..156e51befcd04 100644 --- a/libc/src/math/generic/range_reduction_double_fma.h +++ b/libc/src/__support/math/range_reduction_double_fma.h @@ -16,10 +16,14 @@ #include "src/__support/common.h" #include "src/__support/macros/config.h" #include "src/__support/macros/optimization.h" -#include "src/math/generic/range_reduction_double_common.h" +#include "src/__support/math/range_reduction_double_common.h" namespace LIBC_NAMESPACE_DECL { +namespace math { + +namespace range_reduction_double_internal { + using LIBC_NAMESPACE::fputil::DoubleDouble; LIBC_INLINE unsigned LargeRangeReduction::fast(double x, DoubleDouble &u) { @@ -341,6 +345,10 @@ LIBC_INLINE constexpr DoubleDouble SIN_K_PI_OVER_128[] = { #endif // !LIBC_MATH_HAS_SMALL_TABLES }; +} // namespace range_reduction_double_internal + +} // namespace math + } // namespace LIBC_NAMESPACE_DECL #endif // LLVM_LIBC_SRC_MATH_GENERIC_RANGE_REDUCTION_DOUBLE_FMA_H diff --git a/libc/src/math/generic/range_reduction_double_nofma.h b/libc/src/__support/math/range_reduction_double_nofma.h similarity index 98% rename from libc/src/math/generic/range_reduction_double_nofma.h rename to libc/src/__support/math/range_reduction_double_nofma.h index 9d13d246ce91f..d89cd72e94005 100644 --- a/libc/src/math/generic/range_reduction_double_nofma.h +++ b/libc/src/__support/math/range_reduction_double_nofma.h @@ -16,10 +16,14 @@ #include "src/__support/common.h" #include "src/__support/macros/config.h" #include "src/__support/macros/optimization.h" -#include "src/math/generic/range_reduction_double_common.h" +#include "src/__support/math/range_reduction_double_common.h" namespace LIBC_NAMESPACE_DECL { +namespace math { + +namespace range_reduction_double_internal { + using fputil::DoubleDouble; LIBC_INLINE unsigned LargeRangeReduction::fast(double x, DoubleDouble &u) { @@ -342,6 +346,10 @@ LIBC_INLINE constexpr DoubleDouble SIN_K_PI_OVER_128[] = { #endif // !LIBC_MATH_HAS_SMALL_TABLES }; +} // namespace range_reduction_double_internal + +} // namespace math + } // namespace LIBC_NAMESPACE_DECL #endif // LLVM_LIBC_SRC_MATH_GENERIC_RANGE_REDUCTION_DOUBLE_NOFMA_H diff --git a/libc/src/math/generic/sincos_eval.h b/libc/src/__support/math/sincos_eval.h similarity index 98% rename from libc/src/math/generic/sincos_eval.h rename to libc/src/__support/math/sincos_eval.h index 41a4c75849ff4..fc741af19551f 100644 --- a/libc/src/math/generic/sincos_eval.h +++ b/libc/src/__support/math/sincos_eval.h @@ -18,7 +18,9 @@ namespace LIBC_NAMESPACE_DECL { -namespace generic { +namespace math { + +namespace sincos_eval_internal { using fputil::DoubleDouble; using Float128 = fputil::DyadicFloat<128>; @@ -131,7 +133,9 @@ LIBC_INLINE void sincos_eval(const Float128 &u, Float128 &sin_u, COS_COEFFS[6]); } -} // namespace generic +} // namespace sincos_eval_internal + +} // namespace math } // namespace LIBC_NAMESPACE_DECL diff --git a/libc/src/math/generic/CMakeLists.txt b/libc/src/math/generic/CMakeLists.txt index 44c17cf5b0770..b224808e54e45 100644 --- a/libc/src/math/generic/CMakeLists.txt +++ b/libc/src/math/generic/CMakeLists.txt @@ -275,23 +275,6 @@ add_header_library( libc.src.__support.common ) -add_header_library( - range_reduction_double - HDRS - range_reduction_double_common.h - range_reduction_double_fma.h - range_reduction_double_nofma.h - DEPENDS - libc.src.__support.FPUtil.double_double - libc.src.__support.FPUtil.dyadic_float - libc.src.__support.FPUtil.fp_bits - libc.src.__support.FPUtil.fma - libc.src.__support.FPUtil.multiply_add - libc.src.__support.FPUtil.nearest_integer - libc.src.__support.common - libc.src.__support.integer_literals -) - add_header_library( sincosf_utils HDRS @@ -313,18 +296,6 @@ add_header_library( libc.src.__support.common ) -add_header_library( - sincos_eval - HDRS - sincos_eval.h - DEPENDS - libc.src.__support.FPUtil.double_double - libc.src.__support.FPUtil.dyadic_float - libc.src.__support.FPUtil.multiply_add - libc.src.__support.FPUtil.polyeval - libc.src.__support.integer_literals -) - add_entrypoint_object( cos SRCS @@ -332,16 +303,7 @@ add_entrypoint_object( HDRS ../cos.h DEPENDS - .range_reduction_double - .sincos_eval - libc.hdr.errno_macros - libc.src.errno.errno - libc.src.__support.FPUtil.double_double - libc.src.__support.FPUtil.dyadic_float - libc.src.__support.FPUtil.except_value_utils - libc.src.__support.FPUtil.fenv_impl - libc.src.__support.FPUtil.fp_bits - libc.src.__support.macros.optimization + libc.src.__support.math.cos ) add_entrypoint_object( @@ -420,8 +382,8 @@ add_entrypoint_object( HDRS ../sin.h DEPENDS - .range_reduction_double - .sincos_eval + libc.src.__support.math.range_reduction_double + libc.src.__support.math.sincos_eval libc.hdr.errno_macros libc.src.errno.errno libc.src.__support.FPUtil.double_double @@ -480,8 +442,8 @@ add_entrypoint_object( HDRS ../sincos.h DEPENDS - .range_reduction_double - .sincos_eval + libc.src.__support.math.range_reduction_double + libc.src.__support.math.sincos_eval libc.hdr.errno_macros libc.src.errno.errno libc.src.__support.FPUtil.double_double @@ -553,7 +515,7 @@ add_entrypoint_object( HDRS ../tan.h DEPENDS - .range_reduction_double + libc.src.__support.math.range_reduction_double libc.hdr.errno_macros libc.src.errno.errno libc.src.__support.FPUtil.double_double diff --git a/libc/src/math/generic/cos.cpp b/libc/src/math/generic/cos.cpp index 5da0f86812a89..aabf3bc7edcb0 100644 --- a/libc/src/math/generic/cos.cpp +++ b/libc/src/math/generic/cos.cpp @@ -7,161 +7,10 @@ //===----------------------------------------------------------------------===// #include "src/math/cos.h" -#include "hdr/errno_macros.h" -#include "src/__support/FPUtil/FEnvImpl.h" -#include "src/__support/FPUtil/FPBits.h" -#include "src/__support/FPUtil/double_double.h" -#include "src/__support/FPUtil/dyadic_float.h" -#include "src/__support/FPUtil/except_value_utils.h" -#include "src/__support/common.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/math/generic/range_reduction_double_common.h" -#include "src/math/generic/sincos_eval.h" - -#ifdef LIBC_TARGET_CPU_HAS_FMA_DOUBLE -#include "range_reduction_double_fma.h" -#else -#include "range_reduction_double_nofma.h" -#endif // LIBC_TARGET_CPU_HAS_FMA_DOUBLE +#include "src/__support/math/cos.h" namespace LIBC_NAMESPACE_DECL { -using DoubleDouble = fputil::DoubleDouble; -using Float128 = typename fputil::DyadicFloat<128>; - -LLVM_LIBC_FUNCTION(double, cos, (double x)) { - using FPBits = typename fputil::FPBits<double>; - FPBits xbits(x); - - uint16_t x_e = xbits.get_biased_exponent(); - - DoubleDouble y; - unsigned k; - LargeRangeReduction range_reduction_large{}; - - // |x| < 2^16. - if (LIBC_LIKELY(x_e < FPBits::EXP_BIAS + FAST_PASS_EXPONENT)) { - // |x| < 2^-7 - if (LIBC_UNLIKELY(x_e < FPBits::EXP_BIAS - 7)) { - // |x| < 2^-27 - if (LIBC_UNLIKELY(x_e < FPBits::EXP_BIAS - 27)) { - // Signed zeros. - if (LIBC_UNLIKELY(x == 0.0)) - return 1.0; - - // For |x| < 2^-27, |cos(x) - 1| < |x|^2/2 < 2^-54 = ulp(1 - 2^-53)/2. - return fputil::round_result_slightly_down(1.0); - } - // No range reduction needed. - k = 0; - y.lo = 0.0; - y.hi = x; - } else { - // Small range reduction. - k = range_reduction_small(x, y); - } - } else { - // Inf or NaN - if (LIBC_UNLIKELY(x_e > 2 * FPBits::EXP_BIAS)) { - if (xbits.is_signaling_nan()) { - fputil::raise_except_if_required(FE_INVALID); - return FPBits::quiet_nan().get_val(); - } - // cos(+-Inf) = NaN - if (xbits.get_mantissa() == 0) { - fputil::set_errno_if_required(EDOM); - fputil::raise_except_if_required(FE_INVALID); - } - return x + FPBits::quiet_nan().get_val(); - } - - // Large range reduction. - k = range_reduction_large.fast(x, y); - } - - DoubleDouble sin_y, cos_y; - - [[maybe_unused]] double err = generic::sincos_eval(y, sin_y, cos_y); - - // Look up sin(k * pi/128) and cos(k * pi/128) -#ifdef LIBC_MATH_HAS_SMALL_TABLES - // Memory saving versions. Use 65-entry table. - auto get_idx_dd = [](unsigned kk) -> DoubleDouble { - unsigned idx = (kk & 64) ? 64 - (kk & 63) : (kk & 63); - DoubleDouble ans = SIN_K_PI_OVER_128[idx]; - if (kk & 128) { - ans.hi = -ans.hi; - ans.lo = -ans.lo; - } - return ans; - }; - DoubleDouble msin_k = get_idx_dd(k + 128); - DoubleDouble cos_k = get_idx_dd(k + 64); -#else - // Fast look up version, but needs 256-entry table. - // -sin(k * pi/128) = sin((k + 128) * pi/128) - // cos(k * pi/128) = sin(k * pi/128 + pi/2) = sin((k + 64) * pi/128). - DoubleDouble msin_k = SIN_K_PI_OVER_128[(k + 128) & 255]; - DoubleDouble cos_k = SIN_K_PI_OVER_128[(k + 64) & 255]; -#endif // LIBC_MATH_HAS_SMALL_TABLES - - // After range reduction, k = round(x * 128 / pi) and y = x - k * (pi / 128). - // So k is an integer and -pi / 256 <= y <= pi / 256. - // Then cos(x) = cos((k * pi/128 + y) - // = cos(y) * cos(k*pi/128) - sin(y) * sin(k*pi/128) - DoubleDouble cos_k_cos_y = fputil::quick_mult(cos_y, cos_k); - DoubleDouble msin_k_sin_y = fputil::quick_mult(sin_y, msin_k); - - DoubleDouble rr = fputil::exact_add<false>(cos_k_cos_y.hi, msin_k_sin_y.hi); - rr.lo += msin_k_sin_y.lo + cos_k_cos_y.lo; - -#ifdef LIBC_MATH_HAS_SKIP_ACCURATE_PASS - return rr.hi + rr.lo; -#else - - double rlp = rr.lo + err; - double rlm = rr.lo - err; - - double r_upper = rr.hi + rlp; // (rr.lo + ERR); - double r_lower = rr.hi + rlm; // (rr.lo - ERR); - - // Ziv's rounding test. - if (LIBC_LIKELY(r_upper == r_lower)) - return r_upper; - - Float128 u_f128, sin_u, cos_u; - if (LIBC_LIKELY(x_e < FPBits::EXP_BIAS + FAST_PASS_EXPONENT)) - u_f128 = range_reduction_small_f128(x); - else - u_f128 = range_reduction_large.accurate(); - - generic::sincos_eval(u_f128, sin_u, cos_u); - - auto get_sin_k = [](unsigned kk) -> Float128 { - unsigned idx = (kk & 64) ? 64 - (kk & 63) : (kk & 63); - Float128 ans = SIN_K_PI_OVER_128_F128[idx]; - if (kk & 128) - ans.sign = Sign::NEG; - return ans; - }; - - // -sin(k * pi/128) = sin((k + 128) * pi/128) - // cos(k * pi/128) = sin(k * pi/128 + pi/2) = sin((k + 64) * pi/128). - Float128 msin_k_f128 = get_sin_k(k + 128); - Float128 cos_k_f128 = get_sin_k(k + 64); - - // cos(x) = cos((k * pi/128 + u) - // = cos(u) * cos(k*pi/128) - sin(u) * sin(k*pi/128) - Float128 r = fputil::quick_add(fputil::quick_mul(cos_k_f128, cos_u), - fputil::quick_mul(msin_k_f128, sin_u)); - - // TODO: Add assertion if Ziv's accuracy tests fail in debug mode. - // https://github.com/llvm/llvm-project/issues/96452. - - return static_cast<double>(r); -#endif // !LIBC_MATH_HAS_SKIP_ACCURATE_PASS -} +LLVM_LIBC_FUNCTION(double, cos, (double x)) { return math::cos(x); } } // namespace LIBC_NAMESPACE_DECL diff --git a/libc/src/math/generic/sin.cpp b/libc/src/math/generic/sin.cpp index a614427bd7ee3..1b6310f99366c 100644 --- a/libc/src/math/generic/sin.cpp +++ b/libc/src/math/generic/sin.cpp @@ -18,13 +18,13 @@ #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/math/generic/range_reduction_double_common.h" -#include "src/math/generic/sincos_eval.h" +#include "src/__support/math/range_reduction_double_common.h" +#include "src/__support/math/sincos_eval.h" #ifdef LIBC_TARGET_CPU_HAS_FMA_DOUBLE -#include "range_reduction_double_fma.h" +#include "src/__support/math/range_reduction_double_fma.h" #else -#include "range_reduction_double_nofma.h" +#include "src/__support/math/range_reduction_double_nofma.h" #endif // LIBC_TARGET_CPU_HAS_FMA_DOUBLE namespace LIBC_NAMESPACE_DECL { @@ -33,6 +33,7 @@ using DoubleDouble = fputil::DoubleDouble; using Float128 = typename fputil::DyadicFloat<128>; LLVM_LIBC_FUNCTION(double, sin, (double x)) { + using namespace math::range_reduction_double_internal; using FPBits = typename fputil::FPBits<double>; FPBits xbits(x); @@ -95,7 +96,8 @@ LLVM_LIBC_FUNCTION(double, sin, (double x)) { DoubleDouble sin_y, cos_y; - [[maybe_unused]] double err = generic::sincos_eval(y, sin_y, cos_y); + [[maybe_unused]] double err = + math::sincos_eval_internal::sincos_eval(y, sin_y, cos_y); // Look up sin(k * pi/128) and cos(k * pi/128) #ifdef LIBC_MATH_HAS_SMALL_TABLES @@ -149,7 +151,7 @@ LLVM_LIBC_FUNCTION(double, sin, (double x)) { else u_f128 = range_reduction_large.accurate(); - generic::sincos_eval(u_f128, sin_u, cos_u); + math::sincos_eval_internal::sincos_eval(u_f128, sin_u, cos_u); auto get_sin_k = [](unsigned kk) -> Float128 { unsigned idx = (kk & 64) ? 64 - (kk & 63) : (kk & 63); diff --git a/libc/src/math/generic/sincos.cpp b/libc/src/math/generic/sincos.cpp index 08c8a8298f029..38661deaa886a 100644 --- a/libc/src/math/generic/sincos.cpp +++ b/libc/src/math/generic/sincos.cpp @@ -19,13 +19,13 @@ #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/math/generic/range_reduction_double_common.h" -#include "src/math/generic/sincos_eval.h" +#include "src/__support/math/range_reduction_double_common.h" +#include "src/__support/math/sincos_eval.h" #ifdef LIBC_TARGET_CPU_HAS_FMA_DOUBLE -#include "range_reduction_double_fma.h" +#include "src/__support/math/range_reduction_double_fma.h" #else -#include "range_reduction_double_nofma.h" +#include "src/__support/math/range_reduction_double_nofma.h" #endif // LIBC_TARGET_CPU_HAS_FMA_DOUBLE namespace LIBC_NAMESPACE_DECL { @@ -34,6 +34,7 @@ using DoubleDouble = fputil::DoubleDouble; using Float128 = typename fputil::DyadicFloat<128>; LLVM_LIBC_FUNCTION(void, sincos, (double x, double *sin_x, double *cos_x)) { + using namespace math::range_reduction_double_internal; using FPBits = typename fputil::FPBits<double>; FPBits xbits(x); @@ -106,7 +107,8 @@ LLVM_LIBC_FUNCTION(void, sincos, (double x, double *sin_x, double *cos_x)) { DoubleDouble sin_y, cos_y; - [[maybe_unused]] double err = generic::sincos_eval(y, sin_y, cos_y); + [[maybe_unused]] double err = + math::sincos_eval_internal::sincos_eval(y, sin_y, cos_y); // Look up sin(k * pi/128) and cos(k * pi/128) #ifdef LIBC_MATH_HAS_SMALL_TABLES @@ -179,7 +181,7 @@ LLVM_LIBC_FUNCTION(void, sincos, (double x, double *sin_x, double *cos_x)) { else u_f128 = range_reduction_large.accurate(); - generic::sincos_eval(u_f128, sin_u, cos_u); + math::sincos_eval_internal::sincos_eval(u_f128, sin_u, cos_u); auto get_sin_k = [](unsigned kk) -> Float128 { unsigned idx = (kk & 64) ? 64 - (kk & 63) : (kk & 63); diff --git a/libc/src/math/generic/tan.cpp b/libc/src/math/generic/tan.cpp index 89b812cfc23a0..7ea40c9af75d2 100644 --- a/libc/src/math/generic/tan.cpp +++ b/libc/src/math/generic/tan.cpp @@ -20,12 +20,12 @@ #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/math/generic/range_reduction_double_common.h" +#include "src/__support/math/range_reduction_double_common.h" #ifdef LIBC_TARGET_CPU_HAS_FMA_DOUBLE -#include "range_reduction_double_fma.h" +#include "src/__support/math/range_reduction_double_fma.h" #else -#include "range_reduction_double_nofma.h" +#include "src/__support/math/range_reduction_double_nofma.h" #endif // LIBC_TARGET_CPU_HAS_FMA_DOUBLE namespace LIBC_NAMESPACE_DECL { @@ -121,6 +121,7 @@ LIBC_INLINE double tan_eval(const DoubleDouble &u, DoubleDouble &result) { } // anonymous namespace LLVM_LIBC_FUNCTION(double, tan, (double x)) { + using namespace math::range_reduction_double_internal; using FPBits = typename fputil::FPBits<double>; FPBits xbits(x); diff --git a/libc/test/shared/CMakeLists.txt b/libc/test/shared/CMakeLists.txt index 9685aeaee7d8d..b2264e17ef70c 100644 --- a/libc/test/shared/CMakeLists.txt +++ b/libc/test/shared/CMakeLists.txt @@ -28,6 +28,7 @@ add_fp_unittest( libc.src.__support.math.atanhf16 libc.src.__support.math.cbrt libc.src.__support.math.cbrtf + libc.src.__support.math.cos 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 5e57c49ed23d9..2386b552d181a 100644 --- a/libc/test/shared/shared_math_test.cpp +++ b/libc/test/shared/shared_math_test.cpp @@ -66,8 +66,9 @@ TEST(LlvmLibcSharedMathTest, AllDouble) { EXPECT_FP_EQ(0x1.921fb54442d18p+0, LIBC_NAMESPACE::shared::acos(0.0)); EXPECT_FP_EQ(0x0p+0, LIBC_NAMESPACE::shared::asin(0.0)); EXPECT_FP_EQ(0x0p+0, LIBC_NAMESPACE::shared::atan(0.0)); - EXPECT_FP_EQ(0x0p+0, LIBC_NAMESPACE::shared::cbrt(0.0)); EXPECT_FP_EQ(0x0p+0, LIBC_NAMESPACE::shared::atan2(0.0, 0.0)); + EXPECT_FP_EQ(0x0p+0, LIBC_NAMESPACE::shared::cbrt(0.0)); + EXPECT_FP_EQ(0x1p+0, LIBC_NAMESPACE::shared::cos(0.0)); EXPECT_FP_EQ(0x1p+0, LIBC_NAMESPACE::shared::exp(0.0)); EXPECT_FP_EQ(0x1p+0, LIBC_NAMESPACE::shared::exp10(0.0)); } diff --git a/utils/bazel/llvm-project-overlay/libc/BUILD.bazel b/utils/bazel/llvm-project-overlay/libc/BUILD.bazel index d30c725524d17..40f7c168e75ab 100644 --- a/utils/bazel/llvm-project-overlay/libc/BUILD.bazel +++ b/utils/bazel/llvm-project-overlay/libc/BUILD.bazel @@ -1967,41 +1967,6 @@ libc_support_library( ], ) -libc_support_library( - name = "range_reduction_double", - hdrs = [ - "src/math/generic/range_reduction_double_common.h", - "src/math/generic/range_reduction_double_fma.h", - "src/math/generic/range_reduction_double_nofma.h", - ], - deps = [ - ":__support_common", - ":__support_fputil_double_double", - ":__support_fputil_dyadic_float", - ":__support_fputil_fp_bits", - ":__support_fputil_multiply_add", - ":__support_fputil_nearest_integer", - ":__support_fputil_rounding_mode", - ":__support_integer_literals", - ], -) - -libc_support_library( - name = "sincos_eval", - hdrs = [ - "src/math/generic/sincos_eval.h", - ], - deps = [ - ":__support_common", - ":__support_fputil_double_double", - ":__support_fputil_dyadic_float", - ":__support_fputil_fp_bits", - ":__support_fputil_multiply_add", - ":__support_fputil_polyeval", - ":__support_integer_literals", - ], -) - libc_support_library( name = "sincosf_utils", hdrs = ["src/math/generic/sincosf_utils.h"], @@ -2424,6 +2389,19 @@ libc_support_library( ], ) +libc_support_library( + name = "__support_math_cos", + hdrs = ["src/__support/math/cos.h"], + deps = [ + ":__support_fputil_multiply_add", + ":__support_fputil_except_value_utils", + ":__support_macros_optimization", + ":__support_macros_properties_cpu_features", + ":__support_range_reduction_double", + ":__support_sincos_eval", + ], +) + libc_support_library( name = "__support_math_erff", hdrs = ["src/__support/math/erff.h"], @@ -2678,6 +2656,41 @@ libc_support_library( ], ) +libc_support_library( + name = "__support_range_reduction_double", + hdrs = [ + "src/__support/math/range_reduction_double_common.h", + "src/__support/math/range_reduction_double_fma.h", + "src/__support/math/range_reduction_double_nofma.h", + ], + deps = [ + ":__support_common", + ":__support_fputil_double_double", + ":__support_fputil_dyadic_float", + ":__support_fputil_fp_bits", + ":__support_fputil_multiply_add", + ":__support_fputil_nearest_integer", + ":__support_fputil_rounding_mode", + ":__support_integer_literals", + ], +) + +libc_support_library( + name = "__support_sincos_eval", + hdrs = [ + "src/__support/math/sincos_eval.h", + ], + deps = [ + ":__support_common", + ":__support_fputil_double_double", + ":__support_fputil_dyadic_float", + ":__support_fputil_fp_bits", + ":__support_fputil_multiply_add", + ":__support_fputil_polyeval", + ":__support_integer_literals", + ], +) + ############################### complex targets ################################ libc_function( @@ -3127,11 +3140,7 @@ libc_math_function(name = "copysignf16") libc_math_function( name = "cos", additional_deps = [ - ":__support_fputil_multiply_add", - ":__support_macros_optimization", - ":__support_macros_properties_cpu_features", - ":range_reduction_double", - ":sincos_eval", + ":__support_math_cos", ], ) @@ -4251,8 +4260,8 @@ libc_math_function( ":__support_fputil_multiply_add", ":__support_macros_optimization", ":__support_macros_properties_cpu_features", - ":range_reduction_double", - ":sincos_eval", + ":__support_range_reduction_double", + ":__support_sincos_eval", ], ) @@ -4284,8 +4293,8 @@ libc_math_function( ":__support_fputil_multiply_add", ":__support_macros_optimization", ":__support_macros_properties_cpu_features", - ":range_reduction_double", - ":sincos_eval", + ":__support_range_reduction_double", + ":__support_sincos_eval", ], ) @@ -4379,8 +4388,8 @@ libc_math_function( ":__support_fputil_multiply_add", ":__support_macros_optimization", ":__support_macros_properties_cpu_features", - ":range_reduction_double", - ":sincos_eval", + ":__support_range_reduction_double", + ":__support_sincos_eval", ], ) _______________________________________________ llvm-branch-commits mailing list llvm-branch-commits@lists.llvm.org https://lists.llvm.org/cgi-bin/mailman/listinfo/llvm-branch-commits