move formatfloat out of system (#20195)

* move formatfloat out of system

* fixes doc

* Update changelog.md

* careless

* fixes

* deprecate system/formatfloat

* better handling

1 files changed, 432 insertions, 0 deletions

diff --git a/lib/std/private/schubfach.nim b/lib/std/private/schubfach.nim
new file mode 100644
index 000000000..5b965aaa7
--- /dev/null
+++ b/lib/std/private/schubfach.nim
@@ -0,0 +1,432 @@
+##  Copyright 2020 Alexander Bolz
+##
+##  Distributed under the Boost Software License, Version 1.0.
+##   (See accompanying file LICENSE_1_0.txt or copy at https://www.boost.org/LICENSE_1_0.txt)
+
+# --------------------------------------------------------------------------------------------------
+##  This file contains an implementation of the Schubfach algorithm as described in
+##
+##  [1] Raffaello Giulietti, "The Schubfach way to render doubles",
+##      https://drive.google.com/open?id=1luHhyQF9zKlM8yJ1nebU0OgVYhfC6CBN
+# --------------------------------------------------------------------------------------------------
+
+import std/private/digitsutils
+
+when defined(nimPreviewSlimSystem):
+  import std/assertions
+
+
+template sf_Assert(x: untyped): untyped =
+  assert(x)
+
+# ==================================================================================================
+#
+# ==================================================================================================
+
+type
+  ValueType = float32
+  BitsType = uint32
+  Single {.bycopy.} = object
+    bits: BitsType
+
+const
+  significandSize: int32 = 24
+  MaxExponent = 128
+  exponentBias: int32 = MaxExponent - 1 + (significandSize - 1)
+  maxIeeeExponent: BitsType = BitsType(2 * MaxExponent - 1)
+  hiddenBit: BitsType = BitsType(1) shl (significandSize - 1)
+  significandMask: BitsType = hiddenBit - 1
+  exponentMask: BitsType = maxIeeeExponent shl (significandSize - 1)
+  signMask: BitsType = not (not BitsType(0) shr 1)
+
+proc constructSingle(bits: BitsType): Single {.constructor.} =
+  result.bits = bits
+
+proc constructSingle(value: ValueType): Single {.constructor.} =
+  result.bits = cast[typeof(result.bits)](value)
+
+proc physicalSignificand(this: Single): BitsType {.noSideEffect.} =
+  return this.bits and significandMask
+
+proc physicalExponent(this: Single): BitsType {.noSideEffect.} =
+  return (this.bits and exponentMask) shr (significandSize - 1)
+
+proc isFinite(this: Single): bool {.noSideEffect.} =
+  return (this.bits and exponentMask) != exponentMask
+
+proc isInf(this: Single): bool {.noSideEffect.} =
+  return (this.bits and exponentMask) == exponentMask and
+      (this.bits and significandMask) == 0
+
+proc isNaN(this: Single): bool {.noSideEffect.} =
+  return (this.bits and exponentMask) == exponentMask and
+      (this.bits and significandMask) != 0
+
+proc isZero(this: Single): bool {.noSideEffect.} =
+  return (this.bits and not signMask) == 0
+
+proc signBit(this: Single): int {.noSideEffect.} =
+  return int((this.bits and signMask) != 0)
+
+# ==================================================================================================
+##  Returns floor(x / 2^n).
+##
+##  Technically, right-shift of negative integers is implementation defined...
+##  Should easily be optimized into SAR (or equivalent) instruction.
+
+proc floorDivPow2(x: int32; n: int32): int32 {.inline.} =
+  return x shr n
+
+##  Returns floor(log_10(2^e))
+##  static inline int32_t FloorLog10Pow2(int32_t e)
+##  {
+##      SF_ASSERT(e >= -1500);
+##      SF_ASSERT(e <=  1500);
+##      return FloorDivPow2(e * 1262611, 22);
+##  }
+##  Returns floor(log_10(3/4 2^e))
+##  static inline int32_t FloorLog10ThreeQuartersPow2(int32_t e)
+##  {
+##      SF_ASSERT(e >= -1500);
+##      SF_ASSERT(e <=  1500);
+##      return FloorDivPow2(e * 1262611 - 524031, 22);
+##  }
+##  Returns floor(log_2(10^e))
+
+proc floorLog2Pow10(e: int32): int32 {.inline.} =
+  sf_Assert(e >= -1233)
+  sf_Assert(e <= 1233)
+  return floorDivPow2(e * 1741647, 19)
+
+const
+  kMin: int32 = -31
+  kMax: int32 = 45
+  g: array[kMax - kMin + 1, uint64] = [0x81CEB32C4B43FCF5'u64, 0xA2425FF75E14FC32'u64,
+    0xCAD2F7F5359A3B3F'u64, 0xFD87B5F28300CA0E'u64, 0x9E74D1B791E07E49'u64,
+    0xC612062576589DDB'u64, 0xF79687AED3EEC552'u64, 0x9ABE14CD44753B53'u64,
+    0xC16D9A0095928A28'u64, 0xF1C90080BAF72CB2'u64, 0x971DA05074DA7BEF'u64,
+    0xBCE5086492111AEB'u64, 0xEC1E4A7DB69561A6'u64, 0x9392EE8E921D5D08'u64,
+    0xB877AA3236A4B44A'u64, 0xE69594BEC44DE15C'u64, 0x901D7CF73AB0ACDA'u64,
+    0xB424DC35095CD810'u64, 0xE12E13424BB40E14'u64, 0x8CBCCC096F5088CC'u64,
+    0xAFEBFF0BCB24AAFF'u64, 0xDBE6FECEBDEDD5BF'u64, 0x89705F4136B4A598'u64,
+    0xABCC77118461CEFD'u64, 0xD6BF94D5E57A42BD'u64, 0x8637BD05AF6C69B6'u64,
+    0xA7C5AC471B478424'u64, 0xD1B71758E219652C'u64, 0x83126E978D4FDF3C'u64,
+    0xA3D70A3D70A3D70B'u64, 0xCCCCCCCCCCCCCCCD'u64, 0x8000000000000000'u64,
+    0xA000000000000000'u64, 0xC800000000000000'u64, 0xFA00000000000000'u64,
+    0x9C40000000000000'u64, 0xC350000000000000'u64, 0xF424000000000000'u64,
+    0x9896800000000000'u64, 0xBEBC200000000000'u64, 0xEE6B280000000000'u64,
+    0x9502F90000000000'u64, 0xBA43B74000000000'u64, 0xE8D4A51000000000'u64,
+    0x9184E72A00000000'u64, 0xB5E620F480000000'u64, 0xE35FA931A0000000'u64,
+    0x8E1BC9BF04000000'u64, 0xB1A2BC2EC5000000'u64, 0xDE0B6B3A76400000'u64,
+    0x8AC7230489E80000'u64, 0xAD78EBC5AC620000'u64, 0xD8D726B7177A8000'u64,
+    0x878678326EAC9000'u64, 0xA968163F0A57B400'u64, 0xD3C21BCECCEDA100'u64,
+    0x84595161401484A0'u64, 0xA56FA5B99019A5C8'u64, 0xCECB8F27F4200F3A'u64,
+    0x813F3978F8940985'u64, 0xA18F07D736B90BE6'u64, 0xC9F2C9CD04674EDF'u64,
+    0xFC6F7C4045812297'u64, 0x9DC5ADA82B70B59E'u64, 0xC5371912364CE306'u64,
+    0xF684DF56C3E01BC7'u64, 0x9A130B963A6C115D'u64, 0xC097CE7BC90715B4'u64,
+    0xF0BDC21ABB48DB21'u64, 0x96769950B50D88F5'u64, 0xBC143FA4E250EB32'u64,
+    0xEB194F8E1AE525FE'u64, 0x92EFD1B8D0CF37BF'u64, 0xB7ABC627050305AE'u64,
+    0xE596B7B0C643C71A'u64, 0x8F7E32CE7BEA5C70'u64, 0xB35DBF821AE4F38C'u64]
+
+proc computePow10Single(k: int32): uint64 {.inline.} =
+  ##  There are unique beta and r such that 10^k = beta 2^r and
+  ##  2^63 <= beta < 2^64, namely r = floor(log_2 10^k) - 63 and
+  ##  beta = 2^-r 10^k.
+  ##  Let g = ceil(beta), so (g-1) 2^r < 10^k <= g 2^r, with the latter
+  ##  value being a pretty good overestimate for 10^k.
+  ##  NB: Since for all the required exponents k, we have g < 2^64,
+  ##      all constants can be stored in 128-bit integers.
+  sf_Assert(k >= kMin)
+  sf_Assert(k <= kMax)
+  return g[k - kMin]
+
+proc lo32(x: uint64): uint32 {.inline.} =
+  return cast[uint32](x)
+
+proc hi32(x: uint64): uint32 {.inline.} =
+  return cast[uint32](x shr 32)
+
+when defined(sizeof_Int128):
+  proc roundToOdd(g: uint64; cp: uint32): uint32 {.inline.} =
+    let p: uint128 = uint128(g) * cp
+    let y1: uint32 = lo32(cast[uint64](p shr 64))
+    let y0: uint32 = hi32(cast[uint64](p))
+    return y1 or uint32(y0 > 1)
+
+elif defined(vcc) and defined(cpu64):
+  proc umul128(x, y: uint64, z: ptr uint64): uint64 {.importc: "_umul128", header: "<intrin.h>".}
+  proc roundToOdd(g: uint64; cpHi: uint32): uint32 {.inline.} =
+    var p1: uint64 = 0
+    var p0: uint64 = umul128(g, cpHi, addr(p1))
+    let y1: uint32 = lo32(p1)
+    let y0: uint32 = hi32(p0)
+    return y1 or uint32(y0 > 1)
+
+else:
+  proc roundToOdd(g: uint64; cp: uint32): uint32 {.inline.} =
+    let b01: uint64 = uint64(lo32(g)) * cp
+    let b11: uint64 = uint64(hi32(g)) * cp
+    let hi: uint64 = b11 + hi32(b01)
+    let y1: uint32 = hi32(hi)
+    let y0: uint32 = lo32(hi)
+    return y1 or uint32(y0 > 1)
+
+##  Returns whether value is divisible by 2^e2
+
+proc multipleOfPow2(value: uint32; e2: int32): bool {.inline.} =
+  sf_Assert(e2 >= 0)
+  sf_Assert(e2 <= 31)
+  return (value and ((uint32(1) shl e2) - 1)) == 0
+
+type
+  FloatingDecimal32 {.bycopy.} = object
+    digits: uint32            ##  num_digits <= 9
+    exponent: int32
+
+proc toDecimal32(ieeeSignificand: uint32; ieeeExponent: uint32): FloatingDecimal32 {.
+    inline.} =
+  var c: uint32
+  var q: int32
+  if ieeeExponent != 0:
+    c = hiddenBit or ieeeSignificand
+    q = cast[int32](ieeeExponent) - exponentBias
+    if 0 <= -q and -q < significandSize and multipleOfPow2(c, -q):
+      return FloatingDecimal32(digits: c shr -q, exponent: 0'i32)
+  else:
+    c = ieeeSignificand
+    q = 1 - exponentBias
+  let isEven: bool = (c mod 2 == 0)
+  let lowerBoundaryIsCloser: bool = (ieeeSignificand == 0 and ieeeExponent > 1)
+  ##   const int32_t qb = q - 2;
+  let cbl: uint32 = 4 * c - 2 + uint32(lowerBoundaryIsCloser)
+  let cb: uint32 = 4 * c
+  let cbr: uint32 = 4 * c + 2
+  ##  (q * 1262611         ) >> 22 == floor(log_10(    2^q))
+  ##  (q * 1262611 - 524031) >> 22 == floor(log_10(3/4 2^q))
+  sf_Assert(q >= -1500)
+  sf_Assert(q <= 1500)
+  let k: int32 = floorDivPow2(q * 1262611 - (if lowerBoundaryIsCloser: 524031 else: 0), 22)
+  let h: int32 = q + floorLog2Pow10(-k) + 1
+  sf_Assert(h >= 1)
+  sf_Assert(h <= 4)
+  let pow10: uint64 = computePow10Single(-k)
+  let vbl: uint32 = roundToOdd(pow10, cbl shl h)
+  let vb: uint32 = roundToOdd(pow10, cb shl h)
+  let vbr: uint32 = roundToOdd(pow10, cbr shl h)
+  let lower: uint32 = vbl + uint32(not isEven)
+  let upper: uint32 = vbr - uint32(not isEven)
+  ##  See Figure 4 in [1].
+  ##  And the modifications in Figure 6.
+  let s: uint32 = vb div 4
+  ##  NB: 4 * s == vb & ~3 == vb & -4
+  if s >= 10:
+    let sp: uint32 = s div 10
+    ##  = vb / 40
+    let upInside: bool = lower <= 40 * sp
+    let wpInside: bool = 40 * sp + 40 <= upper
+    ##       if (up_inside || wp_inside) // NB: At most one of u' and w' is in R_v.
+    if upInside != wpInside:
+      return FloatingDecimal32(digits: sp + uint32(wpInside), exponent: k + 1)
+  let uInside: bool = lower <= 4 * s
+  let wInside: bool = 4 * s + 4 <= upper
+  if uInside != wInside:
+    return FloatingDecimal32(digits: s + uint32(wInside), exponent: k)
+  let mid: uint32 = 4 * s + 2
+  ##  = 2(s + t)
+  let roundUp: bool = vb > mid or (vb == mid and (s and 1) != 0)
+  return FloatingDecimal32(digits: s + uint32(roundUp), exponent: k)
+
+## ==================================================================================================
+##  ToChars
+## ==================================================================================================
+
+proc printDecimalDigitsBackwards(buf: var openArray[char]; pos: int; output: uint32): int32 {.inline.} =
+  var output = output
+  var pos = pos
+  var tz: int32 = 0
+  ##  number of trailing zeros removed.
+  var nd: int32 = 0
+  ##  number of decimal digits processed.
+  ##  At most 9 digits remaining
+  if output >= 10000:
+    let q: uint32 = output div 10000
+    let r: uint32 = output mod 10000
+    output = q
+    dec(pos, 4)
+    if r != 0:
+      let rH: uint32 = r div 100
+      let rL: uint32 = r mod 100
+      utoa2Digits(buf, pos, rH)
+      utoa2Digits(buf, pos + 2, rL)
+      tz = trailingZeros2Digits(if rL == 0: rH else: rL) + (if rL == 0: 2 else: 0)
+    else:
+      tz = 4
+    nd = 4
+  if output >= 100:
+    let q: uint32 = output div 100
+    let r: uint32 = output mod 100
+    output = q
+    dec(pos, 2)
+    utoa2Digits(buf, pos, r)
+    if tz == nd:
+      inc(tz, trailingZeros2Digits(r))
+    inc(nd, 2)
+    if output >= 100:
+      let q2: uint32 = output div 100
+      let r2: uint32 = output mod 100
+      output = q2
+      dec(pos, 2)
+      utoa2Digits(buf, pos, r2)
+      if tz == nd:
+        inc(tz, trailingZeros2Digits(r2))
+      inc(nd, 2)
+  sf_Assert(output >= 1)
+  sf_Assert(output <= 99)
+  if output >= 10:
+    let q: uint32 = output
+    dec(pos, 2)
+    utoa2Digits(buf, pos, q)
+    if tz == nd:
+      inc(tz, trailingZeros2Digits(q))
+  else:
+    let q: uint32 = output
+    sf_Assert(q >= 1)
+    sf_Assert(q <= 9)
+    dec(pos)
+    buf[pos] = chr(uint32('0') + q)
+  return tz
+
+proc decimalLength(v: uint32): int32 {.inline.} =
+  sf_Assert(v >= 1)
+  sf_Assert(v <= 999999999'u)
+  if v >= 100000000'u:
+    return 9
+  if v >= 10000000'u:
+    return 8
+  if v >= 1000000'u:
+    return 7
+  if v >= 100000'u:
+    return 6
+  if v >= 10000'u:
+    return 5
+  if v >= 1000'u:
+    return 4
+  if v >= 100'u:
+    return 3
+  if v >= 10'u:
+    return 2
+  return 1
+
+proc formatDigits(buffer: var openArray[char]; pos: int; digits: uint32; decimalExponent: int32;
+                  forceTrailingDotZero: bool = false): int {.inline.} =
+  const
+    minFixedDecimalPoint: int32 = -4
+    maxFixedDecimalPoint: int32 = 9
+  var pos = pos
+  assert(minFixedDecimalPoint <= -1, "internal error")
+  assert(maxFixedDecimalPoint >= 1, "internal error")
+  sf_Assert(digits >= 1)
+  sf_Assert(digits <= 999999999'u)
+  sf_Assert(decimalExponent >= -99)
+  sf_Assert(decimalExponent <= 99)
+  var numDigits: int32 = decimalLength(digits)
+  let decimalPoint: int32 = numDigits + decimalExponent
+  let useFixed: bool = minFixedDecimalPoint <= decimalPoint and
+      decimalPoint <= maxFixedDecimalPoint
+  ##  Prepare the buffer.
+  ##  Avoid calling memset/memcpy with variable arguments below...
+  for i in 0..<32: buffer[pos+i] = '0'
+  assert(minFixedDecimalPoint >= -30, "internal error")
+  assert(maxFixedDecimalPoint <= 32, "internal error")
+  var decimalDigitsPosition: int32
+  if useFixed:
+    if decimalPoint <= 0:
+      ##  0.[000]digits
+      decimalDigitsPosition = 2 - decimalPoint
+    else:
+      ##  dig.its
+      ##  digits[000]
+      decimalDigitsPosition = 0
+  else:
+    ##  dE+123 or d.igitsE+123
+    decimalDigitsPosition = 1
+  var digitsEnd = pos + decimalDigitsPosition + numDigits
+  let tz: int32 = printDecimalDigitsBackwards(buffer, digitsEnd, digits)
+  dec(digitsEnd, tz)
+  dec(numDigits, tz)
+  ##   decimal_exponent += tz; // => decimal_point unchanged.
+  if useFixed:
+    if decimalPoint <= 0:
+      ##  0.[000]digits
+      buffer[pos+1] = '.'
+      pos = digitsEnd
+    elif decimalPoint < numDigits:
+      ##  dig.its
+      for i in countdown(7, 0):
+        buffer[i + decimalPoint + 1] = buffer[i + decimalPoint]
+      buffer[pos+decimalPoint] = '.'
+      pos = digitsEnd + 1
+    else:
+      ##  digits[000]
+      inc(pos, decimalPoint)
+      if forceTrailingDotZero:
+        buffer[pos] = '.'
+        buffer[pos+1] = '0'
+        inc(pos, 2)
+  else:
+    buffer[pos] = buffer[pos+1]
+    if numDigits == 1:
+      ##  dE+123
+      inc(pos)
+    else:
+      ##  d.igitsE+123
+      buffer[pos+1] = '.'
+      pos = digitsEnd
+    let scientificExponent: int32 = decimalPoint - 1
+    ##       SF_ASSERT(scientific_exponent != 0);
+    buffer[pos] = 'e'
+    buffer[pos+1] = if scientificExponent < 0: '-' else: '+'
+    inc(pos, 2)
+    let k: uint32 = cast[uint32](if scientificExponent < 0: -scientificExponent else: scientificExponent)
+    if k < 10:
+      buffer[pos] = chr(uint32('0') + k)
+      inc pos
+    else:
+      utoa2Digits(buffer, pos, k)
+      inc(pos, 2)
+  return pos
+
+proc float32ToChars*(buffer: var openArray[char]; v: float32; forceTrailingDotZero = false): int {.
+    inline.} =
+  let significand: uint32 = physicalSignificand(constructSingle(v))
+  let exponent: uint32 = physicalExponent(constructSingle(v))
+  var pos = 0
+  if exponent != maxIeeeExponent:
+    ##  Finite
+    buffer[pos] = '-'
+    inc(pos, signBit(constructSingle(v)))
+    if exponent != 0 or significand != 0:
+      ##  != 0
+      let dec: auto = toDecimal32(significand, exponent)
+      return formatDigits(buffer, pos, dec.digits, dec.exponent, forceTrailingDotZero)
+    else:
+      buffer[pos] = '0'
+      buffer[pos+1] = '.'
+      buffer[pos+2] = '0'
+      buffer[pos+3] = ' '
+      inc(pos, if forceTrailingDotZero: 3 else: 1)
+      return pos
+  if significand == 0:
+    buffer[pos] = '-'
+    inc(pos, signBit(constructSingle(v)))
+    buffer[pos] = 'i'
+    buffer[pos+1] = 'n'
+    buffer[pos+2] = 'f'
+    buffer[pos+3] = ' '
+    return pos + 3
+  else:
+    buffer[pos] = 'n'
+    buffer[pos+1] = 'a'
+    buffer[pos+2] = 'n'
+    buffer[pos+3] = ' '
+    return pos + 3