Prelim finish of log10 dedicated impl

Author: ckormanyos

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

@@ -32,32 +32,39 @@ constexpr auto log10_impl(T x) noexcept
 
     T result { };
 
-    if (isnan(x))
+    const auto fpc = fpclassify(x);
+
+    if (fpc == FP_ZERO)
+    {
+        result = -std::numeric_limits<T>::infinity();
+    }
+    else if (signbit(x) || (fpc == FP_NAN))
     {
-        result = x;
+        result = std::numeric_limits<T>::quiet_NaN();
     }
-    else if (isinf(x))
+    else if (fpc == FP_INFINITE)
     {
-        result = (!signbit(x)) ? x: std::numeric_limits<T>::quiet_NaN();
+        result = std::numeric_limits<T>::infinity();
     }
     else
     {
         int exp10val { };
 
-        T g { frexp10(x, &exp10val) };
+        auto gn { frexp10(x, &exp10val) };
 
-        while(g > one)
-        {
-          g /= 10;
-
-          ++exp10val;
-        }
+        const auto zeros_removal_result { detail::remove_trailing_zeros(gn) };
 
-        const bool is_pure = (g == one);
+        const bool is_pure { (zeros_removal_result.trimmed_number == one) };
 
         if(is_pure)
         {
-            result = T { exp10val };
+            // Here, a pure power-of-10 argument gets a pure integral result.
+            result =
+                T
+                {
+                      exp10val
+                    + static_cast<int>(zeros_removal_result.number_of_removed_zeros)
+                };
         }
         else
         {
@@ -70,8 +77,6 @@ constexpr auto log10_impl(T x) noexcept
                 }
                 else if ((x == zero) || (-x == zero))
                 {
-                    // Actually, this should be equivalent to -HUGE_VAL.
-
                     result = -std::numeric_limits<T>::infinity();
                 }
                 else
@@ -84,6 +89,11 @@ constexpr auto log10_impl(T x) noexcept
                 // The algorithm for base-10 logarithm is based on Chapter 5, pages 35-36
                 // of Cody and Waite, Software Manual for the Elementary Functions,
                 // Prentice Hall, 1980.
+
+                T g { T { gn } / pow10(std::numeric_limits<T>::digits10) };
+
+                exp10val += std::numeric_limits<T>::digits10;
+
                 constexpr T inv_sqrt10 { UINT64_C(3162277660168379332), -19 };
 
                 const bool reduce_sqrt10 { g < inv_sqrt10 };
@@ -118,6 +128,7 @@ constexpr auto log10_impl(T x) noexcept
             }
         }
     }
+
     return result;
 }
 

test/Jamfile

@@ -106,8 +106,8 @@ run test_legendre.cpp ;
 run test_literals.cpp ;
 run test_lgamma.cpp ;
 run test_log.cpp ;
-run test_log10.cpp ;
 run test_log1p.cpp ;
+run test_log10.cpp ;
 run test_pow.cpp ;
 run test_promotion.cpp ;
 run test_remainder_remquo.cpp ;

test/test_log10.cpp

@@ -22,6 +22,10 @@
 
 #include <boost/core/lightweight_test.hpp>
 
+template<typename DecimalType> auto my_zero() -> DecimalType&;
+template<typename DecimalType> auto my_one () -> DecimalType&;
+template<typename DecimalType> auto my_inf () -> DecimalType&;
+
 namespace local
 {
   template<typename IntegralTimePointType,
@@ -182,6 +186,114 @@ namespace local
     return result_is_ok;
   }
 
