4 Generate powers of ten using William Clinger's ``AlgorithmM`` for use in
5 decimal to floating point conversions.
7 Specifically, computes and outputs (as Rust code) a table of 10^e for some
8 range of exponents e. The output is one array of 64 bit significands and
9 another array of corresponding base two exponents. The approximations are
10 normalized and rounded perfectly, i.e., within 0.5 ULP of the true value.
12 The representation ([u64], [i16]) instead of the more natural [(u64, i16)]
13 is used because (u64, i16) has a ton of padding which would make the table
14 even larger, and it's already uncomfortably large (6 KiB).
16 from __future__ import print_function
17 from math import ceil, log
18 from fractions import Fraction
19 from collections import namedtuple
22 N = 64 # Size of the significand field in bits
23 MIN_SIG = 2 ** (N - 1)
24 MAX_SIG = (2 ** N) - 1
26 # Hand-rolled fp representation without arithmetic or any other operations.
27 # The significand is normalized and always N bit, but the exponent is
28 # unrestricted in range.
29 Fp = namedtuple('Fp', 'sig exp')
32 def algorithm_m(f, e):
52 return ratio_to_float(u, v, k)
55 def ratio_to_float(u, v, k):
71 return Fp(MIN_SIG, z.exp + 1)
73 return Fp(z.sig + 1, z.exp)
77 decimal = f * Fraction(10) ** e
78 binary = z.sig * Fraction(2) ** z.exp
79 abs_err = abs(decimal - binary)
80 # The unit in the last place has value z.exp
81 ulp_err = abs_err / Fraction(2) ** z.exp
86 //! Tables of approximations of powers of ten.
87 //! DO NOT MODIFY: Generated by `src/etc/dec2flt_table.py`
96 print_short_powers(32, 24)
98 print_short_powers(64, 53)
101 def print_proper_powers():
104 e_range = range(MIN_E, MAX_E+1)
107 z = algorithm_m(1, e)
111 print("pub const MIN_E: i16 = {};".format(MIN_E))
112 print("pub const MAX_E: i16 = {};".format(MAX_E))
114 print("#[rustfmt::skip]")
115 typ = "([u64; {0}], [i16; {0}])".format(len(powers))
116 print("pub static POWERS: ", typ, " = (", sep='')
119 print(" 0x{:x},".format(z.sig))
123 print(" {},".format(z.exp))
128 def print_short_powers(num_bits, significand_size):
129 max_sig = 2**significand_size - 1
130 # The fast path bails out for exponents >= ceil(log5(max_sig))
131 max_e = int(ceil(log(max_sig, 5)))
132 e_range = range(max_e)
133 typ = "[f{}; {}]".format(num_bits, len(e_range))
134 print("#[rustfmt::skip]")
135 print("pub const F", num_bits, "_SHORT_POWERS: ", typ, " = [", sep='')
137 print(" 1e{},".format(e))
141 if __name__ == '__main__':