github.com/gnolang/gno@v0.0.0-20240520182011-228e9d0192ce/examples/gno.land/p/demo/json/ryu/ryu64.gno (about)

     1  package ryu
     2  
     3  import (
     4  	"math/bits"
     5  )
     6  
     7  type uint128 struct {
     8  	lo uint64
     9  	hi uint64
    10  }
    11  
    12  // dec64 is a floating decimal type representing m * 10^e.
    13  type dec64 struct {
    14  	m uint64
    15  	e int32
    16  }
    17  
    18  func (d dec64) append(b []byte, neg bool) []byte {
    19  	// Step 5: Print the decimal representation.
    20  	if neg {
    21  		b = append(b, '-')
    22  	}
    23  
    24  	out := d.m
    25  	outLen := decimalLen64(out)
    26  	bufLen := outLen
    27  	if bufLen > 1 {
    28  		bufLen++ // extra space for '.'
    29  	}
    30  
    31  	// Print the decimal digits.
    32  	n := len(b)
    33  	if cap(b)-len(b) >= bufLen {
    34  		// Avoid function call in the common case.
    35  		b = b[:len(b)+bufLen]
    36  	} else {
    37  		b = append(b, make([]byte, bufLen)...)
    38  	}
    39  
    40  	// Avoid expensive 64-bit divisions.
    41  	// We have at most 17 digits, and uint32 can store 9 digits.
    42  	// If the output doesn't fit into a uint32, cut off 8 digits
    43  	// so the rest will fit into a uint32.
    44  	var i int
    45  	if out>>32 > 0 {
    46  		var out32 uint32
    47  		out, out32 = out/1e8, uint32(out%1e8)
    48  		for ; i < 8; i++ {
    49  			b[n+outLen-i] = '0' + byte(out32%10)
    50  			out32 /= 10
    51  		}
    52  	}
    53  	out32 := uint32(out)
    54  	for ; i < outLen-1; i++ {
    55  		b[n+outLen-i] = '0' + byte(out32%10)
    56  		out32 /= 10
    57  	}
    58  	b[n] = '0' + byte(out32%10)
    59  
    60  	// Print the '.' if needed.
    61  	if outLen > 1 {
    62  		b[n+1] = '.'
    63  	}
    64  
    65  	// Print the exponent.
    66  	b = append(b, 'e')
    67  	exp := d.e + int32(outLen) - 1
    68  	if exp < 0 {
    69  		b = append(b, '-')
    70  		exp = -exp
    71  	} else {
    72  		// Unconditionally print a + here to match strconv's formatting.
    73  		b = append(b, '+')
    74  	}
    75  	// Always print at least two digits to match strconv's formatting.
    76  	d2 := exp % 10
    77  	exp /= 10
    78  	d1 := exp % 10
    79  	d0 := exp / 10
    80  	if d0 > 0 {
    81  		b = append(b, '0'+byte(d0))
    82  	}
    83  	b = append(b, '0'+byte(d1), '0'+byte(d2))
    84  
    85  	return b
    86  }
    87  
    88  func float64ToDecimalExactInt(mant, exp uint64) (d dec64, ok bool) {
    89  	e := exp - bias64
    90  	if e > mantBits64 {
    91  		return d, false
    92  	}
    93  	shift := mantBits64 - e
    94  	mant |= 1 << mantBits64 // implicit 1
    95  	d.m = mant >> shift
    96  	if d.m<<shift != mant {
    97  		return d, false
    98  	}
    99  
   100  	for d.m%10 == 0 {
   101  		d.m /= 10
   102  		d.e++
   103  	}
   104  	return d, true
   105  }
   106  
   107  func float64ToDecimal(mant, exp uint64) dec64 {
   108  	var e2 int32
   109  	var m2 uint64
   110  	if exp == 0 {
   111  		// We subtract 2 so that the bounds computation has
   112  		// 2 additional bits.
   113  		e2 = 1 - bias64 - mantBits64 - 2
   114  		m2 = mant
   115  	} else {
   116  		e2 = int32(exp) - bias64 - mantBits64 - 2
   117  		m2 = uint64(1)<<mantBits64 | mant
   118  	}
   119  	even := m2&1 == 0
   120  	acceptBounds := even
   121  
   122  	// Step 2: Determine the interval of valid decimal representations.
   123  	mv := 4 * m2
   124  	mmShift := boolToUint64(mant != 0 || exp <= 1)
   125  	// We would compute mp and mm like this:
   126  	// mp := 4 * m2 + 2;
   127  	// mm := mv - 1 - mmShift;
   128  
   129  	// Step 3: Convert to a decimal power base uing 128-bit arithmetic.
   130  	var (
   131  		vr, vp, vm        uint64
   132  		e10               int32
   133  		vmIsTrailingZeros bool
   134  		vrIsTrailingZeros bool
   135  	)
   136  	if e2 >= 0 {
   137  		// This expression is slightly faster than max(0, log10Pow2(e2) - 1).
   138  		q := log10Pow2(e2) - boolToUint32(e2 > 3)
   139  		e10 = int32(q)
   140  		k := pow5InvNumBits64 + pow5Bits(int32(q)) - 1
   141  		i := -e2 + int32(q) + k
   142  		mul := pow5InvSplit64[q]
   143  		vr = mulShift64(4*m2, mul, i)
   144  		vp = mulShift64(4*m2+2, mul, i)
   145  		vm = mulShift64(4*m2-1-mmShift, mul, i)
   146  		if q <= 21 {
   147  			// This should use q <= 22, but I think 21 is also safe.
   148  			// Smaller values may still be safe, but it's more
   149  			// difficult to reason about them. Only one of mp, mv,
   150  			// and mm can be a multiple of 5, if any.
   151  			if mv%5 == 0 {
   152  				vrIsTrailingZeros = multipleOfPowerOfFive64(mv, q)
   153  			} else if acceptBounds {
   154  				// Same as min(e2 + (^mm & 1), pow5Factor64(mm)) >= q
   155  				// <=> e2 + (^mm & 1) >= q && pow5Factor64(mm) >= q
   156  				// <=> true && pow5Factor64(mm) >= q, since e2 >= q.
   157  				vmIsTrailingZeros = multipleOfPowerOfFive64(mv-1-mmShift, q)
   158  			} else if multipleOfPowerOfFive64(mv+2, q) {
   159  				vp--
   160  			}
   161  		}
   162  	} else {
   163  		// This expression is slightly faster than max(0, log10Pow5(-e2) - 1).
   164  		q := log10Pow5(-e2) - boolToUint32(-e2 > 1)
   165  		e10 = int32(q) + e2
   166  		i := -e2 - int32(q)
   167  		k := pow5Bits(i) - pow5NumBits64
   168  		j := int32(q) - k
   169  		mul := pow5Split64[i]
   170  		vr = mulShift64(4*m2, mul, j)
   171  		vp = mulShift64(4*m2+2, mul, j)
   172  		vm = mulShift64(4*m2-1-mmShift, mul, j)
   173  		if q <= 1 {
   174  			// {vr,vp,vm} is trailing zeros if {mv,mp,mm} has at least q trailing 0 bits.
   175  			// mv = 4 * m2, so it always has at least two trailing 0 bits.
   176  			vrIsTrailingZeros = true
   177  			if acceptBounds {
   178  				// mm = mv - 1 - mmShift, so it has 1 trailing 0 bit iff mmShift == 1.
   179  				vmIsTrailingZeros = mmShift == 1
   180  			} else {
   181  				// mp = mv + 2, so it always has at least one trailing 0 bit.
   182  				vp--
   183  			}
   184  		} else if q < 63 { // TODO(ulfjack/cespare): Use a tighter bound here.
   185  			// We need to compute min(ntz(mv), pow5Factor64(mv) - e2) >= q - 1
   186  			// <=> ntz(mv) >= q - 1 && pow5Factor64(mv) - e2 >= q - 1
   187  			// <=> ntz(mv) >= q - 1 (e2 is negative and -e2 >= q)
   188  			// <=> (mv & ((1 << (q - 1)) - 1)) == 0
   189  			// We also need to make sure that the left shift does not overflow.
   190  			vrIsTrailingZeros = multipleOfPowerOfTwo64(mv, q-1)
   191  		}
   192  	}
   193  
   194  	// Step 4: Find the shortest decimal representation
   195  	// in the interval of valid representations.
   196  	var removed int32
   197  	var lastRemovedDigit uint8
   198  	var out uint64
   199  	// On average, we remove ~2 digits.
   200  	if vmIsTrailingZeros || vrIsTrailingZeros {
   201  		// General case, which happens rarely (~0.7%).
   202  		for {
   203  			vpDiv10 := vp / 10
   204  			vmDiv10 := vm / 10
   205  			if vpDiv10 <= vmDiv10 {
   206  				break
   207  			}
   208  			vmMod10 := vm % 10
   209  			vrDiv10 := vr / 10
   210  			vrMod10 := vr % 10
   211  			vmIsTrailingZeros = vmIsTrailingZeros && vmMod10 == 0
   212  			vrIsTrailingZeros = vrIsTrailingZeros && lastRemovedDigit == 0
   213  			lastRemovedDigit = uint8(vrMod10)
   214  			vr = vrDiv10
   215  			vp = vpDiv10
   216  			vm = vmDiv10
   217  			removed++
   218  		}
   219  		if vmIsTrailingZeros {
   220  			for {
   221  				vmDiv10 := vm / 10
   222  				vmMod10 := vm % 10
   223  				if vmMod10 != 0 {
   224  					break
   225  				}
   226  				vpDiv10 := vp / 10
   227  				vrDiv10 := vr / 10
   228  				vrMod10 := vr % 10
   229  				vrIsTrailingZeros = vrIsTrailingZeros && lastRemovedDigit == 0
   230  				lastRemovedDigit = uint8(vrMod10)
   231  				vr = vrDiv10
   232  				vp = vpDiv10
   233  				vm = vmDiv10
   234  				removed++
   235  			}
   236  		}
   237  		if vrIsTrailingZeros && lastRemovedDigit == 5 && vr%2 == 0 {
   238  			// Round even if the exact number is .....50..0.
   239  			lastRemovedDigit = 4
   240  		}
   241  		out = vr
   242  		// We need to take vr + 1 if vr is outside bounds
   243  		// or we need to round up.
   244  		if (vr == vm && (!acceptBounds || !vmIsTrailingZeros)) || lastRemovedDigit >= 5 {
   245  			out++
   246  		}
   247  	} else {
   248  		// Specialized for the common case (~99.3%).
   249  		// Percentages below are relative to this.
   250  		roundUp := false
   251  		for vp/100 > vm/100 {
   252  			// Optimization: remove two digits at a time (~86.2%).
   253  			roundUp = vr%100 >= 50
   254  			vr /= 100
   255  			vp /= 100
   256  			vm /= 100
   257  			removed += 2
   258  		}
   259  		// Loop iterations below (approximately), without optimization above:
   260  		// 0: 0.03%, 1: 13.8%, 2: 70.6%, 3: 14.0%, 4: 1.40%, 5: 0.14%, 6+: 0.02%
   261  		// Loop iterations below (approximately), with optimization above:
   262  		// 0: 70.6%, 1: 27.8%, 2: 1.40%, 3: 0.14%, 4+: 0.02%
   263  		for vp/10 > vm/10 {
   264  			roundUp = vr%10 >= 5
   265  			vr /= 10
   266  			vp /= 10
   267  			vm /= 10
   268  			removed++
   269  		}
   270  		// We need to take vr + 1 if vr is outside bounds
   271  		// or we need to round up.
   272  		out = vr + boolToUint64(vr == vm || roundUp)
   273  	}
   274  
   275  	return dec64{m: out, e: e10 + removed}
   276  }
   277  
   278  var powersOf10 = [...]uint64{
   279  	1e0,
   280  	1e1,
   281  	1e2,
   282  	1e3,
   283  	1e4,
   284  	1e5,
   285  	1e6,
   286  	1e7,
   287  	1e8,
   288  	1e9,
   289  	1e10,
   290  	1e11,
   291  	1e12,
   292  	1e13,
   293  	1e14,
   294  	1e15,
   295  	1e16,
   296  	1e17,
   297  	// We only need to find the length of at most 17 digit numbers.
   298  }
   299  
   300  func decimalLen64(u uint64) int {
   301  	// http://graphics.stanford.edu/~seander/bithacks.html#IntegerLog10
   302  	log2 := 64 - bits.LeadingZeros64(u) - 1
   303  	t := (log2 + 1) * 1233 >> 12
   304  	return t - boolToInt(u < powersOf10[t]) + 1
   305  }
   306  
   307  func mulShift64(m uint64, mul uint128, shift int32) uint64 {
   308  	hihi, hilo := bits.Mul64(m, mul.hi)
   309  	lohi, _ := bits.Mul64(m, mul.lo)
   310  	sum := uint128{hi: hihi, lo: lohi + hilo}
   311  	if sum.lo < lohi {
   312  		sum.hi++ // overflow
   313  	}
   314  	return shiftRight128(sum, shift-64)
   315  }
   316  
   317  func shiftRight128(v uint128, shift int32) uint64 {
   318  	// The shift value is always modulo 64.
   319  	// In the current implementation of the 64-bit version
   320  	// of Ryu, the shift value is always < 64.
   321  	// (It is in the range [2, 59].)
   322  	// Check this here in case a future change requires larger shift
   323  	// values. In this case this function needs to be adjusted.
   324  	assert(shift < 64, "shift < 64")
   325  	return (v.hi << uint64(64-shift)) | (v.lo >> uint(shift))
   326  }
   327  
   328  func pow5Factor64(v uint64) uint32 {
   329  	for n := uint32(0); ; n++ {
   330  		q, r := v/5, v%5
   331  		if r != 0 {
   332  			return n
   333  		}
   334  		v = q
   335  	}
   336  }
   337  
   338  func multipleOfPowerOfFive64(v uint64, p uint32) bool {
   339  	return pow5Factor64(v) >= p
   340  }
   341  
   342  func multipleOfPowerOfTwo64(v uint64, p uint32) bool {
   343  	return uint32(bits.TrailingZeros64(v)) >= p
   344  }