Add Eisel-Lemire decimal parsing fast path

NOTICE

@@ -42,6 +42,10 @@ OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
 
 File random/random.mbt is adapted from Golang's [`math/rand/v2`](https://pkg.go.dev/math/rand/v2) package.
 
+Files `internal/strconv/strconv_eisel_lemire.mbt` and
+`internal/strconv/strconv_eisel_lemire_table.mbt` are adapted from Go 1.26.2's
+`internal/strconv/atofeisel.go` and generated `internal/strconv/pow10tab.go`.
+
 License from Golang:
 Copyright 2009 The Go Authors.
 

internal/strconv/eisel_lemire_quickcheck_test.mbt

@@ -0,0 +1,79 @@
+// Copyright 2026 International Digital Economy Academy
+//
+// Licensed under the Apache License, Version 2.0 (the "License");
+// you may not use this file except in compliance with the License.
+// You may obtain a copy of the License at
+//
+//     http://www.apache.org/licenses/LICENSE-2.0
+//
+// Unless required by applicable law or agreed to in writing, software
+// distributed under the License is distributed on an "AS IS" BASIS,
+// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+// See the License for the specific language governing permissions and
+// limitations under the License.
+
+// Property tests for the Eisel-Lemire fast path.
+//
+// `try_eisel_lemire64` is free to reject any input, but every value it
+// accepts must be the correctly rounded double — bit for bit what the
+// arbitrary-precision Decimal conversion produces. The oracle reaches that
+// exact path through the public parser: appending 21 mantissa zeros (with
+// the exponent compensating) pushes the digit count past the 19 the fast
+// paths tolerate, so `parse_double` is forced onto the Decimal slow path.
+
+///|
+fn eisel_lemire_is_exact(
+  mantissa : UInt64,
+  exponent : Int,
+  negative : Bool,
+) -> Bool {
+  let fast = @strconv.try_eisel_lemire64(
+    mantissa,
+    exponent.to_int64(),
+    negative,
+  )
+  if fast.is_nan() {
+    // Rejection is always allowed: those inputs stay on the exact fallback.
+    return true
+  }
+  let sign = if negative { "-" } else { "" }
+  let padded = "\{sign}\{mantissa}\{repeat_char('0', 21)}e\{exponent - 21}"
+  parsed(padded) is Some(exact) && same_double(fast, exact)
+}
+
+///|
+/// Whatever the fast path accepts must match the exact conversion.
+test "quickcheck: accepted Eisel-Lemire values are correctly rounded" {
+  @quickcheck.check(
+    (input : (UInt64, Int, Int, Bool)) => {
+      let (raw, shift_code, exp_code, negative) = input
+      // The shift spreads mantissas across every magnitude, including the
+      // small values a uniform UInt64 almost never produces.
+      let mantissa = raw >> wrap_index(shift_code, 64)
+      // Exponents overshoot the table range (-348..=347) on both sides so
+      // the bounds check is exercised alongside the conversion.
+      let exponent = wrap_index(exp_code, 723) - 361
+      eisel_lemire_is_exact(mantissa, exponent, negative)
+    },
+    count=20000,
+  )
+}
+
+///|
+/// Mantissas next to powers of two sit on binade boundaries, where rounding
+/// carries into the next exponent and halfway cases cluster; the ambiguity
+/// rejection has to fire exactly there.
+test "quickcheck: Eisel-Lemire stays exact near binade boundaries" {
+  @quickcheck.check(
+    (input : (Int, Int, Int, Bool)) => {
+      let (bit_code, delta_code, exp_code, negative) = input
+      let base = 1UL << wrap_index(bit_code, 64)
+      let delta = (wrap_index(delta_code, 9) - 4).to_int64()
+      // Wrapping addition is fine: any UInt64 is a valid mantissa.
+      let mantissa = base + delta.reinterpret_as_uint64()
+      let exponent = wrap_index(exp_code, 723) - 361
+      eisel_lemire_is_exact(mantissa, exponent, negative)
+    },
+    count=20000,
+  )
+}

internal/strconv/parse_double_bench_test.mbt

@@ -27,6 +27,12 @@ let parse_double_underscore_bench_inputs : FixedArray[String] = [
   "123_456_789_012_345e-2", "876_543_210_987_654e-3", "1_234_567_890_123e+2", "7_654_321_098_765e-1",
 ]
 
+///|
+let parse_double_long_mantissa_bench_inputs : FixedArray[String] = [
+  "-65.613616999999977", "43.420273000000009", "-65.619720000000029", "43.418052999999986",
+  "-65.625000000000000", "43.412101000000000", "-65.630279999999994", "43.406101000000010",
+]
+
 ///|
 fn parse_double_bench_sum(inputs : FixedArray[String]) -> Double {
   let mut sum = 0.0
@@ -51,3 +57,10 @@ test "bench parse_double underscores n=4096" (it : @bench.T) {
     it.keep(parse_double_bench_sum(parse_double_underscore_bench_inputs))
   })
 }
+
+///|
+test "bench parse_double long mantissa n=4096" (it : @bench.T) {
+  it.bench(fn() {
+    it.keep(parse_double_bench_sum(parse_double_long_mantissa_bench_inputs))
+  })
+}

internal/strconv/strconv_double.mbt

