Heapify and heapsort: building a heap in O(n)

Bottom-up heapify, published by Robert Floyd in 1964, turns any list into a heap in O(n) time: it sifts down every parent of the tree, starting with the last parent, at index n // 2 − 1, and ending at the root. Heapsort, published by J. W. J. Williams the same year, then sorts the list in place. Build a max-heap, swap its root, the largest item, with the last item of the heap, shrink the heap by one and sift the new root down; after n − 1 such steps the list is in ascending order. Heapsort takes O(n log n) time on every input and O(1) extra memory, but it is not stable, and on large arrays it is usually slower than quicksort and merge sort, because it jumps around in memory.

Heapify: sift down every parent, last one first

When all the items are already in a list, there is a cheaper way than the one push per item of the previous lesson. Leave the items where they are and repair the heap property from the bottom up, as Python’s heapq.heapify does “in-place, in linear time”, in the words of its documentation:

1. The indices n // 2 to n − 1 are leaves, and a single item is already a heap.
2. For i from n // 2 − 1 down to 0, sift a[i] down: swap it with its larger child (for a max-heap) while that child is larger.
3. When the loop has handled index 0, the whole list is a heap.

The order matters. Sifting down at index i assumes that both subtrees below i are heaps already, so that the larger child is the largest item of its subtree. Going from the last parent towards the root makes that true at every step.

The cost is linear because most items live near the bottom of the tree, where there is little room to sink. Half the items are leaves and never move. The quarter just above the leaves can sink one level each, the eighth above them two, and every level higher holds half as many items with one more level of room. With pushes it is the other way round: the many items at the bottom are the ones that may climb the whole height.

Heapsort: the largest item to the end, n − 1 times

Heapsort uses a max-heap, so that the largest item is at index 0 and can move straight to the end of the list:

1. Heapify the whole list into a max-heap.
2. For end from n − 1 down to 1: swap a[0] with a[end], then sift the new root down inside a[0:end] only.

After each swap the slots from end onwards hold the largest items in ascending order and are never touched again, and the heap part shrinks by one. The code is short:

"""Heapsort on eight numbers: the bottom-up heapify, then every move of the largest item to the sorted tail."""

from heap_tools import sift_down

