heapq.py 22 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601
  1. """Heap queue algorithm (a.k.a. priority queue).
  2. Heaps are arrays for which a[k] <= a[2*k+1] and a[k] <= a[2*k+2] for
  3. all k, counting elements from 0. For the sake of comparison,
  4. non-existing elements are considered to be infinite. The interesting
  5. property of a heap is that a[0] is always its smallest element.
  6. Usage:
  7. heap = [] # creates an empty heap
  8. heappush(heap, item) # pushes a new item on the heap
  9. item = heappop(heap) # pops the smallest item from the heap
  10. item = heap[0] # smallest item on the heap without popping it
  11. heapify(x) # transforms list into a heap, in-place, in linear time
  12. item = heapreplace(heap, item) # pops and returns smallest item, and adds
  13. # new item; the heap size is unchanged
  14. Our API differs from textbook heap algorithms as follows:
  15. - We use 0-based indexing. This makes the relationship between the
  16. index for a node and the indexes for its children slightly less
  17. obvious, but is more suitable since Python uses 0-based indexing.
  18. - Our heappop() method returns the smallest item, not the largest.
  19. These two make it possible to view the heap as a regular Python list
  20. without surprises: heap[0] is the smallest item, and heap.sort()
  21. maintains the heap invariant!
  22. """
  23. # Original code by Kevin O'Connor, augmented by Tim Peters and Raymond Hettinger
  24. __about__ = """Heap queues
  25. [explanation by François Pinard]
  26. Heaps are arrays for which a[k] <= a[2*k+1] and a[k] <= a[2*k+2] for
  27. all k, counting elements from 0. For the sake of comparison,
  28. non-existing elements are considered to be infinite. The interesting
  29. property of a heap is that a[0] is always its smallest element.
  30. The strange invariant above is meant to be an efficient memory
  31. representation for a tournament. The numbers below are `k', not a[k]:
  32. 0
  33. 1 2
  34. 3 4 5 6
  35. 7 8 9 10 11 12 13 14
  36. 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30
  37. In the tree above, each cell `k' is topping `2*k+1' and `2*k+2'. In
  38. a usual binary tournament we see in sports, each cell is the winner
  39. over the two cells it tops, and we can trace the winner down the tree
  40. to see all opponents s/he had. However, in many computer applications
  41. of such tournaments, we do not need to trace the history of a winner.
  42. To be more memory efficient, when a winner is promoted, we try to
  43. replace it by something else at a lower level, and the rule becomes
  44. that a cell and the two cells it tops contain three different items,
  45. but the top cell "wins" over the two topped cells.
  46. If this heap invariant is protected at all time, index 0 is clearly
  47. the overall winner. The simplest algorithmic way to remove it and
  48. find the "next" winner is to move some loser (let's say cell 30 in the
  49. diagram above) into the 0 position, and then percolate this new 0 down
  50. the tree, exchanging values, until the invariant is re-established.
  51. This is clearly logarithmic on the total number of items in the tree.
  52. By iterating over all items, you get an O(n ln n) sort.
  53. A nice feature of this sort is that you can efficiently insert new
  54. items while the sort is going on, provided that the inserted items are
  55. not "better" than the last 0'th element you extracted. This is
  56. especially useful in simulation contexts, where the tree holds all
  57. incoming events, and the "win" condition means the smallest scheduled
  58. time. When an event schedule other events for execution, they are
  59. scheduled into the future, so they can easily go into the heap. So, a
  60. heap is a good structure for implementing schedulers (this is what I
  61. used for my MIDI sequencer :-).
  62. Various structures for implementing schedulers have been extensively
  63. studied, and heaps are good for this, as they are reasonably speedy,
  64. the speed is almost constant, and the worst case is not much different
  65. than the average case. However, there are other representations which
  66. are more efficient overall, yet the worst cases might be terrible.
  67. Heaps are also very useful in big disk sorts. You most probably all
  68. know that a big sort implies producing "runs" (which are pre-sorted
  69. sequences, which size is usually related to the amount of CPU memory),
  70. followed by a merging passes for these runs, which merging is often
  71. very cleverly organised[1]. It is very important that the initial
  72. sort produces the longest runs possible. Tournaments are a good way
  73. to that. If, using all the memory available to hold a tournament, you
  74. replace and percolate items that happen to fit the current run, you'll
  75. produce runs which are twice the size of the memory for random input,
  76. and much better for input fuzzily ordered.
  77. Moreover, if you output the 0'th item on disk and get an input which
  78. may not fit in the current tournament (because the value "wins" over
  79. the last output value), it cannot fit in the heap, so the size of the
  80. heap decreases. The freed memory could be cleverly reused immediately
  81. for progressively building a second heap, which grows at exactly the
  82. same rate the first heap is melting. When the first heap completely
  83. vanishes, you switch heaps and start a new run. Clever and quite
  84. effective!
  85. In a word, heaps are useful memory structures to know. I use them in
  86. a few applications, and I think it is good to keep a `heap' module
  87. around. :-)
  88. --------------------
  89. [1] The disk balancing algorithms which are current, nowadays, are
  90. more annoying than clever, and this is a consequence of the seeking
  91. capabilities of the disks. On devices which cannot seek, like big
  92. tape drives, the story was quite different, and one had to be very
  93. clever to ensure (far in advance) that each tape movement will be the
  94. most effective possible (that is, will best participate at
  95. "progressing" the merge). Some tapes were even able to read
  96. backwards, and this was also used to avoid the rewinding time.
  97. Believe me, real good tape sorts were quite spectacular to watch!
  98. From all times, sorting has always been a Great Art! :-)
  99. """
  100. __all__ = ['heappush', 'heappop', 'heapify', 'heapreplace', 'merge',
  101. 'nlargest', 'nsmallest', 'heappushpop']
  102. def heappush(heap, item):
  103. """Push item onto heap, maintaining the heap invariant."""
  104. heap.append(item)
  105. _siftdown(heap, 0, len(heap)-1)
  106. def heappop(heap):
  107. """Pop the smallest item off the heap, maintaining the heap invariant."""
  108. lastelt = heap.pop() # raises appropriate IndexError if heap is empty
  109. if heap:
  110. returnitem = heap[0]
  111. heap[0] = lastelt
  112. _siftup(heap, 0)
  113. return returnitem
  114. return lastelt
  115. def heapreplace(heap, item):
  116. """Pop and return the current smallest value, and add the new item.
  117. This is more efficient than heappop() followed by heappush(), and can be
  118. more appropriate when using a fixed-size heap. Note that the value
  119. returned may be larger than item! That constrains reasonable uses of
  120. this routine unless written as part of a conditional replacement:
  121. if item > heap[0]:
  122. item = heapreplace(heap, item)
  123. """
  124. returnitem = heap[0] # raises appropriate IndexError if heap is empty
  125. heap[0] = item
  126. _siftup(heap, 0)
  127. return returnitem
  128. def heappushpop(heap, item):
  129. """Fast version of a heappush followed by a heappop."""
  130. if heap and heap[0] < item:
  131. item, heap[0] = heap[0], item
  132. _siftup(heap, 0)
  133. return item
  134. def heapify(x):
  135. """Transform list into a heap, in-place, in O(len(x)) time."""
  136. n = len(x)
  137. # Transform bottom-up. The largest index there's any point to looking at
  138. # is the largest with a child index in-range, so must have 2*i + 1 < n,
  139. # or i < (n-1)/2. If n is even = 2*j, this is (2*j-1)/2 = j-1/2 so
  140. # j-1 is the largest, which is n//2 - 1. If n is odd = 2*j+1, this is
  141. # (2*j+1-1)/2 = j so j-1 is the largest, and that's again n//2-1.
  142. for i in reversed(range(n//2)):
  143. _siftup(x, i)
  144. def _heappop_max(heap):
  145. """Maxheap version of a heappop."""
  146. lastelt = heap.pop() # raises appropriate IndexError if heap is empty
  147. if heap:
  148. returnitem = heap[0]
  149. heap[0] = lastelt
  150. _siftup_max(heap, 0)
  151. return returnitem
  152. return lastelt
  153. def _heapreplace_max(heap, item):
  154. """Maxheap version of a heappop followed by a heappush."""
  155. returnitem = heap[0] # raises appropriate IndexError if heap is empty
  156. heap[0] = item
  157. _siftup_max(heap, 0)
  158. return returnitem
  159. def _heapify_max(x):
  160. """Transform list into a maxheap, in-place, in O(len(x)) time."""
  161. n = len(x)
  162. for i in reversed(range(n//2)):
  163. _siftup_max(x, i)
  164. # 'heap' is a heap at all indices >= startpos, except possibly for pos. pos
  165. # is the index of a leaf with a possibly out-of-order value. Restore the
  166. # heap invariant.
  167. def _siftdown(heap, startpos, pos):
  168. newitem = heap[pos]
  169. # Follow the path to the root, moving parents down until finding a place
  170. # newitem fits.
  171. while pos > startpos:
  172. parentpos = (pos - 1) >> 1
  173. parent = heap[parentpos]
  174. if newitem < parent:
  175. heap[pos] = parent
  176. pos = parentpos
  177. continue
  178. break
  179. heap[pos] = newitem
  180. # The child indices of heap index pos are already heaps, and we want to make
  181. # a heap at index pos too. We do this by bubbling the smaller child of
  182. # pos up (and so on with that child's children, etc) until hitting a leaf,
  183. # then using _siftdown to move the oddball originally at index pos into place.
  184. #
  185. # We *could* break out of the loop as soon as we find a pos where newitem <=
  186. # both its children, but turns out that's not a good idea, and despite that
  187. # many books write the algorithm that way. During a heap pop, the last array
  188. # element is sifted in, and that tends to be large, so that comparing it
  189. # against values starting from the root usually doesn't pay (= usually doesn't
  190. # get us out of the loop early). See Knuth, Volume 3, where this is
  191. # explained and quantified in an exercise.
  192. #
  193. # Cutting the # of comparisons is important, since these routines have no
  194. # way to extract "the priority" from an array element, so that intelligence
  195. # is likely to be hiding in custom comparison methods, or in array elements
  196. # storing (priority, record) tuples. Comparisons are thus potentially
  197. # expensive.
  198. #
  199. # On random arrays of length 1000, making this change cut the number of
  200. # comparisons made by heapify() a little, and those made by exhaustive
  201. # heappop() a lot, in accord with theory. Here are typical results from 3
  202. # runs (3 just to demonstrate how small the variance is):
  203. #
  204. # Compares needed by heapify Compares needed by 1000 heappops
  205. # -------------------------- --------------------------------
  206. # 1837 cut to 1663 14996 cut to 8680
  207. # 1855 cut to 1659 14966 cut to 8678
  208. # 1847 cut to 1660 15024 cut to 8703
  209. #
  210. # Building the heap by using heappush() 1000 times instead required
  211. # 2198, 2148, and 2219 compares: heapify() is more efficient, when
  212. # you can use it.
  213. #
  214. # The total compares needed by list.sort() on the same lists were 8627,
  215. # 8627, and 8632 (this should be compared to the sum of heapify() and
  216. # heappop() compares): list.sort() is (unsurprisingly!) more efficient
  217. # for sorting.
  218. def _siftup(heap, pos):
  219. endpos = len(heap)
  220. startpos = pos
  221. newitem = heap[pos]
  222. # Bubble up the smaller child until hitting a leaf.
  223. childpos = 2*pos + 1 # leftmost child position
  224. while childpos < endpos:
  225. # Set childpos to index of smaller child.
  226. rightpos = childpos + 1
  227. if rightpos < endpos and not heap[childpos] < heap[rightpos]:
  228. childpos = rightpos
  229. # Move the smaller child up.
  230. heap[pos] = heap[childpos]
  231. pos = childpos
  232. childpos = 2*pos + 1
  233. # The leaf at pos is empty now. Put newitem there, and bubble it up
  234. # to its final resting place (by sifting its parents down).
  235. heap[pos] = newitem
  236. _siftdown(heap, startpos, pos)
  237. def _siftdown_max(heap, startpos, pos):
  238. 'Maxheap variant of _siftdown'
  239. newitem = heap[pos]
  240. # Follow the path to the root, moving parents down until finding a place
  241. # newitem fits.
  242. while pos > startpos:
  243. parentpos = (pos - 1) >> 1
  244. parent = heap[parentpos]
  245. if parent < newitem:
  246. heap[pos] = parent
  247. pos = parentpos
  248. continue
  249. break
  250. heap[pos] = newitem
  251. def _siftup_max(heap, pos):
  252. 'Maxheap variant of _siftup'
  253. endpos = len(heap)
  254. startpos = pos
  255. newitem = heap[pos]
  256. # Bubble up the larger child until hitting a leaf.
  257. childpos = 2*pos + 1 # leftmost child position
  258. while childpos < endpos:
  259. # Set childpos to index of larger child.
  260. rightpos = childpos + 1
  261. if rightpos < endpos and not heap[rightpos] < heap[childpos]:
  262. childpos = rightpos
  263. # Move the larger child up.
  264. heap[pos] = heap[childpos]
  265. pos = childpos
  266. childpos = 2*pos + 1
  267. # The leaf at pos is empty now. Put newitem there, and bubble it up
  268. # to its final resting place (by sifting its parents down).
  269. heap[pos] = newitem
  270. _siftdown_max(heap, startpos, pos)
  271. def merge(*iterables, key=None, reverse=False):
  272. '''Merge multiple sorted inputs into a single sorted output.
  273. Similar to sorted(itertools.chain(*iterables)) but returns a generator,
  274. does not pull the data into memory all at once, and assumes that each of
  275. the input streams is already sorted (smallest to largest).
  276. >>> list(merge([1,3,5,7], [0,2,4,8], [5,10,15,20], [], [25]))
  277. [0, 1, 2, 3, 4, 5, 5, 7, 8, 10, 15, 20, 25]
  278. If *key* is not None, applies a key function to each element to determine
  279. its sort order.
  280. >>> list(merge(['dog', 'horse'], ['cat', 'fish', 'kangaroo'], key=len))
  281. ['dog', 'cat', 'fish', 'horse', 'kangaroo']
  282. '''
  283. h = []
  284. h_append = h.append
  285. if reverse:
  286. _heapify = _heapify_max
  287. _heappop = _heappop_max
  288. _heapreplace = _heapreplace_max
  289. direction = -1
  290. else:
  291. _heapify = heapify
  292. _heappop = heappop
  293. _heapreplace = heapreplace
  294. direction = 1
  295. if key is None:
  296. for order, it in enumerate(map(iter, iterables)):
  297. try:
  298. next = it.__next__
  299. h_append([next(), order * direction, next])
  300. except StopIteration:
  301. pass
  302. _heapify(h)
  303. while len(h) > 1:
  304. try:
  305. while True:
  306. value, order, next = s = h[0]
  307. yield value
  308. s[0] = next() # raises StopIteration when exhausted
  309. _heapreplace(h, s) # restore heap condition
  310. except StopIteration:
  311. _heappop(h) # remove empty iterator
  312. if h:
  313. # fast case when only a single iterator remains
  314. value, order, next = h[0]
  315. yield value
  316. yield from next.__self__
  317. return
  318. for order, it in enumerate(map(iter, iterables)):
  319. try:
  320. next = it.__next__
  321. value = next()
  322. h_append([key(value), order * direction, value, next])
  323. except StopIteration:
  324. pass
  325. _heapify(h)
  326. while len(h) > 1:
  327. try:
  328. while True:
  329. key_value, order, value, next = s = h[0]
  330. yield value
  331. value = next()
  332. s[0] = key(value)
  333. s[2] = value
  334. _heapreplace(h, s)
  335. except StopIteration:
  336. _heappop(h)
  337. if h:
  338. key_value, order, value, next = h[0]
  339. yield value
  340. yield from next.__self__
  341. # Algorithm notes for nlargest() and nsmallest()
  342. # ==============================================
  343. #
  344. # Make a single pass over the data while keeping the k most extreme values
  345. # in a heap. Memory consumption is limited to keeping k values in a list.
  346. #
  347. # Measured performance for random inputs:
  348. #
  349. # number of comparisons
  350. # n inputs k-extreme values (average of 5 trials) % more than min()
  351. # ------------- ---------------- --------------------- -----------------
  352. # 1,000 100 3,317 231.7%
  353. # 10,000 100 14,046 40.5%
  354. # 100,000 100 105,749 5.7%
  355. # 1,000,000 100 1,007,751 0.8%
  356. # 10,000,000 100 10,009,401 0.1%
  357. #
  358. # Theoretical number of comparisons for k smallest of n random inputs:
  359. #
  360. # Step Comparisons Action
  361. # ---- -------------------------- ---------------------------
  362. # 1 1.66 * k heapify the first k-inputs
  363. # 2 n - k compare remaining elements to top of heap
  364. # 3 k * (1 + lg2(k)) * ln(n/k) replace the topmost value on the heap
  365. # 4 k * lg2(k) - (k/2) final sort of the k most extreme values
  366. #
  367. # Combining and simplifying for a rough estimate gives:
  368. #
  369. # comparisons = n + k * (log(k, 2) * log(n/k) + log(k, 2) + log(n/k))
  370. #
  371. # Computing the number of comparisons for step 3:
  372. # -----------------------------------------------
  373. # * For the i-th new value from the iterable, the probability of being in the
  374. # k most extreme values is k/i. For example, the probability of the 101st
  375. # value seen being in the 100 most extreme values is 100/101.
  376. # * If the value is a new extreme value, the cost of inserting it into the
  377. # heap is 1 + log(k, 2).
  378. # * The probability times the cost gives:
  379. # (k/i) * (1 + log(k, 2))
  380. # * Summing across the remaining n-k elements gives:
  381. # sum((k/i) * (1 + log(k, 2)) for i in range(k+1, n+1))
  382. # * This reduces to:
  383. # (H(n) - H(k)) * k * (1 + log(k, 2))
  384. # * Where H(n) is the n-th harmonic number estimated by:
  385. # gamma = 0.5772156649
  386. # H(n) = log(n, e) + gamma + 1 / (2 * n)
  387. # http://en.wikipedia.org/wiki/Harmonic_series_(mathematics)#Rate_of_divergence
  388. # * Substituting the H(n) formula:
  389. # comparisons = k * (1 + log(k, 2)) * (log(n/k, e) + (1/n - 1/k) / 2)
  390. #
  391. # Worst-case for step 3:
  392. # ----------------------
  393. # In the worst case, the input data is reversed sorted so that every new element
  394. # must be inserted in the heap:
  395. #
  396. # comparisons = 1.66 * k + log(k, 2) * (n - k)
  397. #
  398. # Alternative Algorithms
  399. # ----------------------
  400. # Other algorithms were not used because they:
  401. # 1) Took much more auxiliary memory,
  402. # 2) Made multiple passes over the data.
  403. # 3) Made more comparisons in common cases (small k, large n, semi-random input).
  404. # See the more detailed comparison of approach at:
  405. # http://code.activestate.com/recipes/577573-compare-algorithms-for-heapqsmallest
  406. def nsmallest(n, iterable, key=None):
  407. """Find the n smallest elements in a dataset.
  408. Equivalent to: sorted(iterable, key=key)[:n]
  409. """
  410. # Short-cut for n==1 is to use min()
  411. if n == 1:
  412. it = iter(iterable)
  413. sentinel = object()
  414. result = min(it, default=sentinel, key=key)
  415. return [] if result is sentinel else [result]
  416. # When n>=size, it's faster to use sorted()
  417. try:
  418. size = len(iterable)
  419. except (TypeError, AttributeError):
  420. pass
  421. else:
  422. if n >= size:
  423. return sorted(iterable, key=key)[:n]
  424. # When key is none, use simpler decoration
  425. if key is None:
  426. it = iter(iterable)
  427. # put the range(n) first so that zip() doesn't
  428. # consume one too many elements from the iterator
  429. result = [(elem, i) for i, elem in zip(range(n), it)]
  430. if not result:
  431. return result
  432. _heapify_max(result)
  433. top = result[0][0]
  434. order = n
  435. _heapreplace = _heapreplace_max
  436. for elem in it:
  437. if elem < top:
  438. _heapreplace(result, (elem, order))
  439. top, _order = result[0]
  440. order += 1
  441. result.sort()
  442. return [elem for (elem, order) in result]
  443. # General case, slowest method
  444. it = iter(iterable)
  445. result = [(key(elem), i, elem) for i, elem in zip(range(n), it)]
  446. if not result:
  447. return result
  448. _heapify_max(result)
  449. top = result[0][0]
  450. order = n
  451. _heapreplace = _heapreplace_max
  452. for elem in it:
  453. k = key(elem)
  454. if k < top:
  455. _heapreplace(result, (k, order, elem))
  456. top, _order, _elem = result[0]
  457. order += 1
  458. result.sort()
  459. return [elem for (k, order, elem) in result]
  460. def nlargest(n, iterable, key=None):
  461. """Find the n largest elements in a dataset.
  462. Equivalent to: sorted(iterable, key=key, reverse=True)[:n]
  463. """
  464. # Short-cut for n==1 is to use max()
  465. if n == 1:
  466. it = iter(iterable)
  467. sentinel = object()
  468. result = max(it, default=sentinel, key=key)
  469. return [] if result is sentinel else [result]
  470. # When n>=size, it's faster to use sorted()
  471. try:
  472. size = len(iterable)
  473. except (TypeError, AttributeError):
  474. pass
  475. else:
  476. if n >= size:
  477. return sorted(iterable, key=key, reverse=True)[:n]
  478. # When key is none, use simpler decoration
  479. if key is None:
  480. it = iter(iterable)
  481. result = [(elem, i) for i, elem in zip(range(0, -n, -1), it)]
  482. if not result:
  483. return result
  484. heapify(result)
  485. top = result[0][0]
  486. order = -n
  487. _heapreplace = heapreplace
  488. for elem in it:
  489. if top < elem:
  490. _heapreplace(result, (elem, order))
  491. top, _order = result[0]
  492. order -= 1
  493. result.sort(reverse=True)
  494. return [elem for (elem, order) in result]
  495. # General case, slowest method
  496. it = iter(iterable)
  497. result = [(key(elem), i, elem) for i, elem in zip(range(0, -n, -1), it)]
  498. if not result:
  499. return result
  500. heapify(result)
  501. top = result[0][0]
  502. order = -n
  503. _heapreplace = heapreplace
  504. for elem in it:
  505. k = key(elem)
  506. if top < k:
  507. _heapreplace(result, (k, order, elem))
  508. top, _order, _elem = result[0]
  509. order -= 1
  510. result.sort(reverse=True)
  511. return [elem for (k, order, elem) in result]
  512. # If available, use C implementation
  513. try:
  514. from _heapq import *
  515. except ImportError:
  516. pass
  517. try:
  518. from _heapq import _heapreplace_max
  519. except ImportError:
  520. pass
  521. try:
  522. from _heapq import _heapify_max
  523. except ImportError:
  524. pass
  525. try:
  526. from _heapq import _heappop_max
  527. except ImportError:
  528. pass
  529. if __name__ == "__main__":
  530. import doctest # pragma: no cover
  531. print(doctest.testmod()) # pragma: no cover