Implement correct bigint to double conversion.

runtime/vm/bigint_operations.cc

Index: runtime/vm/bigint_operations.cc
diff --git a/runtime/vm/bigint_operations.cc b/runtime/vm/bigint_operations.cc
index d46e78513244b24b6067289a067e70dde5a3a93a..1c1b94e27172b7d9a4ada9159e6b75a9c84381ec 100644
--- a/runtime/vm/bigint_operations.cc
+++ b/runtime/vm/bigint_operations.cc
@@ -398,14 +398,132 @@ RawSmi* BigintOperations::ToSmi(const Bigint& bigint) {
 
 
 RawDouble* BigintOperations::ToDouble(const Bigint& bigint) {
-  // TODO(floitsch/benl): This is a quick and dirty implementation to unblock
-  // other areas of the code. It does not handle all bit-twiddling correctly.
-  const double shift_value = (1 << kDigitBitSize);
-  double value = 0.0;
-  for (int i = bigint.Length() - 1; i >= 0; i--) {
-    value *= shift_value;
-    value += static_cast<double>(bigint.GetChunkAt(i));
+  if (bigint.IsZero()) {

[cshapiro, 2012/03/07 21:49:18]

This fast path should be expanded to covert values that are in the signed 64-bit
integer range as it is possible to have bigint objects whose values in the mint
and smi range.

[floitsch, 2012/03/20 04:06:48]

Done.

+    return Double::New(0.0);
+  }
+

[cshapiro, 2012/03/07 21:49:18]

How about an early exit here if the value is too large to fit into a double?  We
can quickly reject a large class of out-of-range values up front.

[floitsch, 2012/03/20 04:06:48]

Done.

+  static const int kPhysicalSignificandSize = 52;
+  // The significand size has an additional hidden bit.
+  static const int kSignificandSize = kPhysicalSignificandSize + 1;
+  static const int kExponentBias = 0x3FF + kPhysicalSignificandSize;
+  static const int kMaxExponent = 0x7FF - kExponentBias;
+  static const uint64_t kOne64 = 1;
+
+  // A double is composed of an exponent e and a significand s. Its value equals
+  // s * 2^e. The significand has 53 bits of which the first one must always be
+  // 1 (at least for then numbers we are working with here) and is therefore
+  // omitted. The physical size of the significand is thus 52 bits.
+  // The exponent has 11 bits and is biased by 0x3FF + 52. For example an
+  // exponent e = 10 is written as 0x3FF + 52 + 10 (in the 11 bits that are
+  // reserved for the exponent).
+  // When converting the given bignum to a double we have to pay attention to
+  // the rounding. In particular we have to decide which double to pick if an
+  // input lies exactly between two doubles. As usual with double operations
+  // we pick the double with an even significand in such cases.
+  //
+  // General approach of this algorithm: Get 54 bits (one more than the
+  // significand size) of the bigint. If the last bit is then 1, then (without
+  // knowledge of the remaining bits) we could have a half-way number.
+  // If the second-to-last bit is odd then we know that we have to round up:
+  // if the remaining bits are not zero then the input lies closer to the higher
+  // double. If the remaining bits are zero then we have a half-way case and
+  // we need to round up too (rounding to the even double).
+  // If the second-to-last bit is even then we need to look at the remaining
+  // bits to determine if any of them is not zero. If that's the case then the
+  // number lies closer to the next-higher double. Otherwise we round the
+  // half-way case down to even.
+
+  uint64_t significand_bits = 0;
+  bool remaining_bits_are_zero = true;
+  intptr_t exponent = 0;
+

[cshapiro, 2012/03/07 21:49:18]

Why not use the ShiftRight method to extract the upper 54 bits?

[floitsch, 2012/03/20 04:06:48]

We could, but it seems overkill. I hope the new version is simple enough.

+  intptr_t digit_index = bigint.Length() - 1;
+  // In order to round correctly we need to look at half-way cases. Therefore we
+  // get kSignificandSize + 1 bits. If the last bit is 1 then we have to look
+  // at the remaining bits to know if we have to round up.
+  int needed_bits = kSignificandSize + 1;
+  while ((digit_index) >= 0 && (needed_bits > 0)) {
+    Chunk digit = bigint.GetChunkAt(digit_index--);
+    if (significand_bits == 0) {
+      significand_bits = digit;
+      while (digit != 0) {
+        needed_bits--;
+        digit >>= 1;
+      }
+    } else if (needed_bits >= kDigitBitSize) {
+      significand_bits <<= kDigitBitSize;
+      significand_bits |= digit;
+      needed_bits -= kDigitBitSize;
+    } else {
+      significand_bits <<= needed_bits;
+      int discarded_bits_count = kDigitBitSize - needed_bits;
+      significand_bits |= digit >> discarded_bits_count;
+      uint64_t discarded_bits_mask = (kOne64 << discarded_bits_count) - 1;
+      remaining_bits_are_zero = ((digit & discarded_bits_mask) == 0);
+      exponent = discarded_bits_count;
+      needed_bits = 0;
+    }
+  }
+  ASSERT(significand_bits != 0);
+  ASSERT(needed_bits > 0 || ((significand_bits >> kSignificandSize) == 1));
+
+  exponent += (digit_index + 1) * kDigitBitSize;
+
+  if (needed_bits > 0) {
+    // The value cannot be a half-way case. Just shift the significand to the
+    // correct position and update the exponent accordingly.
+    significand_bits <<= (needed_bits - 1);
+    exponent -= (needed_bits - 1);
+  } else if ((significand_bits & 1) == 1) {
+    bool round_up = false;
+
+    if ((significand_bits & 2) != 0 || !remaining_bits_are_zero) {
+      // Even if the remaining bits are zero we still need to round up since we
+      // want to round to even for half-way cases.
+      round_up = true;
+    } else {
+      // Could be a half-way case. See if the remaining bits are non-zero.
+      for (intptr_t i = 0; i <= digit_index; i++) {
+        if (!remaining_bits_are_zero) break;
+        remaining_bits_are_zero = (bigint.GetChunkAt(i) == 0);
+      }
+      round_up = !remaining_bits_are_zero;
+    }
+    // Remove the bit we used to determine half-way cases.
+    significand_bits >>= 1;
+    exponent++;
+    if (round_up) {
+      significand_bits += 1;
+      // It might be that we just went from 53 bits to 54 bits.
+      // Example: After adding 1 to 1FFF..FF (with 53 bits set to 1) we have
+      // 2000..00 (= 2 ^ 54). In this case we adjust the exponent and shift to
+      // the right yielding again an 53 bit significand: 1000...00 ('1' followed
+      // by 52 '0' digits).
+      if (significand_bits == (kOne64 << kSignificandSize)) {
+        significand_bits >>= 1;
+        exponent++;
+      }
+    }
+  } else {
+    // Simply round down.
+    significand_bits >>= 1;
+    exponent++;
   }
+
+  if (exponent >= kMaxExponent) {
+    // Infinity.
+    exponent = kMaxExponent;
+    significand_bits = 0;
+  } else {
+    // Clear the hidden bit in the significand.
+    ASSERT((significand_bits >> (kSignificandSize - 1)) == 1);
+    significand_bits &= ((kOne64 << kPhysicalSignificandSize) - 1);
+  }
+  uint64_t biased_exponent = exponent + kExponentBias;
+  uint64_t double_bits =
+      (biased_exponent << kPhysicalSignificandSize) | significand_bits;
+
+  double value = bit_cast<double>(double_bits);
   if (bigint.IsNegative()) {
     value = -value;

[cshapiro, 2012/03/07 21:49:18]

Why not put the sign bit in when you are constructing "value" above?

[floitsch, 2012/03/20 04:06:48]

It seems cleaner not to duplicate the bit_cast expression. It wouldn't fit on
one line anyways.

   }