{{ message }}
Repository navigation
Expand file tree
/
Copy path_pylong.py
More file actions
329 lines (275 loc) · 10.5 KB
/
Copy path_pylong.py
File metadata and controls
329 lines (275 loc) · 10.5 KB
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
"""Python implementations of some algorithms for use by longobject.c.
The goal is to provide asymptotically faster algorithms that can be
used for operations on integers with many digits. In those cases, the
performance overhead of the Python implementation is not significant
since the asymptotic behavior is what dominates runtime. Functions
provided by this module should be considered private and not part of any
public API.
Note: for ease of maintainability, please prefer clear code and avoid
"micro-optimizations". This module will only be imported and used for
integers with a huge number of digits. Saving a few microseconds with
tricky or non-obvious code is not worth it. For people looking for
maximum performance, they should use something like gmpy2."""
import re
import decimal
try:
import _decimal
except ImportError:
_decimal = None
def int_to_decimal(n):
"""Asymptotically fast conversion of an 'int' to Decimal."""
# Function due to Tim Peters. See GH issue #90716 for details.
# https://github.com/python/cpython/issues/90716
#
# The implementation in longobject.c of base conversion algorithms
# between power-of-2 and non-power-of-2 bases are quadratic time.
# This function implements a divide-and-conquer algorithm that is
# faster for large numbers. Builds an equal decimal.Decimal in a
# "clever" recursive way. If we want a string representation, we
# apply str to _that_.
D = decimal.Decimal
D2 = D(2)
BITLIM = 128
mem = {}
def w2pow(w):
"""Return D(2)**w and store the result. Also possibly save some
intermediate results. In context, these are likely to be reused
across various levels of the conversion to Decimal."""
if (result := mem.get(w)) is None:
if w <= BITLIM:
result = D2**w
elif w - 1 in mem:
result = (t := mem[w - 1]) + t
else:
w2 = w >> 1
# If w happens to be odd, w-w2 is one larger then w2
# now. Recurse on the smaller first (w2), so that it's
# in the cache and the larger (w-w2) can be handled by
# the cheaper `w-1 in mem` branch instead.
result = w2pow(w2) * w2pow(w - w2)
mem[w] = result
return result
def inner(n, w):
if w <= BITLIM:
return D(n)
w2 = w >> 1
hi = n >> w2
lo = n - (hi << w2)
return inner(lo, w2) + inner(hi, w - w2) * w2pow(w2)
with decimal.localcontext() as ctx:
ctx.prec = decimal.MAX_PREC
ctx.Emax = decimal.MAX_EMAX
ctx.Emin = decimal.MIN_EMIN
ctx.traps[decimal.Inexact] = 1
if n < 0:
negate = True
n = -n
else:
negate = False
result = inner(n, n.bit_length())
if negate:
result = -result
return result
def int_to_decimal_string(n):
"""Asymptotically fast conversion of an 'int' to a decimal string."""
w = n.bit_length()
if w > 450_000 and _decimal is not None:
# It is only usable with the C decimal implementation.
# _pydecimal.py calls str() on very large integers, which in its
# turn calls int_to_decimal_string(), causing very deep recursion.
return str(int_to_decimal(n))
# Fallback algorithm for the case when the C decimal module isn't
# available. This algorithm is asymptotically worse than the algorithm
# using the decimal module, but better than the quadratic time
# implementation in longobject.c.
def inner(n, w):
if w <= 1000:
return str(n)
w2 = w >> 1
d = pow10_cache.get(w2)
if d is None:
d = pow10_cache[w2] = 5**w2 << w2 # 10**i = (5*2)**i = 5**i * 2**i
hi, lo = divmod(n, d)
return inner(hi, w - w2) + inner(lo, w2).zfill(w2)
# The estimation of the number of decimal digits.
# There is no harm in small error. If we guess too large, there may
# be leading 0's that need to be stripped. If we guess too small, we
# may need to call str() recursively for the remaining highest digits,
# which can still potentially be a large integer. This is manifested
# only if the number has way more than 10**15 digits, that exceeds
# the 52-bit physical address limit in both Intel64 and AMD64.
w = int(w * 0.3010299956639812 + 1) # log10(2)
pow10_cache = {}
if n < 0:
n = -n
sign = '-'
else:
sign = ''
s = inner(n, w)
if s[0] == '0' and n:
# If our guess of w is too large, there may be leading 0's that
# need to be stripped.
s = s.lstrip('0')
return sign + s
def _str_to_int_inner(s):
"""Asymptotically fast conversion of a 'str' to an 'int'."""
# Function due to Bjorn Martinsson. See GH issue #90716 for details.
# https://github.com/python/cpython/issues/90716
#
# The implementation in longobject.c of base conversion algorithms
# between power-of-2 and non-power-of-2 bases are quadratic time.
# This function implements a divide-and-conquer algorithm making use
# of Python's built in big int multiplication. Since Python uses the
# Karatsuba algorithm for multiplication, the time complexity
# of this function is O(len(s)**1.58).
DIGLIM = 2048
mem = {}
def w5pow(w):
"""Return 5**w and store the result.
Also possibly save some intermediate results. In context, these
are likely to be reused across various levels of the conversion
to 'int'.
"""
if (result := mem.get(w)) is None:
if w <= DIGLIM:
result = 5**w
elif w - 1 in mem:
result = mem[w - 1] * 5
else:
w2 = w >> 1
# If w happens to be odd, w-w2 is one larger then w2
# now. Recurse on the smaller first (w2), so that it's
# in the cache and the larger (w-w2) can be handled by
# the cheaper `w-1 in mem` branch instead.
result = w5pow(w2) * w5pow(w - w2)
mem[w] = result
return result
def inner(a, b):
if b - a <= DIGLIM:
return int(s[a:b])
mid = (a + b + 1) >> 1
return inner(mid, b) + ((inner(a, mid) * w5pow(b - mid)) << (b - mid))
return inner(0, len(s))
def int_from_string(s):
"""Asymptotically fast version of PyLong_FromString(), conversion
of a string of decimal digits into an 'int'."""
# PyLong_FromString() has already removed leading +/-, checked for invalid
# use of underscore characters, checked that string consists of only digits
# and underscores, and stripped leading whitespace. The input can still
# contain underscores and have trailing whitespace.
s = s.rstrip().replace('_', '')
return _str_to_int_inner(s)
def str_to_int(s):
"""Asymptotically fast version of decimal string to 'int' conversion."""
# FIXME: this doesn't support the full syntax that int() supports.
m = re.match(r'\s*([+-]?)([0-9_]+)\s*', s)
if not m:
raise ValueError('invalid literal for int() with base 10')
v = int_from_string(m.group(2))
if m.group(1) == '-':
v = -v
return v
# Fast integer division, based on code from Mark Dickinson, fast_div.py
# GH-47701. Additional refinements and optimizations by Bjorn Martinsson. The
# algorithm is due to Burnikel and Ziegler, in their paper "Fast Recursive
# Division".
_DIV_LIMIT = 4000
def _div2n1n(a, b, n):
"""Divide a 2n-bit nonnegative integer a by an n-bit positive integer
b, using a recursive divide-and-conquer algorithm.
Inputs:
n is a positive integer
b is a positive integer with exactly n bits
a is a nonnegative integer such that a < 2**n * b
Output:
(q, r) such that a = b*q+r and 0 <= r < b.
"""
if a.bit_length() - n <= _DIV_LIMIT:
return divmod(a, b)
pad = n & 1
if pad:
a <<= 1
b <<= 1
n += 1
half_n = n >> 1
mask = (1 << half_n) - 1
b1, b2 = b >> half_n, b & mask
q1, r = _div3n2n(a >> n, (a >> half_n) & mask, b, b1, b2, half_n)
q2, r = _div3n2n(r, a & mask, b, b1, b2, half_n)
if pad:
r >>= 1
return q1 << half_n | q2, r
def _div3n2n(a12, a3, b, b1, b2, n):
"""Helper function for _div2n1n; not intended to be called directly."""
if a12 >> n == b1:
q, r = (1 << n) - 1, a12 - (b1 << n) + b1
else:
q, r = _div2n1n(a12, b1, n)
r = (r << n | a3) - q * b2
while r < 0:
q -= 1
r += b
return q, r
def _int2digits(a, n):
"""Decompose non-negative int a into base 2**n
Input:
a is a non-negative integer
Output:
List of the digits of a in base 2**n in little-endian order,
meaning the most significant digit is last. The most
significant digit is guaranteed to be non-zero.
If a is 0 then the output is an empty list.
"""
a_digits = [0] * ((a.bit_length() + n - 1) // n)
def inner(x, L, R):
if L + 1 == R:
a_digits[L] = x
return
mid = (L + R) >> 1
shift = (mid - L) * n
upper = x >> shift
lower = x ^ (upper << shift)
inner(lower, L, mid)
inner(upper, mid, R)
if a:
inner(a, 0, len(a_digits))
return a_digits
def _digits2int(digits, n):
"""Combine base-2**n digits into an int. This function is the
inverse of `_int2digits`. For more details, see _int2digits.
"""
def inner(L, R):
if L + 1 == R:
return digits[L]
mid = (L + R) >> 1
shift = (mid - L) * n
return (inner(mid, R) << shift) + inner(L, mid)
return inner(0, len(digits)) if digits else 0
def _divmod_pos(a, b):
"""Divide a non-negative integer a by a positive integer b, giving
quotient and remainder."""
# Use grade-school algorithm in base 2**n, n = nbits(b)
n = b.bit_length()
a_digits = _int2digits(a, n)
r = 0
q_digits = []
for a_digit in reversed(a_digits):
q_digit, r = _div2n1n((r << n) + a_digit, b, n)
q_digits.append(q_digit)
q_digits.reverse()
q = _digits2int(q_digits, n)
return q, r
def int_divmod(a, b):
"""Asymptotically fast replacement for divmod, for 'int'.
Its time complexity is O(n**1.58), where n = #bits(a) + #bits(b).
"""
if b == 0:
raise ZeroDivisionError
elif b < 0:
q, r = int_divmod(-a, -b)
return q, -r
elif a < 0:
q, r = int_divmod(~a, b)
return ~q, b + ~r
else:
return _divmod_pos(a, b)
You can’t perform that action at this time.