+  template<typename DecimalType, typename FloatType>
+  auto test_log10_edge() -> bool
+  {
+    using decimal_type = DecimalType;
+    using float_type   = FloatType;
+
+    std::mt19937_64 gen;
+
+    gen.seed(time_point<typename std::mt19937_64::result_type>());
+
+    std::uniform_real_distribution<float_type> dist(1.0F, 2.0F);
+
+    volatile auto result_is_ok = true;
+
+    for(auto index = static_cast<unsigned>(UINT8_C(0)); index < static_cast<unsigned>(UINT8_C(4)); ++index)
+    {
+      static_cast<void>(index);
+
+      const auto log_zero = log10(::my_zero<decimal_type>() * static_cast<decimal_type>(dist(gen)));
+
+      const volatile auto result_log_zero_is_ok = (isinf(log_zero) && (log_zero < ::my_zero<decimal_type>()));
+
+      BOOST_TEST(result_log_zero_is_ok);
+
+      result_is_ok = (result_log_zero_is_ok && result_is_ok);
+    }
+
+    for(auto index = static_cast<unsigned>(UINT8_C(0)); index < static_cast<unsigned>(UINT8_C(4)); ++index)
+    {
+      static_cast<void>(index);
+
+      const auto log_zero_minus = log10(-::my_zero<decimal_type>() * static_cast<decimal_type>(dist(gen)));
+
+      const volatile auto result_log_zero_minus_is_ok = (isinf(log_zero_minus) && (log_zero_minus < ::my_zero<decimal_type>()));
+
+      BOOST_TEST(result_log_zero_minus_is_ok);
+
+      result_is_ok = (result_log_zero_minus_is_ok && result_is_ok);
+    }
+
+    for(auto index = static_cast<unsigned>(UINT8_C(0)); index < static_cast<unsigned>(UINT8_C(4)); ++index)
+    {
+      static_cast<void>(index);
+
+      const auto log_one = log10(::my_one<decimal_type>());
+
+      const volatile auto result_log_one_is_ok = (log_one == ::my_zero<decimal_type>() * static_cast<decimal_type>(dist(gen)));
+
+      BOOST_TEST(result_log_one_is_ok);
+
+      result_is_ok = (result_log_one_is_ok && result_is_ok);
+    }
+
+    for(auto index = static_cast<unsigned>(UINT8_C(0)); index < static_cast<unsigned>(UINT8_C(4)); ++index)
+    {
+      static_cast<void>(index);
+
+      const auto log_one_minus = log10(-::my_one<decimal_type>());
+
+      const volatile auto result_log_one_minus_is_ok = isnan(log_one_minus);
+
+      BOOST_TEST(result_log_one_minus_is_ok);
+
+      result_is_ok = (result_log_one_minus_is_ok && result_is_ok);
+    }
+
+    for(auto index = static_cast<unsigned>(UINT8_C(0)); index < static_cast<unsigned>(UINT8_C(4)); ++index)
+    {
+      static_cast<void>(index);
+
+      const auto log_inf = log10(::my_inf<decimal_type>() * static_cast<decimal_type>(dist(gen)));
+
+      const volatile auto result_log_inf_is_ok = isinf(log_inf);
+
+      BOOST_TEST(result_log_inf_is_ok);
+
+      result_is_ok = (result_log_inf_is_ok && result_is_ok);
+    }
+
+    for(auto index = static_cast<unsigned>(UINT8_C(0)); index < static_cast<unsigned>(UINT8_C(4)); ++index)
+    {
+      static_cast<void>(index);
+
+      const auto log_inf_minus = log10(-::my_inf<decimal_type>() * static_cast<decimal_type>(dist(gen)));
+
+      const volatile auto result_log_inf_minus_is_ok = isnan(log_inf_minus);
+
+      BOOST_TEST(result_log_inf_minus_is_ok);
+
+      result_is_ok = (result_log_inf_minus_is_ok && result_is_ok);
+    }
+
+    for(auto index = static_cast<unsigned>(UINT8_C(0)); index < static_cast<unsigned>(UINT8_C(4)); ++index)
+    {
+      static_cast<void>(index);
+
+      const auto log_nan = log10(std::numeric_limits<decimal_type>::quiet_NaN() * static_cast<decimal_type>(dist(gen)));
+
+      const volatile auto result_log_nan_is_ok = isnan(log_nan);
+
+      BOOST_TEST(result_log_nan_is_ok);
+
+      result_is_ok = (result_log_nan_is_ok && result_is_ok);
+    }
+
+    return result_is_ok;
+  }
+
 } // namespace local
 
 auto main() -> int
@@ -194,6 +306,8 @@ auto main() -> int
 
     const auto test_log10_is_ok = local::test_log10<decimal_type, float_type>(128);
 
+    BOOST_TEST(test_log10_is_ok);
+
     result_is_ok = (test_log10_is_ok && result_is_ok);
   }
 
@@ -203,10 +317,53 @@ auto main() -> int
 
     const auto test_log10_pow10_is_ok = local::test_log10_pow10<decimal_type, float_type>();
 
+    BOOST_TEST(test_log10_pow10_is_ok);
+
     result_is_ok = (test_log10_pow10_is_ok && result_is_ok);
   }
 
+  {
+    using decimal_type = boost::decimal::decimal32;
+    using float_type   = float;
+
+    const auto test_log10_edge_is_ok = local::test_log10_edge<decimal_type, float_type>();
+
+    BOOST_TEST(test_log10_edge_is_ok);
+
+    result_is_ok = (test_log10_edge_is_ok && result_is_ok);
+  }
+
   result_is_ok = ((boost::report_errors() == 0) && result_is_ok);
 
   return (result_is_ok ? 0 : -1);
 }
+
+template<typename DecimalType>
+auto my_zero() -> DecimalType&
+{
+  using decimal_type = DecimalType;
+
+  static decimal_type val_zero { 0, 0 };
+
+  return val_zero;
+}
+
+template<typename DecimalType>
+auto my_one() -> DecimalType&
+{
+  using decimal_type = DecimalType;
+
+  static decimal_type val_one { 1, 0 };
+
+  return val_one;
+}
+
+template<typename DecimalType>
+auto my_inf() -> DecimalType&
+{
+  using decimal_type = DecimalType;
+
+  static decimal_type val_inf { std::numeric_limits<decimal_type>::infinity() };
+
+  return val_inf;
+}