a = [4, 10, 3, 5, 1, 8, 7, 2]
print(f"start:       {a}")
n = len(a)
for i in range(n // 2 - 1, -1, -1):
    item = a[i]
    swaps = sift_down(a, i, n)
    print(f"sift {item:>2} at {i}: {a}  ({swaps} swap{'' if swaps == 1 else 's'})")
print(f"max-heap:    {a}")
for end in range(n - 1, 0, -1):
    a[0], a[end] = a[end], a[0]
    swaps = sift_down(a, 0, end)
    print(f"{a[end]:>2} to {end}:     {a[:end]} | {a[end:]}  ({swaps} swap{'' if swaps == 1 else 's'})")
print(f"sorted:      {a}")

start:       [4, 10, 3, 5, 1, 8, 7, 2]
sift  5 at 3: [4, 10, 3, 5, 1, 8, 7, 2]  (0 swaps)
sift  3 at 2: [4, 10, 8, 5, 1, 3, 7, 2]  (1 swap)
sift 10 at 1: [4, 10, 8, 5, 1, 3, 7, 2]  (0 swaps)
sift  4 at 0: [10, 5, 8, 4, 1, 3, 7, 2]  (2 swaps)
max-heap:    [10, 5, 8, 4, 1, 3, 7, 2]
10 to 7:     [8, 5, 7, 4, 1, 3, 2] | [10]  (2 swaps)
 8 to 6:     [7, 5, 3, 4, 1, 2] | [8, 10]  (2 swaps)
 7 to 5:     [5, 4, 3, 2, 1] | [7, 8, 10]  (2 swaps)
 5 to 4:     [4, 2, 3, 1] | [5, 7, 8, 10]  (2 swaps)
 4 to 3:     [3, 2, 1] | [4, 5, 7, 8, 10]  (1 swap)
 3 to 2:     [2, 1] | [3, 4, 5, 7, 8, 10]  (1 swap)
 2 to 1:     [1] | [2, 3, 4, 5, 7, 8, 10]  (0 swaps)
sorted:      [1, 2, 3, 4, 5, 7, 8, 10]

"""Bottom-up heapify and heapsort on a Python list, with a max-heap: no child is larger than its parent.

Indices start at 0, so the children of i are 2i + 1 and 2i + 2, and only a[0:n] belongs to the heap.
"""


def sift_down(a, i, n):
    """Moves a[i] down inside a[0:n] until no child is larger. Returns the number of swaps it made."""
    swaps = 0
    while True:
        child = 2 * i + 1
        if child >= n:
            return swaps  # a[i] is a leaf
        if child + 1 < n and a[child + 1] > a[child]:
            child += 1  # the larger child
        if a[child] <= a[i]:
            return swaps
        a[i], a[child] = a[child], a[i]
        swaps += 1
        i = child


def heapify(a):
    """Turns the whole list into a max-heap in place: sift down every parent, the last parent first."""
    return sum(sift_down(a, i, len(a)) for i in range(len(a) // 2 - 1, -1, -1))


def heapsort(a):
    """Sorts a in place in ascending order: heapify, then swap the maximum to the end of the heap n - 1 times."""
    heapify(a)
    for end in range(len(a) - 1, 0, -1):
        a[0], a[end] = a[end], a[0]  # the largest item left joins the sorted tail
        sift_down(a, 0, end)  # the heap is now a[0:end]

Heapify starts at index 3, the last parent of 8 items. 5 is already larger than its only child, 3 sinks below its larger child 8, 10 stays, and 4 sinks twice, past 10 and then past 5: three swaps build the heap. Each later line moves the root into the slot just left of the sorted tail, 10 to index 7, then 8 to index 6, and sifts the item that took its place down again.

Heapify against n pushes, counted

This run builds max-heaps of 2¹⁰ to 2¹⁶ items both ways, from shuffled input and from ascending input, and counts the swaps. The ratio column compares each size with the one before, which is half as large:

"""Swaps made by bottom-up heapify against n pushes into an empty max-heap, for n = 2**10 to 2**16."""

import random

from heap_tools import heapify


def build_by_pushes(items):
    """Pushes the items one at a time (append, then sift up). Returns the heap and the swaps."""
    heap, swaps = [], 0
    for x in items:
        heap.append(x)
        i = len(heap) - 1
        while i > 0 and heap[(i - 1) // 2] < heap[i]:
            heap[i], heap[(i - 1) // 2] = heap[(i - 1) // 2], heap[i]
            swaps += 1
            i = (i - 1) // 2
    return heap, swaps


rng = random.Random(5)
for order in ("random", "ascending"):
    print(f"{order} input: swaps, and the ratio to the size before")
    previous = None
    for k in range(10, 17):
        n = 2**k
        items = list(range(n))
        if order == "random":
            rng.shuffle(items)
        floyd = items[:]
        h = heapify(floyd)
        _, p = build_by_pushes(items)
        ratio = f"  x{h / previous[0]:.2f} and x{p / previous[1]:.2f}" if previous else ""
        print(f"  n = 2^{k} = {n:>6,}: heapify {h:>7,} ({h / n:.2f} n), pushes {p:>9,} ({p / n:.2f} n){ratio}")
        previous = (h, p)

random input: swaps, and the ratio to the size before
  n = 2^10 =  1,024: heapify     770 (0.75 n), pushes     1,294 (1.26 n)
  n = 2^11 =  2,048: heapify   1,537 (0.75 n), pushes     2,640 (1.29 n)  x2.00 and x2.04
  n = 2^12 =  4,096: heapify   3,069 (0.75 n), pushes     5,223 (1.28 n)  x2.00 and x1.98
  n = 2^13 =  8,192: heapify   6,111 (0.75 n), pushes    10,543 (1.29 n)  x1.99 and x2.02
  n = 2^14 = 16,384: heapify  12,254 (0.75 n), pushes    21,074 (1.29 n)  x2.01 and x2.00
  n = 2^15 = 32,768: heapify  24,447 (0.75 n), pushes    42,187 (1.29 n)  x2.00 and x2.00
  n = 2^16 = 65,536: heapify  48,552 (0.74 n), pushes    83,895 (1.28 n)  x1.99 and x1.99
ascending input: swaps, and the ratio to the size before
  n = 2^10 =  1,024: heapify   1,015 (0.99 n), pushes     8,204 (8.01 n)
  n = 2^11 =  2,048: heapify   2,038 (1.00 n), pushes    18,445 (9.01 n)  x2.01 and x2.25
  n = 2^12 =  4,096: heapify   4,085 (1.00 n), pushes    40,974 (10.00 n)  x2.00 and x2.22
  n = 2^13 =  8,192: heapify   8,180 (1.00 n), pushes    90,127 (11.00 n)  x2.00 and x2.20
  n = 2^14 = 16,384: heapify  16,371 (1.00 n), pushes   196,624 (12.00 n)  x2.00 and x2.18
  n = 2^15 = 32,768: heapify  32,754 (1.00 n), pushes   426,001 (13.00 n)  x2.00 and x2.17
  n = 2^16 = 65,536: heapify  65,521 (1.00 n), pushes   917,522 (14.00 n)  x2.00 and x2.15

On shuffled input both stay close to linear: about 0.75 swaps per item for heapify and 1.29 for the pushes. Ascending input is the worst case for pushes into a max-heap, since every new item is the largest so far and climbs to the root: 8 swaps per item at 1,024 items and 14 at 65,536, each doubling costing 2.15 to 2.25 times as much, the signature of n log n. Heapify stayed at 1.00 swap per item. When all the items are there at the start, heapify is the way to build a heap; pushes are for items that arrive one at a time.

Heapsort, quicksort and merge sort compared

Heapsort is the only one of the three that guarantees O(n log n) time with O(1) extra memory, which is why libraries keep it as a safety net inside their quicksort.

Exercise and self-check

Write heapsort yourself, in place, without the built-in sort. The sample tests also count the comparisons, so a quadratic sort fails on 4,096 items.

def heapsort(a):
    """Sorts the list a in ascending order in place, with a max-heap kept in a itself. Returns nothing."""
    # Replace this line: heapify a bottom up, then move the largest item to the end n - 1 times.
    raise NotImplementedError

import math
import random

from heapsort import heapsort


class NoSort(list):
    """A list whose built-in sort refuses to run."""

    def sort(self, *args, **kwargs):
        raise AssertionError("write heapsort yourself instead of calling sort()")


def test_small_lists():
    """sorts the lesson's eight numbers and small edge cases in place"""
    for values in ([4, 10, 3, 5, 1, 8, 7, 2], [], [7], [2, 1], [3, 1, 3, 1, 2], [-5, 0, -5, 9], [1, 2, 3, 4], [9, 7, 5, 3, 1]):
        a = NoSort(values)
        assert heapsort(a) is None, "heapsort returns nothing: it sorts the list it was given"
        assert list(a) == sorted(values), f"{values}: got {list(a)}"


def test_random_lists():
    """200 random lists of up to 60 numbers match sorted()"""
    rng = random.Random(4)
    for _ in range(200):
        values = [rng.randint(-20, 20) for _ in range(rng.randint(0, 60))]
        a = NoSort(values)
        heapsort(a)
        assert list(a) == sorted(values), f"{values}: got {list(a)}"


class Counted:
    """A number that counts every comparison made with it."""

    comparisons = 0

    def __init__(self, v):
        self.v = v

    def __lt__(self, other):
        Counted.comparisons += 1
        return self.v < other.v

    def __le__(self, other):
        Counted.comparisons += 1
        return self.v <= other.v

    def __gt__(self, other):
        Counted.comparisons += 1
        return self.v > other.v

    def __ge__(self, other):
        Counted.comparisons += 1
        return self.v >= other.v


def test_comparisons_on_4096_items():
    """4,096 random numbers take fewer than 2n log2 n + 2n comparisons"""
    rng = random.Random(8)
    n = 4096
    values = [rng.random() for _ in range(n)]
    a = NoSort(Counted(v) for v in values)
    Counted.comparisons = 0
    heapsort(a)
    limit = 2 * n * math.log2(n) + 2 * n
    assert Counted.comparisons < limit, f"{Counted.comparisons:,} comparisons; the limit is {limit:,.0f}"
    assert [x.v for x in a] == sorted(values)

/** Sorts the array a in ascending order in place, with a max-heap kept in a itself. Returns nothing. */
export function heapsort(a) {
  // Replace this line: heapify a bottom up, then move the largest item to the end n - 1 times.
  throw new Error('not written yet');
}

import { test, assert } from 'mysmartcopilot:test';
import { heapsort } from './heapsort.mjs';

const ascending = (a, b) => a - b;
function noSort(values) {
  const a = [...values];
  a.sort = () => { throw new Error('write heapsort yourself instead of calling sort()'); };
  return a;
}

test("sorts the lesson's eight numbers and small edge cases in place", () => {
  for (const values of [[4, 10, 3, 5, 1, 8, 7, 2], [], [7], [2, 1], [3, 1, 3, 1, 2], [-5, 0, -5, 9], [1, 2, 3, 4], [9, 7, 5, 3, 1]]) {
    const a = noSort(values);
    assert.equal(heapsort(a), undefined, 'heapsort returns nothing: it sorts the array it was given');
    assert.deepEqual([...a], [...values].sort(ascending), JSON.stringify(values));
  }
});

test('200 random arrays of up to 60 numbers match a sorted copy', () => {
  let seed = 4;
  const random = () => (seed = (seed * 48271) % 2147483647) / 2147483647;
  for (let k = 0; k < 200; k++) {
    const values = Array.from({ length: Math.floor(random() * 61) }, () => Math.floor(random() * 41) - 20);
    const a = noSort(values);
    heapsort(a);
    assert.deepEqual([...a], [...values].sort(ascending), JSON.stringify(values));
  }
});

test('4,096 random numbers take fewer than 2n log2 n + 2n comparisons', () => {
  // Comparing two of these numbers with <, <=, > or >= reads both values, so the reads are twice the comparisons.
  const box = { reads: 0 };
  let seed = 8;
  const random = () => (seed = (seed * 48271) % 2147483647) / 2147483647;
  const n = 4096;
  const values = Array.from({ length: n }, random);
  const a = noSort(values.map((v) => ({ v, valueOf() { box.reads++; return v; } })));
  heapsort(a);
  const comparisons = box.reads / 2;
  const limit = 2 * n * Math.log2(n) + 2 * n;
  assert.ok(comparisons < limit, `${comparisons} comparisons; the limit is ${limit}`);
  assert.deepEqual(a.map((x) => x.v), [...values].sort(ascending));
});

Why heapify is linear: the invariant, the exact bound and the comparisons Pro: opens with a Pro, Premium or Ultimate pass

Heapsort in depth: its proof, its instability and its speed, measured Pro: opens with a Pro, Premium or Ultimate pass

Heapsort in five languages, the slips that break it, and a practice set Pro: opens with a Pro, Premium or Ultimate pass

Interview questions

Key takeaways

• Bottom-up heapify sifts down every parent from index n // 2 − 1 to 0 and builds a heap in O(n): most items are near the bottom and can barely sink.
• Pushes are for items that arrive one at a time; on ascending input they cost O(n log n), 14 swaps per item at 65,536 items against heapify’s 1.
• Heapsort heapifies into a max-heap, then swaps the root to the end and sifts down n − 1 times: O(n log n) on every input, O(1) extra memory, not stable.
• Among the three classic O(n log n) sorts, heapsort is the one that needs no extra memory and has no bad case, which is why libraries keep it as a fallback for quicksort.