From 468c243ca81b9eae579af1bd8eadd27ce8d543a5 Mon Sep 17 00:00:00 2001
From: Victor Zverovich <victor.zverovich@gmail.com>
Date: Sun, 27 May 2018 10:57:26 -0700
Subject: [PATCH] Add a function to get cached power of 10

---
 include/fmt/format-inl.h | 13 +++++++++++++
 include/fmt/format.h     | 12 +++++-------
 test/util-test.cc        | 12 ++++++++++++
 3 files changed, 30 insertions(+), 7 deletions(-)

diff --git a/include/fmt/format-inl.h b/include/fmt/format-inl.h
index b8f9a2b0..00c97c3f 100644
--- a/include/fmt/format-inl.h
+++ b/include/fmt/format-inl.h
@@ -342,6 +342,19 @@ FMT_FUNC fp operator*(fp x, fp y) {
   uint64_t mid = (bd >> 32) + (ad & mask) + (bc & mask) + (1U << 31);
   return fp(ac + (ad >> 32) + (bc >> 32) + (mid >> 32), x.e + y.e + 64);
 }
+
+FMT_FUNC fp get_cached_power(int min_exponent, int &pow10_exponent) {
+  const double one_over_log2_10 = 0.30102999566398114;  // 1 / log2(10)
+  int index = static_cast<int>(std::ceil(
+        (min_exponent + fp::fp_significand_size - 1) * one_over_log2_10));
+  // Decimal exponent of the first (smallest) cached power of 10.
+  const int first_dec_exp = -348;
+  // Difference between two consecutive decimal exponents in cached powers of 10.
+  const int dec_exp_step = 8;
+  index = (index - first_dec_exp - 1) / dec_exp_step + 1;
+  pow10_exponent = first_dec_exp + index * dec_exp_step;
+  return fp(data::POW10_SIGNIFICANDS[index], data::POW10_EXPONENTS[index]);
+}
 }  // namespace internal
 
 #if FMT_USE_WINDOWS_H
diff --git a/include/fmt/format.h b/include/fmt/format.h
index 590ce83d..a9111d4d 100644
--- a/include/fmt/format.h
+++ b/include/fmt/format.h
@@ -272,6 +272,8 @@ class fp {
   uint64_t f;
   int e;
 
+  static constexpr int fp_significand_size = sizeof(f) * char_size;
+
   fp(uint64_t f, int e): f(f), e(e) {}
 
   // Constructs fp from an IEEE754 double. It is a template to prevent compile
@@ -305,7 +307,6 @@ class fp {
       f <<= 1;
       --e;
     }
-    auto fp_significand_size = sizeof(f) * char_size;
     // Subtract 1 to account for hidden bit.
     auto offset = fp_significand_size - double_significand_size - SHIFT - 1;
     f <<= offset;
@@ -323,12 +324,9 @@ inline fp operator-(fp x, fp y) {
 // with half-up tie breaking, r.e = x.e + y.e + 32. Result may not be normalized.
 fp operator*(fp x, fp y);
 
-// Compute k such that its cached power c_k = c_k.f * pow(2, c_k.e) satisfies
-// min_exponent <= c_k.e + e <= min_exponent + 3.
-inline int compute_cached_power_index(int e, int min_exponent) {
-  const double one_over_log2_10 = 0.30102999566398114;  // 1 / log2(10)
-  return static_cast<int>(std::ceil((min_exponent - e + 63) * one_over_log2_10));
-}
+// Returns cached power (of 10) c_k = c_k.f * pow(2, c_k.e) such that its
+// (binary) exponent satisfies min_exponent <= c_k.e <= min_exponent + 3.
+fp get_cached_power(int min_exponent, int &pow10_exponent);
 
 template <typename Allocator>
 typename Allocator::value_type *allocate(Allocator& alloc, std::size_t n) {
diff --git a/test/util-test.cc b/test/util-test.cc
index 62fba708..101f3e01 100644
--- a/test/util-test.cc
+++ b/test/util-test.cc
@@ -918,3 +918,15 @@ TEST(FPTest, Multiply) {
   EXPECT_EQ(v.f, (123 * 567 + 1u) / 2);
   EXPECT_EQ(v.e, 4 + 8 + 64);
 }
+
+TEST(FPTest, GetCachedPower) {
+  typedef std::numeric_limits<double> limits;
+  for (auto exp = limits::min_exponent; exp <= limits::max_exponent; ++exp) {
+    int dec_exp = 0;
+    auto fp = fmt::internal::get_cached_power(exp, dec_exp);
+    EXPECT_LE(exp, fp.e);
+    int dec_exp_step = 8;
+    EXPECT_LE(fp.e, exp + dec_exp_step * log2(10));
+    EXPECT_EQ(pow(10, dec_exp), ldexp(fp.f, fp.e));
+  }
+}