internal/strconv: extract fixed-precision ftoa from ftoaryu.go

author Russ Cox <rsc@golang.org>

Sat, 1 Nov 2025 13:41:40 +0000 (09:41 -0400)

committer Gopher Robot <gobot@golang.org>

Tue, 4 Nov 2025 04:29:54 +0000 (20:29 -0800)
author Russ Cox <rsc@golang.org>
Sat, 1 Nov 2025 13:41:40 +0000 (09:41 -0400)
committer Gopher Robot <gobot@golang.org>
Tue, 4 Nov 2025 04:29:54 +0000 (20:29 -0800)
diff --git a/src/internal/strconv/atoi.go b/src/internal/strconv/atoi.go

index 5bc259e7e55e83d8906cb7790b23397fa099a46a..4bbcb4f5da7aeb88c2fb75000d8e4c9231633091 100644 (file)
--- a/src/internal/strconv/atoi.go
+++ b/src/internal/strconv/atoi.go
@@ -41,8 +41,6 @@ const intSize = 32 << (^uint(0) >> 63)
  // IntSize is the size in bits of an int or uint value.
  const IntSize = intSize
  
-const maxUint64 = 1<<64 - 1
-
  // ParseUint is like [ParseInt] but for unsigned numbers.
  //
  // A sign prefix is not permitted.
diff --git a/src/internal/strconv/export_test.go b/src/internal/strconv/export_test.go