@@ -79,7 +79,18 @@ pub fn parse_double(str : StringView) -> Double raise {
       // Clinger's fast path (How to read floating point numbers accurately)[https://doi.org/10.1145/989393.989430]
       match num.try_fast_path() {
         Some(value) => value
-        None => parse_decimal_priv(str).to_double_priv() // fallback to slow path
+        None => {
+          let fast = if num.many_digits {
+            @double.not_a_number
+          } else {
+            try_eisel_lemire64(num.mantissa, num.exponent, num.negative)
+          }
+          if fast.is_nan() {
+            parse_decimal_priv(str).to_double_priv() // fallback to slow path
+          } else {
+            fast
+          }
+        }
       }
   }
 }

internal/strconv/strconv_eisel_lemire.mbt

@@ -0,0 +1,128 @@
+// Copyright 2026 International Digital Economy Academy
+//
+// Licensed under the Apache License, Version 2.0 (the "License");
+// you may not use this file except in compliance with the License.
+// You may obtain a copy of the License at
+//
+//     http://www.apache.org/licenses/LICENSE-2.0
+//
+// Unless required by applicable law or agreed to in writing, software
+// distributed under the License is distributed on an "AS IS" BASIS,
+// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+// See the License for the specific language governing permissions and
+// limitations under the License.
+
+///|
+#valtype
+priv struct EiselProduct {
+  lo : UInt64
+  hi : UInt64
+}
+
+///|
+fn eisel_umul128(a : UInt64, b : UInt64) -> EiselProduct {
+  let a_lo = a & 0xffffffffUL
+  let a_hi = a >> 32
+  let b_lo = b & 0xffffffffUL
+  let b_hi = b >> 32
+  let x = a_lo * b_lo
+  let y = a_hi * b_lo + (x >> 32)
+  let z = a_lo * b_hi + (y & 0xffffffffUL)
+  let hi = a_hi * b_hi + (y >> 32) + (z >> 32)
+  { lo: a * b, hi }
+}
+
+///|
+fn eisel_mul_log2_10(exponent : Int) -> Int {
+  // floor(exponent * log2(10)) for -500 <= exponent <= 500.
+  (exponent * 108853) >> 15
+}
+
+///|
+/// Attempts Eisel-Lemire conversion of `mantissa * 10^exponent`.
+///
+/// A NaN result is a private failure sentinel. The algorithm deliberately
+/// rejects values whose correct rounding cannot be certified; callers must
+/// retain an exact Decimal fallback for those inputs.
+#doc(hidden)
+pub fn try_eisel_lemire64(
+  mantissa : UInt64,
+  exponent : Int64,
+  negative : Bool,
+) -> Double {
+  if mantissa == 0UL {
+    return if negative {
+      0x8000000000000000UL.reinterpret_as_double()
+    } else {
+      0.0
+    }
+  }
+  if exponent < EISEL_LEMIRE_POW10_MIN.to_int64() ||
+    exponent > EISEL_LEMIRE_POW10_MAX.to_int64() {
+    return @double.not_a_number
+  }
+  let exponent = exponent.to_int()
+  let table_index = (exponent - EISEL_LEMIRE_POW10_MIN) * 2
+  let pow_hi = eisel_lemire_pow10_table[table_index]
+  let pow_lo = eisel_lemire_pow10_table[table_index + 1]
+  let pow_exp2 = 1 + eisel_mul_log2_10(exponent)
+
+  // Normalize the decimal mantissa so its most significant bit is set.
+  let leading_zeros = mantissa.clz()
+  let normalized = mantissa << leading_zeros
+  let mut result_exp2 = pow_exp2 + 63 + 1023 - leading_zeros
+
+  let product = eisel_umul128(normalized, pow_hi)
+  let mut product_hi = product.hi
+  let mut product_lo = product.lo
+
+  // Use the low limb of the cached power when the first product does not
+  // contain enough information to determine the rounded result.
+  if (product_hi & 0x1ffUL) == 0x1ffUL && product_lo + normalized < normalized {
+    let wider = eisel_umul128(normalized, pow_lo)
+    let mut merged_hi = product_hi
+    let merged_lo = product_lo + wider.hi
+    if merged_lo < product_lo {
+      merged_hi += 1UL
+    }
+    if (merged_hi & 0x1ffUL) == 0x1ffUL &&
+      merged_lo + 1UL == 0UL &&
+      wider.lo + normalized < normalized {
+      return @double.not_a_number
+    }
+    product_hi = merged_hi
+    product_lo = merged_lo
+  }
+
+  // Keep 54 significant bits, then round down to the binary64 precision.
+  let top_bit = (product_hi >> 63).to_int()
+  let mut result_mantissa = product_hi >> (top_bit + 9)
+  result_exp2 -= 1 - top_bit
+
+  // An exact halfway case needs the Decimal fallback to resolve ties safely.
+  if product_lo == 0UL &&
+    (product_hi & 0x1ffUL) == 0UL &&
+    (result_mantissa & 3UL) == 1UL {
+    return @double.not_a_number
+  }
+
+  result_mantissa += result_mantissa & 1UL
+  result_mantissa = result_mantissa >> 1
+  if result_mantissa >> 53 > 0UL {
+    result_mantissa = result_mantissa >> 1
+    result_exp2 += 1
+  }
+
+  // Subnormal, overflow, and special-value boundaries remain on the exact
+  // fallback path.
+  if result_exp2 <= 0 || result_exp2 >= 0x7ff {
+    return @double.not_a_number
+  }
+  let exponent_bits = UInt64::extend_uint(result_exp2.reinterpret_as_uint()) <<
+    52
+  let mut result_bits = exponent_bits | (result_mantissa & 0x000fffffffffffffUL)
+  if negative {
+    result_bits = result_bits | 0x8000000000000000UL
+  }
+  result_bits.reinterpret_as_double()
+}