Merge pull request #1145 from cppalliance/1137

Handle Git issue 1137

Author: mborland

include/boost/decimal/detail/cmath/asin.hpp

@@ -1,4 +1,5 @@
-// Copyright 2024 Matt Borland
+// Copyright 2024 - 2025 Matt Borland
+// Copyright 2024 - 2025 Christopher Kormanyos
 // Distributed under the Boost Software License, Version 1.0.
 // https://www.boost.org/LICENSE_1_0.txt
 
@@ -41,29 +42,42 @@ constexpr auto asin_impl(const T x) noexcept
         return x;
     }
 
-    const auto absx {fabs(x)};
-    T result {};
+    const auto absx { fabs(x) };
 
-    if (absx <= std::numeric_limits<T>::epsilon())
+    T result { };
+
+    constexpr T cbrt_eps { cbrt(std::numeric_limits<T>::epsilon()) };
+
+    constexpr T one { 1 };
+
+    if (absx <= cbrt_eps)
     {
-        result = absx;
+        result = absx * (one + (absx / 6) * absx);
     }
-    else if (absx <= T{5, -1})
+    else if (absx <= T { 5, -1 })
     {
         result = asin_series(absx);
     }
-    else if (absx <= T{1, 0})
-    {
-        constexpr T half_pi {numbers::pi_v<T> / 2};
-        result = half_pi - 2 * asin_series(sqrt((1 - absx) / 2));
-    }
     else
     {
-        #ifndef BOOST_DECIMAL_FAST_MATH
-        result = std::numeric_limits<T>::quiet_NaN();
-        #else
-        result = T{0};
-        #endif
+        constexpr T half_pi { numbers::pi_v<T> / 2 };
+
+        if (absx < one)
+        {
+            result = half_pi - 2 * asin_series(sqrt((1 - absx) / 2));
+        }
+        else if (absx > one)
+        {
+            #ifndef BOOST_DECIMAL_FAST_MATH
+            result = std::numeric_limits<T>::quiet_NaN();
+            #else
+            result = T{0};
+            #endif
+        }
+        else
+        {
+            result = half_pi;
+        }
     }
 
     // arcsin(-x) == -arcsin(x)

test/test_asin.cpp

@@ -1,4 +1,5 @@
-// Copyright 2024 Matt Borland
+// Copyright 2024 - 2025 Matt Borland
+// Copyright 2024 - 2025 Christopher Kormanyos
 // Distributed under the Boost Software License, Version 1.0.
 // https://www.boost.org/LICENSE_1_0.txt
 
@@ -14,7 +15,7 @@
 #  pragma clang diagnostic ignored "-Wconversion"
 #  pragma clang diagnostic ignored "-Wsign-conversion"
 #  pragma clang diagnostic ignored "-Wfloat-equal"
-#  if __clang_major__ >= 20
+#  if (__clang_major__ >= 20)
 #    pragma clang diagnostic ignored "-Wfortify-source"
 #  endif
 #elif defined(__GNUC__)
@@ -28,9 +29,9 @@
 
 #include <boost/math/special_functions/next.hpp>
 #include <boost/core/lightweight_test.hpp>
-#include <iostream>
-#include <random>
+
 #include <cmath>
+#include <random>
 
 #if !defined(BOOST_DECIMAL_REDUCE_TEST_DEPTH)
 static constexpr auto N = static_cast<std::size_t>(128U); // Number of trials
@@ -42,6 +43,9 @@ static std::mt19937_64 rng(42);
 
 using namespace boost::decimal;
 
+template<typename T> auto my_zero() -> T;
+template<typename T> auto my_one () -> T;
+
 template <typename Dec>
 void test_asin()
 {
@@ -137,6 +141,100 @@ void print_value(T value, const char* str)
               << "\nExp: " << ptr << "\n" << std::endl;
 }
 