index 86435f66cf84c646da1ab76b6cd6fc263f1484a6..c879f24480a450bb8d78ee729fcaebf81ce775da 100644 (file)
--- a/src/internal/strconv/export_test.go
+++ b/src/internal/strconv/export_test.go
@@ -18,6 +18,9 @@ var (
         Pow10            = pow10
         Umul128          = umul128
         Umul192          = umul192
+       Div5Tab          = div5Tab
+       DivisiblePow5    = divisiblePow5
+       TrimZeros        = trimZeros
  )
  
  func NewDecimal(i uint64) *decimal {
diff --git a/src/internal/strconv/ftoa.go b/src/internal/strconv/ftoa.go

index 1aec5447ece8b8a789985cee762a4a74f57bae1d..fd30f28289ae9d4efc6bbede5f2966f4477f0f0b 100644 (file)
--- a/src/internal/strconv/ftoa.go
+++ b/src/internal/strconv/ftoa.go
@@ -123,16 +123,17 @@ func genericFtoa(dst []byte, val float64, fmt byte, prec, bitSize int) []byte {
                 return bigFtoa(dst, prec, fmt, neg, mant, exp, flt)
         }
  
-       var digs decimalSlice
-       ok := false
         // Negative precision means "only as much as needed to be exact."
         shortest := prec < 0
+       var digs decimalSlice
+       if mant == 0 {
+               return formatDigits(dst, shortest, neg, digs, prec, fmt)
+       }
         if shortest {
                 // Use Ryu algorithm.
                 var buf [32]byte
                 digs.d = buf[:]
                 ryuFtoaShortest(&digs, mant, exp-int(flt.mantbits), flt)
-               ok = true
                 // Precision for shortest representation mode.
                 switch fmt {
                 case 'e', 'E':
@@ -142,7 +143,11 @@ func genericFtoa(dst []byte, val float64, fmt byte, prec, bitSize int) []byte {
                 case 'g', 'G':
                         prec = digs.nd
                 }
-       } else if fmt != 'f' {
+               return formatDigits(dst, shortest, neg, digs, prec, fmt)
+       }
+
+       // TODO figure out when we can use fast code for f
+       if fmt != 'f' {
                 // Fixed number of digits.
                 digits := prec
                 switch fmt {
@@ -157,21 +162,15 @@ func genericFtoa(dst []byte, val float64, fmt byte, prec, bitSize int) []byte {
                         // Invalid mode.
                         digits = 1
                 }
-               var buf [24]byte
-               if bitSize == 32 && digits <= 9 {
+               if digits <= 18 {
+                       var buf [24]byte
                         digs.d = buf[:]
-                       ryuFtoaFixed32(&digs, uint32(mant), exp-int(flt.mantbits), digits)
-                       ok = true
-               } else if digits <= 18 {
-                       digs.d = buf[:]
-                       ryuFtoaFixed64(&digs, mant, exp-int(flt.mantbits), digits)
-                       ok = true
+                       fixedFtoa(&digs, mant, exp-int(flt.mantbits), digits)
+                       return formatDigits(dst, false, neg, digs, prec, fmt)
                 }
         }
-       if !ok {
-               return bigFtoa(dst, prec, fmt, neg, mant, exp, flt)
-       }
-       return formatDigits(dst, shortest, neg, digs, prec, fmt)
+
+       return bigFtoa(dst, prec, fmt, neg, mant, exp, flt)
  }
  
  // bigFtoa uses multiprecision computations to format a float.
diff --git a/src/internal/strconv/ftoa_test.go b/src/internal/strconv/ftoa_test.go

index d510629537548bbae7ff529d9e7fcdc2d8dcf09c..4e6f4629288c35276f14fc33fd13e8237ddcb063 100644 (file)
--- a/src/internal/strconv/ftoa_test.go
+++ b/src/internal/strconv/ftoa_test.go
@@ -177,6 +177,16 @@ var ftoatests = []ftoaTest{
         {1.801439850948199e+16, 'g', -1, "1.801439850948199e+16"},
         {5.960464477539063e-08, 'g', -1, "5.960464477539063e-08"},
         {1.012e-320, 'g', -1, "1.012e-320"},
+
+       // Cases from TestFtoaRandom that caught bugs in fixedFtoa.
+       {8177880169308380. * (1 << 1), 'e', 14, "1.63557603386168e+16"},
+       {8393378656576888. * (1 << 1), 'e', 15, "1.678675731315378e+16"},
+       {8738676561280626. * (1 << 4), 'e', 16, "1.3981882498049002e+17"},
+       {8291032395191335. / (1 << 30), 'e', 5, "7.72163e+06"},
+
+       // Exercise divisiblePow5 case in fixedFtoa
+       {2384185791015625. * (1 << 12), 'e', 5, "9.76562e+18"},
+       {2384185791015625. * (1 << 13), 'e', 5, "1.95312e+19"},
  }
  
  func TestFtoa(t *testing.T) {
@@ -253,7 +263,7 @@ func TestFtoaRandom(t *testing.T) {
                 shortSlow = FormatFloat(x, 'e', prec, 64)
                 SetOptimize(true)
                 if shortSlow != shortFast {
-                       t.Errorf("%b printed as %s, want %s", x, shortFast, shortSlow)
+                       t.Errorf("%b printed with %%.%de as %s, want %s", x, prec, shortFast, shortSlow)
                 }
         }
  }
diff --git a/src/internal/strconv/ftoafixed.go b/src/internal/strconv/ftoafixed.go

new file mode 100644 (file)

index 0000000..f3542d1
--- /dev/null
+++ b/src/internal/strconv/ftoafixed.go
@@ -0,0 +1,156 @@
+// Copyright 2025 The Go Authors. All rights reserved.
+// Use of this source code is governed by a BSD-style
+// license that can be found in the LICENSE file.
+
+package strconv
+
+import "math/bits"
+
+var uint64pow10 = [...]uint64{
+       1, 1e1, 1e2, 1e3, 1e4, 1e5, 1e6, 1e7, 1e8, 1e9,
+       1e10, 1e11, 1e12, 1e13, 1e14, 1e15, 1e16, 1e17, 1e18, 1e19,
+}
+
+// fixedFtoa formats a number of decimal digits of mant*(2^exp) into d,
+// where mant > 0 and 1 ≤ digits ≤ 18.
+func fixedFtoa(d *decimalSlice, mant uint64, exp, digits int) {
+       // The strategy here is to multiply (mant * 2^exp) by a power of 10
+       // to make the resulting integer be the number of digits we want.
+       //
+       // Adams proved in the Ryu paper that 128-bit precision in the
+       // power-of-10 constant is sufficient to produce correctly
+       // rounded output for all float64s, up to 18 digits.
+       // https://dl.acm.org/doi/10.1145/3192366.3192369
+       //
+       // TODO(rsc): The paper is not focused on, nor terribly clear about,
+       // this fact in this context, and the proof seems too complicated.
+       // Post a shorter, more direct proof and link to it here.
+
+       if digits > 18 {
+               panic("fixedFtoa called with digits > 18")
+       }
+
+       // Shift mantissa to have 64 bits,
+       // so that the 192-bit product below will
+       // have at least 63 bits in its top word.
+       b := 64 - bits.Len64(mant)
+       mant <<= b
+       exp -= b
+
+       // We have f = mant * 2^exp ≥ 2^(63+exp)
+       // and we want to multiply it by some 10^p
+       // to make it have the number of digits plus one rounding bit:
+       //
+       //      2 * 10^(digits-1) ≤ f * 10^p < ~2 * 10^digits
+       //
+       // The lower bound is required, but the upper bound is approximate:
+       // we must not have too few digits, but we can round away extra ones.
+       //
+       //      f * 10^p ≥ 2 * 10^(digits-1)
+       //      10^p ≥ 2 * 10^(digits-1) / f                         [dividing by f]
+       //      p ≥ (log₁₀ 2) + (digits-1) - log₁₀ f                 [taking log₁₀]
+       //      p ≥ (log₁₀ 2) + (digits-1) - log₁₀ (mant * 2^exp)    [expanding f]
+       //      p ≥ (log₁₀ 2) + (digits-1) - (log₁₀ 2) * (64 + exp)  [mant < 2⁶⁴]
+       //      p ≥ (digits - 1) - (log₁₀ 2) * (63 + exp)            [refactoring]
+       //
+       // Once we have p, we can compute the scaled value:
+       //
+       //      dm * 2^de = mant * 2^exp * 10^p
+       //                = mant * 2^exp * pow/2^128 * 2^exp2.
+       //                = (mant * pow/2^128) * 2^(exp+exp2).
+       p := (digits - 1) - mulLog10_2(63+exp)
+       pow, exp2, ok := pow10(p)
+       if !ok {
+               // This never happens due to the range of float32/float64 exponent
+               panic("fixedFtoa: pow10 out of range")
+       }
+       if -22 <= p && p < 0 {
+               // Special case: Let q=-p. q is in [1,22]. We are dividing by 10^q
+               // and the mantissa may be a multiple of 5^q (5^22 < 2^53),
+               // in which case the division must be computed exactly and
+               // recorded as exact for correct rounding. Our normal computation is:
+               //
+               //      dm = floor(mant * floor(10^p * 2^s))
+               //
+               // for some scaling shift s. To make this an exact division,
+               // it suffices to change the inner floor to a ceil:
+               //
+               //      dm = floor(mant * ceil(10^p * 2^s))
+               //
+               // In the range of values we are using, the floor and ceil
+               // cancel each other out and the high 64 bits of the product
+               // come out exactly right.
+               // (This is the same trick compilers use for division by constants.
+               // See Hacker's Delight, 2nd ed., Chapter 10.)
+               pow.Lo++
+       }
+       dm, lo1, lo0 := umul192(mant, pow)
+       de := exp + exp2
+
+       // Check whether any bits have been truncated from dm.
+       // If so, set dt != 0. If not, leave dt == 0 (meaning dm is exact).
+       var dt uint
+       switch {
+       default:
+               // Most powers of 10 use a truncated constant,
+               // meaning the result is also truncated.
+               dt = 1
+       case 0 <= p && p <= 55:
+               // Small positive powers of 10 (up to 10⁵⁵) can be represented
+               // precisely in a 128-bit mantissa (5⁵⁵ ≤ 2¹²⁸), so the only truncation
+               // comes from discarding the low bits of the 192-bit product.
+               //
+               // TODO(rsc): The new proof mentioned above should also
+               // prove that we can't have lo1 == 0 and lo0 != 0.
+               // After proving that, drop computation and use of lo0 here.
+               dt = bool2uint(lo1|lo0 != 0)
+       case -22 <= p && p < 0 && divisiblePow5(mant, -p):
+               // If the original mantissa was a multiple of 5^p,
+               // the result is exact. (See comment above for pow.Lo++.)
+               dt = 0
+       }
+
+       // The value we want to format is dm * 2^de, where de < 0.
+       // Multply by 2^de by shifting, but leave one extra bit for rounding.
+       // After the shift, the "integer part" of dm is dm>>1,
+       // the "rounding bit" (the first fractional bit) is dm&1,
+       // and the "truncated bit" (have any bits been discarded?) is dt.
+       shift := -de - 1
+       dt |= bool2uint(dm&(1<<shift-1) != 0)
+       dm >>= shift
+
+       // Set decimal point in eventual formatted digits,
+       // so we can update it as we adjust the digits.
+       d.dp = digits - p
+
+       // Trim excess digit if any, updating truncation and decimal point.
+       // The << 1 is leaving room for the rounding bit.
+       max := uint64pow10[digits] << 1
+       if dm >= max {
+               var r uint
+               dm, r = dm/10, uint(dm%10)
+               dt |= bool2uint(r != 0)
+               d.dp++
+       }
+
+       // Round and shift away rounding bit.
+       // We want to round up when
+       // (a) the fractional part is > 0.5 (dm&1 != 0 and dt == 1)
+       // (b) or the fractional part is ≥ 0.5 and the integer part is odd
+       //     (dm&1 != 0 and dm&2 != 0).
+       // The bitwise expression encodes that logic.
+       dm += uint64(uint(dm) & (dt | uint(dm)>>1) & 1)
+       dm >>= 1
+       if dm == max>>1 {
+               // 999... rolled over to 1000...
+               dm = uint64pow10[digits-1]
+               d.dp++
+       }
+
+       // Format digits into d.
+       formatBase10(d.d[:digits], dm)
+       d.nd = digits
+       for d.d[d.nd-1] == '0' {
+               d.nd--
+       }
+}
diff --git a/src/internal/strconv/ftoaryu.go b/src/internal/strconv/ftoaryu.go

index 999af515029204041081db6fce2ab8c56514e801..9407bfec445680fddb237ea1a6d7866270ed1205 100644 (file)
--- a/src/internal/strconv/ftoaryu.go
+++ b/src/internal/strconv/ftoaryu.go
@@ -4,203 +4,11 @@
  
  package strconv
  
-import (
-       "math/bits"
-)
+import "math/bits"
  
  // binary to decimal conversion using the Ryū algorithm.
  //
  // See Ulf Adams, "Ryū: Fast Float-to-String Conversion" (doi:10.1145/3192366.3192369)
-//
-// Fixed precision formatting is a variant of the original paper's
-// algorithm, where a single multiplication by 10^k is required,
-// sharing the same rounding guarantees.
-
-// ryuFtoaFixed32 formats mant*(2^exp) with prec decimal digits.
-func ryuFtoaFixed32(d *decimalSlice, mant uint32, exp int, prec int) {
-       if prec < 0 {
-               panic("ryuFtoaFixed32 called with negative prec")
-       }
-       if prec > 9 {
-               panic("ryuFtoaFixed32 called with prec > 9")
-       }
-       // Zero input.
-       if mant == 0 {
-               d.nd, d.dp = 0, 0
-               return
-       }
-       // Renormalize to a 25-bit mantissa.
-       e2 := exp
-       if b := bits.Len32(mant); b < 25 {
-               mant <<= uint(25 - b)
-               e2 += b - 25
-       }
-       // Choose an exponent such that rounded mant*(2^e2)*(10^q) has
-       // at least prec decimal digits, i.e
-       //     mant*(2^e2)*(10^q) >= 10^(prec-1)
-       // Because mant >= 2^24, it is enough to choose:
-       //     2^(e2+24) >= 10^(-q+prec-1)
-       // or q = -mulLog10_2(e2+24) + prec - 1
-       q := -mulLog10_2(e2+24) + prec - 1
-
-       // Now compute mant*(2^e2)*(10^q).
-       // Is it an exact computation?
-       // Only small positive powers of 10 are exact (5^28 has 66 bits).
-       exact := q <= 27 && q >= 0
-
-       di, dexp2, d0 := mult64bitPow10(mant, e2, q)
-       if dexp2 >= 0 {
-               panic("not enough significant bits after mult64bitPow10")
-       }
-       // As a special case, computation might still be exact, if exponent
-       // was negative and if it amounts to computing an exact division.
-       // In that case, we ignore all lower bits.
-       // Note that division by 10^11 cannot be exact as 5^11 has 26 bits.
-       if q < 0 && q >= -10 && divisibleByPower5(uint64(mant), -q) {
-               exact = true
-               d0 = true
-       }
-       // Remove extra lower bits and keep rounding info.
-       extra := uint(-dexp2)
-       extraMask := uint32(1<<extra - 1)
-
-       di, dfrac := di>>extra, di&extraMask
-       roundUp := false
-       if exact {
-               // If we computed an exact product, d + 1/2
-               // should round to d+1 if 'd' is odd.
-               roundUp = dfrac > 1<<(extra-1) ||
-                       (dfrac == 1<<(extra-1) && !d0) ||
-                       (dfrac == 1<<(extra-1) && d0 && di&1 == 1)
-       } else {
-               // otherwise, d+1/2 always rounds up because
-               // we truncated below.
-               roundUp = dfrac>>(extra-1) == 1
-       }
-       if dfrac != 0 {
-               d0 = false
-       }
-       // Proceed to the requested number of digits
-       formatDecimal(d, uint64(di), !d0, roundUp, prec)
-       // Adjust exponent
-       d.dp -= q
-}
-
-// ryuFtoaFixed64 formats mant*(2^exp) with prec decimal digits.
-func ryuFtoaFixed64(d *decimalSlice, mant uint64, exp int, prec int) {
-       if prec > 18 {
-               panic("ryuFtoaFixed64 called with prec > 18")
-       }
-       // Zero input.
-       if mant == 0 {
-               d.nd, d.dp = 0, 0
-               return
-       }
-       // Renormalize to a 55-bit mantissa.
-       e2 := exp
-       if b := bits.Len64(mant); b < 55 {
-               mant = mant << uint(55-b)
-               e2 += b - 55
-       }
-       // Choose an exponent such that rounded mant*(2^e2)*(10^q) has
-       // at least prec decimal digits, i.e
-       //     mant*(2^e2)*(10^q) >= 10^(prec-1)
-       // Because mant >= 2^54, it is enough to choose:
-       //     2^(e2+54) >= 10^(-q+prec-1)
-       // or q = -mulLog10_2(e2+54) + prec - 1
-       //
-       // The minimal required exponent is -mulLog10_2(1025)+18 = -291
-       // The maximal required exponent is mulLog10_2(1074)+18 = 342
-       q := -mulLog10_2(e2+54) + prec - 1
-
-       // Now compute mant*(2^e2)*(10^q).
-       // Is it an exact computation?
-       // Only small positive powers of 10 are exact (5^55 has 128 bits).
-       exact := q <= 55 && q >= 0
-
-       di, dexp2, d0 := mult128bitPow10(mant, e2, q)
-       if dexp2 >= 0 {
-               panic("not enough significant bits after mult128bitPow10")
-       }
-       // As a special case, computation might still be exact, if exponent
-       // was negative and if it amounts to computing an exact division.
-       // In that case, we ignore all lower bits.
-       // Note that division by 10^23 cannot be exact as 5^23 has 54 bits.
-       if q < 0 && q >= -22 && divisibleByPower5(mant, -q) {
-               exact = true
-               d0 = true
-       }
-       // Remove extra lower bits and keep rounding info.
-       extra := uint(-dexp2)
-       extraMask := uint64(1<<extra - 1)
-
-       di, dfrac := di>>extra, di&extraMask
-       roundUp := false
-       if exact {
-               // If we computed an exact product, d + 1/2
-               // should round to d+1 if 'd' is odd.
-               roundUp = dfrac > 1<<(extra-1) ||
-                       (dfrac == 1<<(extra-1) && !d0) ||
-                       (dfrac == 1<<(extra-1) && d0 && di&1 == 1)
-       } else {
-               // otherwise, d+1/2 always rounds up because
-               // we truncated below.
-               roundUp = dfrac>>(extra-1) == 1
-       }
-       if dfrac != 0 {
-               d0 = false
-       }
-       // Proceed to the requested number of digits
-       formatDecimal(d, di, !d0, roundUp, prec)
-       // Adjust exponent
-       d.dp -= q
-}
-
-var uint64pow10 = [...]uint64{
-       1, 1e1, 1e2, 1e3, 1e4, 1e5, 1e6, 1e7, 1e8, 1e9,
-       1e10, 1e11, 1e12, 1e13, 1e14, 1e15, 1e16, 1e17, 1e18, 1e19,
-}
-
-// formatDecimal fills d with at most prec decimal digits
-// of mantissa m. The boolean trunc indicates whether m
-// is truncated compared to the original number being formatted.
-func formatDecimal(d *decimalSlice, m uint64, trunc bool, roundUp bool, prec int) {
-       max := uint64pow10[prec]
-       trimmed := 0
-       for m >= max {
-               a, b := m/10, m%10
-               m = a
-               trimmed++
-               if b > 5 {
-                       roundUp = true
-               } else if b < 5 {
-                       roundUp = false
-               } else { // b == 5
-                       // round up if there are trailing digits,
-                       // or if the new value of m is odd (round-to-even convention)
-                       roundUp = trunc || m&1 == 1
-               }
-               if b != 0 {
-                       trunc = true
-               }
-       }
-       if roundUp {
-               m++
-       }
-       if m >= max {
-               // Happens if di was originally 99999....xx
-               m /= 10
-               trimmed++
-       }
-       // render digits
-       formatBase10(d.d[:prec], m)
-       d.nd = prec
-       for d.d[d.nd-1] == '0' {
-               d.nd--
-               trimmed++
-       }
-       d.dp = d.nd + trimmed
-}
  
  // ryuFtoaShortest formats mant*2^exp with prec decimal digits.
  func ryuFtoaShortest(d *decimalSlice, mant uint64, exp int, flt *floatInfo) {
@@ -249,13 +57,13 @@ func ryuFtoaShortest(d *decimalSlice, mant uint64, exp int, flt *floatInfo) {
         if q < 0 && q >= -24 {
                 // Division by a power of ten may be exact.
                 // (note that 5^25 is a 59-bit number so division by 5^25 is never exact).
-               if divisibleByPower5(ml, -q) {
+               if divisiblePow5(ml, -q) {
                         dl0 = true
                 }
-               if divisibleByPower5(mc, -q) {
+               if divisiblePow5(mc, -q) {
                         dc0 = true
                 }
-               if divisibleByPower5(mu, -q) {
+               if divisiblePow5(mu, -q) {
                         du0 = true
                 }
         }
@@ -497,16 +305,3 @@ func mult128bitPow10(m uint64, e2, q int) (resM uint64, resE int, exact bool) {
         hi, mid, lo := umul192(m, pow)
         return hi<<9 | mid>>55, e2, mid<<9 == 0 && lo == 0
  }
-
-func divisibleByPower5(m uint64, k int) bool {
-       if m == 0 {
-               return true
-       }
-       for i := 0; i < k; i++ {
-               if m%5 != 0 {
-                       return false
-               }
-               m /= 5
-       }
-       return true
-}
diff --git a/src/internal/strconv/import_test.go b/src/internal/strconv/import_test.go

index ed1015ee5d214d12a051a20f51325eef5b5e4441..3dab2bf9e56e2b98e5372a69aea35529faea723a 100644 (file)
--- a/src/internal/strconv/import_test.go
+++ b/src/internal/strconv/import_test.go
@@ -20,4 +20,7 @@ var (
         pow10            = Pow10
         umul128          = Umul128
         umul192          = Umul192
+       div5Tab          = Div5Tab
+       divisiblePow5    = DivisiblePow5
+       trimZeros        = TrimZeros
  )
diff --git a/src/internal/strconv/math.go b/src/internal/strconv/math.go

index 37303d76dbbaf52a9ac37150f43c96da7b5c5869..3b884e846a62224d780f8c88d552b2fbf431460e 100644 (file)
--- a/src/internal/strconv/math.go
+++ b/src/internal/strconv/math.go
@@ -56,3 +56,124 @@ func mulLog2_10(x int) int {
         // log(10)/log(2) ≈ 3.32192809489 ≈ 108853 / 2^15
         return (x * 108853) >> 15
  }
+
+func bool2uint(b bool) uint {
+       if b {
+               return 1
+       }
+       return 0
+}
+
+// Exact Division and Remainder Checking
+//
+// An exact division x/c (exact means x%c == 0)
+// can be implemented by x*m where m is the multiplicative inverse of c (m*c == 1).
+//
+// Since c is also the multiplicative inverse of m, x*m is lossless,
+// and all the exact multiples of c map to all of [0, maxUint64/c].
+// The non-multiples are forced to map to larger values.
+// This also gives a quick test for whether x is an exact multiple of c:
+// compute the exact division and check whether it's at most maxUint64/c:
+//     x%c == 0 => x*m <= maxUint64/c.
+//
+// Only odd c have multiplicative inverses mod powers of two.
+// To do an exact divide x / (c<<s) we can use (x/c)>>s instead.
+// And to check for remainder, we need to check that those low s
+// bits are all zero before we shift them away. We can merge that
+// with the <= for the exact odd remainder check by rotating the
+// shifted bits into the high part instead:
+//     x%(c<<s) == 0 => bits.RotateLeft64(x*m, -s) <= maxUint64/c.
+//
+// The compiler does this transformation automatically in general,
+// but we apply it here by hand in a few ways that the compiler can't help with.
+//
+// For a more detailed explanation, see
+// Henry S. Warren, Jr., Hacker's Delight, 2nd ed., sections 10-16 and 10-17.
+
+// divisiblePow5 reports whether x is divisible by 5^p.
+// It returns false for p not in [1, 22],
+// because we only care about float64 mantissas, and 5^23 > 2^53.
+func divisiblePow5(x uint64, p int) bool {
+       return 1 <= p && p <= 22 && x*div5Tab[p-1][0] <= div5Tab[p-1][1]
+}
+
+const maxUint64 = 1<<64 - 1
+
+// div5Tab[p-1] is the multiplicative inverse of 5^p and maxUint64/5^p.
+var div5Tab = [22][2]uint64{
+       {0xcccccccccccccccd, maxUint64 / 5},
+       {0x8f5c28f5c28f5c29, maxUint64 / 5 / 5},
+       {0x1cac083126e978d5, maxUint64 / 5 / 5 / 5},
+       {0xd288ce703afb7e91, maxUint64 / 5 / 5 / 5 / 5},
+       {0x5d4e8fb00bcbe61d, maxUint64 / 5 / 5 / 5 / 5 / 5},
+       {0x790fb65668c26139, maxUint64 / 5 / 5 / 5 / 5 / 5 / 5},
+       {0xe5032477ae8d46a5, maxUint64 / 5 / 5 / 5 / 5 / 5 / 5 / 5},
+       {0xc767074b22e90e21, maxUint64 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5},
+       {0x8e47ce423a2e9c6d, maxUint64 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5},
+       {0x4fa7f60d3ed61f49, maxUint64 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5},
+       {0x0fee64690c913975, maxUint64 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5},
+       {0x3662e0e1cf503eb1, maxUint64 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5},
+       {0xa47a2cf9f6433fbd, maxUint64 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5},
+       {0x54186f653140a659, maxUint64 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5},
+       {0x7738164770402145, maxUint64 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5},
+       {0xe4a4d1417cd9a041, maxUint64 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5},
+       {0xc75429d9e5c5200d, maxUint64 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5},
+       {0xc1773b91fac10669, maxUint64 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5},
+       {0x26b172506559ce15, maxUint64 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5},
+       {0xd489e3a9addec2d1, maxUint64 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5},
+       {0x90e860bb892c8d5d, maxUint64 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5},
+       {0x502e79bf1b6f4f79, maxUint64 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5 / 5},
+}
+
+// trimZeros trims trailing zeros from x.
+// It finds the largest p such that x % 10^p == 0
+// and then returns x / 10^p, p.
+//
+// This is here for reference and tested, because it is an optimization
+// used by other ftoa algorithms, but in our implementations it has
+// never been benchmarked to be faster than trimming zeros after
+// formatting into decimal bytes.
+func trimZeros(x uint64) (uint64, int) {
+       const (
+               div1e8m  = 0xc767074b22e90e21
+               div1e8le = maxUint64 / 100000000
+
+               div1e4m  = 0xd288ce703afb7e91
+               div1e4le = maxUint64 / 10000
+
+               div1e2m  = 0x8f5c28f5c28f5c29
+               div1e2le = maxUint64 / 100
+
+               div1e1m  = 0xcccccccccccccccd
+               div1e1le = maxUint64 / 10
+       )
+
+       // _ = assert[x - y] asserts at compile time that x == y.
+       // Assert that the multiplicative inverses are correct
+       // by checking that (div1eNm * 5^N) % 1<<64 == 1.
+       var assert [1]struct{}
+       _ = assert[(div1e8m*5*5*5*5*5*5*5*5)%(1<<64)-1]
+       _ = assert[(div1e4m*5*5*5*5)%(1<<64)-1]
+       _ = assert[(div1e2m*5*5)%(1<<64)-1]
+       _ = assert[(div1e1m*5)%(1<<64)-1]
+
+       // Cut 8 zeros, then 4, then 2, then 1.
+       p := 0
+       for d := bits.RotateLeft64(x*div1e8m, -8); d <= div1e8le; d = bits.RotateLeft64(x*div1e8m, -8) {
+               x = d
+               p += 8
+       }
+       if d := bits.RotateLeft64(x*div1e4m, -4); d <= div1e4le {
+               x = d
+               p += 4
+       }
+       if d := bits.RotateLeft64(x*div1e2m, -2); d <= div1e2le {
+               x = d
+               p += 2
+       }
+       if d := bits.RotateLeft64(x*div1e1m, -1); d <= div1e1le {
+               x = d
+               p += 1
+       }
+       return x, p
+}
diff --git a/src/internal/strconv/math_test.go b/src/internal/strconv/math_test.go

index 3a1ff3400c06bbba982f42755ad38dd881338b85..55e25f98cfee284141057a15fe6073e4131460db 100644 (file)
--- a/src/internal/strconv/math_test.go
+++ b/src/internal/strconv/math_test.go
@@ -93,3 +93,73 @@ func TestMulLog2_10(t *testing.T) {
                 }
         }
  }
+
+func pow5(p int) uint64 {
+       x := uint64(1)
+       for range p {
+               x *= 5
+       }
+       return x
+}
+
+func TestDivisiblePow5(t *testing.T) {
+       for p := 1; p <= 22; p++ {
+               x := pow5(p)
+               if divisiblePow5(1, p) {
+                       t.Errorf("divisiblePow5(1, %d) = true, want, false", p)
+               }
+               if divisiblePow5(x-1, p) {
+                       t.Errorf("divisiblePow5(%d, %d) = true, want false", x-1, p)
+               }
+               if divisiblePow5(x+1, p) {
+                       t.Errorf("divisiblePow5(%d, %d) = true, want false", x-1, p)
+               }
+               if divisiblePow5(x/5, p) {
+                       t.Errorf("divisiblePow5(%d, %d) = true, want false", x/5, p)
+               }
+               if !divisiblePow5(0, p) {
+                       t.Errorf("divisiblePow5(0, %d) = false, want true", p)
+               }
+               if !divisiblePow5(x, p) {
+                       t.Errorf("divisiblePow5(%d, %d) = false, want true", x, p)
+               }
+               if 2*x > x && !divisiblePow5(2*x, p) {
+                       t.Errorf("divisiblePow5(%d, %d) = false, want true", 2*x, p)
+               }
+       }
+}
+
+func TestDiv5Tab(t *testing.T) {
+       for p := 1; p <= 22; p++ {
+               m := div5Tab[p-1][0]
+               le := div5Tab[p-1][1]
+
+               // See comment in math.go on div5Tab.
+               // m needs to be multiplicative inverse of pow5(p).
+               if m*pow5(p) != 1 {
+                       t.Errorf("pow5Tab[%d-1][0] = %#x, but %#x * (5**%d) = %d, want 1", p, m, m, p, m*pow5(p))
+               }
+
+               // le needs to be ⌊(1<<64 - 1) / 5^p⌋.
+               want := (1<<64 - 1) / pow5(p)
+               if le != want {
+                       t.Errorf("pow5Tab[%d-1][1] = %#x, want %#x", p, le, want)
+               }
+       }
+}
+
+func TestTrimZeros(t *testing.T) {
+       for _, x := range []uint64{1, 2, 3, 4, 101, 123} {
+               want := x
+               for p := range 20 {
+                       haveX, haveP := trimZeros(x)
+                       if haveX != want || haveP != p {
+                               t.Errorf("trimZeros(%d) = %d, %d, want %d, %d", x, haveX, haveP, want, p)
+                       }
+                       if x >= (1<<64-1)/10 {
+                               break
+                       }
+                       x *= 10
+               }
+       }
+}
author	Russ Cox <rsc@golang.org>
	Sat, 1 Nov 2025 13:41:40 +0000 (09:41 -0400)
committer	Gopher Robot <gobot@golang.org>
	Tue, 4 Nov 2025 04:29:54 +0000 (20:29 -0800)
src/internal/strconv/atoi.go		patch \| blob \| history
src/internal/strconv/export_test.go		patch \| blob \| history
src/internal/strconv/ftoa.go		patch \| blob \| history
src/internal/strconv/ftoa_test.go		patch \| blob \| history
src/internal/strconv/ftoafixed.go	[new file with mode: 0644]	patch \| blob
src/internal/strconv/ftoaryu.go		patch \| blob \| history
src/internal/strconv/import_test.go		patch \| blob \| history
src/internal/strconv/math.go		patch \| blob \| history
src/internal/strconv/math_test.go		patch \| blob \| history