]>
Commit | Line | Data |
---|---|---|
ba9703b0 | 1 | #!/usr/bin/env python3 |
e9174d1e SL |
2 | |
3 | """ | |
136023e0 | 4 | Generate powers of five using Daniel Lemire's ``Eisel-Lemire algorithm`` for use in |
e9174d1e SL |
5 | decimal to floating point conversions. |
6 | ||
7 | Specifically, computes and outputs (as Rust code) a table of 10^e for some | |
136023e0 XL |
8 | range of exponents e. The output is one array of 128 bit significands. |
9 | The base two exponents can be inferred using a logarithmic slope | |
10 | of the decimal exponent. The approximations are normalized and rounded perfectly, | |
11 | i.e., within 0.5 ULP of the true value. | |
e9174d1e | 12 | |
136023e0 XL |
13 | Adapted from Daniel Lemire's fast_float ``table_generation.py``, |
14 | available here: <https://github.com/fastfloat/fast_float/blob/main/script/table_generation.py>. | |
e9174d1e SL |
15 | """ |
16 | from __future__ import print_function | |
136023e0 | 17 | from math import ceil, floor, log, log2 |
e9174d1e | 18 | from fractions import Fraction |
136023e0 | 19 | from collections import deque |
74b04a01 | 20 | |
9cc50fc6 | 21 | HEADER = """ |
136023e0 XL |
22 | //! Pre-computed tables powers-of-5 for extended-precision representations. |
23 | //! | |
24 | //! These tables enable fast scaling of the significant digits | |
25 | //! of a float to the decimal exponent, with minimal rounding | |
26 | //! errors, in a 128 or 192-bit representation. | |
27 | //! | |
9cc50fc6 | 28 | //! DO NOT MODIFY: Generated by `src/etc/dec2flt_table.py` |
e9174d1e SL |
29 | """ |
30 | ||
136023e0 XL |
31 | STATIC_WARNING = """ |
32 | // Use static to avoid long compile times: Rust compiler errors | |
33 | // can have the entire table compiled multiple times, and then | |
34 | // emit code multiple times, even if it's stripped out in | |
35 | // the final binary. | |
36 | """ | |
9cc50fc6 | 37 | |
e9174d1e | 38 | def main(): |
136023e0 XL |
39 | min_exp = minimum_exponent(10) |
40 | max_exp = maximum_exponent(10) | |
41 | bias = -minimum_exponent(5) | |
42 | ||
9cc50fc6 SL |
43 | print(HEADER.strip()) |
44 | print() | |
136023e0 XL |
45 | print('pub const SMALLEST_POWER_OF_FIVE: i32 = {};'.format(min_exp)) |
46 | print('pub const LARGEST_POWER_OF_FIVE: i32 = {};'.format(max_exp)) | |
47 | print('pub const N_POWERS_OF_FIVE: usize = ', end='') | |
48 | print('(LARGEST_POWER_OF_FIVE - SMALLEST_POWER_OF_FIVE + 1) as usize;') | |
9cc50fc6 | 49 | print() |
136023e0 XL |
50 | print_proper_powers(min_exp, max_exp, bias) |
51 | ||
52 | ||
53 | def minimum_exponent(base): | |
54 | return ceil(log(5e-324, base) - log(0xFFFFFFFFFFFFFFFF, base)) | |
55 | ||
9cc50fc6 | 56 | |
136023e0 XL |
57 | def maximum_exponent(base): |
58 | return floor(log(1.7976931348623157e+308, base)) | |
9cc50fc6 | 59 | |
136023e0 XL |
60 | |
61 | def print_proper_powers(min_exp, max_exp, bias): | |
62 | powers = deque() | |
63 | ||
64 | # Add negative exponents. | |
65 | # 2^(2b)/(5^−q) with b=64 + int(math.ceil(log2(5^−q))) | |
e9174d1e | 66 | powers = [] |
136023e0 XL |
67 | for q in range(min_exp, 0): |
68 | power5 = 5 ** -q | |
69 | z = 0 | |
70 | while (1 << z) < power5: | |
71 | z += 1 | |
72 | if q >= -27: | |
73 | b = z + 127 | |
74 | c = 2 ** b // power5 + 1 | |
75 | powers.append((c, q)) | |
76 | else: | |
77 | b = 2 * z + 2 * 64 | |
78 | c = 2 ** b // power5 + 1 | |
79 | # truncate | |
80 | while c >= (1<<128): | |
81 | c //= 2 | |
82 | powers.append((c, q)) | |
83 | ||
84 | # Add positive exponents | |
85 | for q in range(0, max_exp + 1): | |
86 | power5 = 5 ** q | |
87 | # move the most significant bit in position | |
88 | while power5 < (1<<127): | |
89 | power5 *= 2 | |
90 | # *truncate* | |
91 | while power5 >= (1<<128): | |
92 | power5 //= 2 | |
93 | powers.append((power5, q)) | |
94 | ||
95 | # Print the powers. | |
96 | print(STATIC_WARNING.strip()) | |
97 | print('#[rustfmt::skip]') | |
98 | typ = '[(u64, u64); N_POWERS_OF_FIVE]' | |
99 | print('pub static POWER_OF_FIVE_128: {} = ['.format(typ)) | |
100 | lo_mask = (1 << 64) - 1 | |
101 | for c, exp in powers: | |
102 | hi = '0x{:x}'.format(c // (1 << 64)) | |
103 | lo = '0x{:x}'.format(c % (1 << 64)) | |
104 | value = ' ({}, {}), '.format(hi, lo) | |
105 | comment = '// {}^{}'.format(5, exp) | |
106 | print(value.ljust(46, ' ') + comment) | |
107 | print('];') | |
9cc50fc6 SL |
108 | |
109 | ||
e9174d1e SL |
110 | if __name__ == '__main__': |
111 | main() |