+template<typename T>
+auto test_asin_edge() -> void
+{
+    using nl = std::numeric_limits<T>;
+
+    const T tiny0 { nl::epsilon() * 999 / 1000 };
+    const T tiny1 { nl::epsilon() };
+    const T tiny2 { nl::epsilon() * 1000 / 999 };
+
+    const T asin_tiny0 { asin(tiny0) };
+    const T asin_tiny1 { asin(tiny1) };
+    const T asin_tiny2 { asin(tiny2) };
+
+    // tiny1: 1
+    // tiny2: 1.001001
+    // tiny0: 0.999
+    // tiny1: 1
+    // tiny2: 1.001001001001001
+    // tiny0: 0.999
+    // tiny1: 1
+    // tiny2: 1.001001001001001001001001001001001
+
+    constexpr T ctrl_tiny2
+    {
+          std::numeric_limits<T>::digits10 < 10 ? T("1.001001")
+        : std::numeric_limits<T>::digits10 < 20 ? T("1.001001001001001")
+        :                                         T("1.001001001001001001001001001001001")
+    };
+
+    BOOST_TEST_EQ(asin_tiny0 / nl::epsilon(), T(999, -3));
+    BOOST_TEST_EQ(asin_tiny1 / nl::epsilon(), T(1));
+    BOOST_TEST_EQ(asin_tiny2 / nl::epsilon(), ctrl_tiny2);
+
+    constexpr T half_pi { numbers::pi_v<T> / 2 };
+
+    BOOST_TEST_EQ(asin(my_zero<T>() + my_one<T>()), half_pi);
+    BOOST_TEST_EQ(asin(my_zero<T>() - my_one<T>()), -half_pi);
+}
+
+template<typename T>
+void test_asin_1137()
+{
+    using nl = std::numeric_limits<T>;
+
+    const T tiny0 { nl::epsilon() * 999/1000 };
+    const T tiny1 { nl::epsilon() };
+    const T tiny2 { nl::epsilon() * 1000/999 };
+
+    BOOST_TEST(tiny0 != tiny1);
+    BOOST_TEST(tiny1 != tiny2);
+
+    std::stringstream strm { };
+
+    BOOST_TEST_EQ(tiny0, asin(tiny0));
+    BOOST_TEST_EQ(tiny1, asin(tiny1));
+    BOOST_TEST_EQ(tiny2, asin(tiny2));
+
+    const T sqrt_tiny0 { sqrt(nl::epsilon() * 999/1000) };
+    const T sqrt_tiny1 { sqrt(nl::epsilon()) };
+    const T sqrt_tiny2 { sqrt(nl::epsilon() * 1000/999) };
+
+    BOOST_TEST_EQ(sqrt_tiny0, asin(sqrt_tiny0));
+    BOOST_TEST_EQ(sqrt_tiny1, asin(sqrt_tiny1));
+    BOOST_TEST_EQ(sqrt_tiny2, asin(sqrt_tiny2));
+
+    const T cbrt_tiny0 { cbrt(nl::epsilon() * 999/1000) };
+    const T cbrt_tiny1 { cbrt(nl::epsilon()) };
+    const T cbrt_tiny2 { cbrt(nl::epsilon() * 1000/999) };
+    const T cbrt_tiny3 { cbrt(nl::epsilon() * 1004/999) };
+
+    auto mini_series
+    {
+        [](const T eps)
+        {
+            return eps * (1 + (eps / 6) * eps);
+        }
+    };
+
+    auto is_close
+    {
+        [](const T a, const T b)
+        {
+            const T delta { fabs(a - b) };
+
+            return (delta < (std::numeric_limits<T>::epsilon() * 4));
+        }
+    };
+
+    BOOST_TEST(is_close(asin(cbrt_tiny0), mini_series(cbrt_tiny0)));
+    BOOST_TEST(is_close(asin(cbrt_tiny1), mini_series(cbrt_tiny1)));
+    BOOST_TEST(is_close(asin(cbrt_tiny2), mini_series(cbrt_tiny2)));
+    BOOST_TEST(is_close(asin(cbrt_tiny3), mini_series(cbrt_tiny3)));
+}
+
 int main()
 {
     #ifdef BOOST_DECIMAL_GENERATE_CONSTANT_SIGS
@@ -187,12 +285,22 @@ int main()
 
     test_asin<decimal32_t>();
     test_asin<decimal64_t>();
-
     #if !defined(BOOST_DECIMAL_REDUCE_TEST_DEPTH)
     test_asin<decimal128_t>();
     #endif
 
     test_asin<decimal_fast32_t>();
 
+    test_asin_edge<decimal32_t>();
+    test_asin_edge<decimal64_t>();
+    test_asin_edge<decimal128_t>();
+
+    test_asin_1137<decimal32_t>();
+    test_asin_1137<decimal64_t>();
+    test_asin_1137<decimal128_t>();
+
     return boost::report_errors();
 }
+
+template<typename T> auto my_zero() -> T { T zero { 0 }; return zero; }
+template<typename T> auto my_one () -> T { T one  { 1 }; return one; }