02 - Sorting Algorithms: Bubble, Merge, Quick, Heap

You need to sort a list and the library already has a function that does it, so the sensible question is not "how do you write a sort" but "what does that function do, when does it let me down, and why do textbooks describe four of them". Bubble, merge, quick and heap sort answer the same request with different trade-offs: worst-case time, additional memory, behavior with duplicates. To choose well you have to see the numbers, not memorize a table of acronyms.

This is the second article in the Algorithms and Data Structures series and it picks up the method of the first: every complexity comes with a count of comparisons. I implement the four algorithms in Python 3.14 (python3 -I, standard library only), verify them against sorted and against all permutations up to six elements, and build the worst cases instead of just naming them. The outputs are the ones obtained by running the listings twice with an identical result.

What is measured and how correctness is checked

A comparison sort has a single reference cost: the number of times it compares two keys. Every listing counts that number. The elements are (key, label) pairs and only the key is compared: that way, with repeated keys, the label reveals whether the original order of the duplicates was preserved, that is, whether the algorithm is stable. Python's documentation guarantees that sorted is stable, and that is why I use it as the reference for stability too.

Correctness is not declared, it is proved in the only way that scales: all permutations from 0 to 6 elements (874 inputs in total, covering every possible relative order for small lists) and then 3,000 random lists with many duplicates and a fixed seed. A property that holds on all small inputs is not a proof, but it catches nearly all the index and partition errors you make when writing these algorithms by hand.

Bubble sort: the worst case is the reversed case

Bubble sort walks the list swapping every adjacent pair that is out of order; at each pass the maximum "rises" to the end. The version with early exit stops at the first pass without swaps. The listing verifies it and compares the count on the already sorted case and on the reversed one. The independent check is the number of inversions, the pairs (i, j) with i < j and a greater key, computed by brute force.

# Bubble sort with early exit. Items are (key, label) pairs: only the key is compared.
import itertools, random

def bubble_sort(a):
    a, comparisons, swaps = a[:], 0, 0
    for end in range(len(a) - 1, 0, -1):
        swapped = False
        for i in range(end):
            comparisons += 1
            if a[i][0] > a[i + 1][0]:
                a[i], a[i + 1] = a[i + 1], a[i]
                swaps += 1
                swapped = True
        if not swapped:                      # a pass without swaps: already sorted
            break
    return a, comparisons, swaps

def inversions(keys):                        # independent reference, brute force
    return sum(1 for i in range(len(keys)) for j in range(i + 1, len(keys)) if keys[i] > keys[j])

# correctness: ALL permutations up to n = 6, then 3000 random cases with duplicates (fixed seed)
for n in range(0, 7):
    for p in itertools.permutations(range(n)):
        out, _, sw = bubble_sort([(k, k) for k in p])
        assert [k for k, _ in out] == sorted(p) and sw == inversions(p)
rng = random.Random(1)
for _ in range(3000):
    a = [(rng.randrange(8), t) for t in range(rng.randrange(0, 40))]
    out, _, sw = bubble_sort(a)
    assert out == sorted(a, key=lambda x: x[0])          # sorted() is stable: the labels are compared too
    assert sw == inversions([k for k, _ in a])
print("bubble sort: correct, stable, swaps == inversions (permutations n<=6 and 3000 cases)")

print(f"{'n':>5} {'sorted':>9} {'reversed':>9} {'n(n-1)/2':>9} {'swaps reversed':>15}")
for n in (50, 100, 200, 400):
    _, c0, _ = bubble_sort([(k, k) for k in range(n)])
    _, c1, s1 = bubble_sort([(k, k) for k in range(n, 0, -1)])
    print(f"{n:>5} {c0:>9} {c1:>9} {n * (n - 1) // 2:>9} {s1:>15}")

# Output:
# bubble sort: correct, stable, swaps == inversions (permutations n<=6 and 3000 cases)
#     n    sorted  reversed  n(n-1)/2  swaps reversed
#    50        49      1225      1225            1225
#   100        99      4950      4950            4950
#   200       199     19900     19900           19900
#   400       399     79800     79800           79800

Each swap removes exactly one inversion, which is why the check "swaps == inversions" passes on every input: it is the property that governs the cost. On the reversed case the comparisons are n(n − 1)/2 and the swaps are just as many (400 elements: 79,800 and 79,800); on the sorted case with early exit the comparisons drop to n − 1. Without the early exit the sorted case would still cost quadratic time, and it is the typical mistake of whoever implements it from the textbook. Bubble sort remains a teaching tool: it swaps more than necessary, and the two following algorithms beat it already at a few hundred elements.

Merge sort: divide, sort, merge

Merge sort splits the list in half, sorts the two halves and merges them in a single pass. The merge is the delicate point: it compares the heads of the two halves and takes the smaller. In the listing the condition src[j][0] < src[i][0] uses the strict less-than, so on equal keys the left element wins, and it is this line that makes the algorithm stable. If you replace it with <= the result is still sorted but the order of the duplicates is reversed, and the comparison with sorted flags it.

# Merge sort with an auxiliary buffer: comparison count, BUILT worst case, lower bound
import itertools, math, random

def merge_sort(a):
    comparisons = 0
    def sort(lo, hi, src, buf):              # sorts src[lo:hi] using buf as scratch space
        nonlocal comparisons
        if hi - lo < 2:
            return
        mid = (lo + hi) // 2
        sort(lo, mid, src, buf)
        sort(mid, hi, src, buf)
        i, j, k = lo, mid, lo
        while i < mid and j < hi:
            comparisons += 1
            if src[j][0] < src[i][0]:        # strict '<': on ties the left item wins => stable
                buf[k] = src[j]; j += 1
            else:
                buf[k] = src[i]; i += 1
            k += 1
        buf[k:hi] = src[i:mid] if i < mid else src[j:hi]
        src[lo:hi] = buf[lo:hi]
    out = a[:]
    sort(0, len(out), out, [None] * len(out))
    return out, comparisons

def worst_input(keys):                       # "un-merges": left = even positions, right = odd
    if len(keys) < 2:
        return keys
    return worst_input(keys[0::2]) + worst_input(keys[1::2])

for n in range(0, 7):
    for p in itertools.permutations(range(n)):
        assert [k for k, _ in merge_sort([(k, k) for k in p])[0]] == sorted(p)
rng = random.Random(2)
for _ in range(3000):
    a = [(rng.randrange(8), t) for t in range(rng.randrange(0, 60))]
    assert merge_sort(a)[0] == sorted(a, key=lambda x: x[0])      # includes stability
print("merge sort: correct and stable (permutations n<=6 and 3000 cases)")

print(f"{'n':>5} {'random (mean of 20)':>19} {'built worst case':>19} {'n lg n - n + 1':>15} {'ceil(lg n!)':>12}")
for n in (64, 256, 1024, 2048):
    avg = sum(merge_sort([(k, k) for k in rng.sample(range(n), n)])[1] for _ in range(20)) / 20
    worst = merge_sort([(k, k) for k in worst_input(list(range(n)))])[1]
    bound = math.ceil(math.log2(math.factorial(n)))
    print(f"{n:>5} {avg:>19.0f} {worst:>19} {n * int(math.log2(n)) - n + 1:>15} {bound:>12}")

# Output:
# merge sort: correct and stable (permutations n<=6 and 3000 cases)
#     n random (mean of 20)    built worst case  n lg n - n + 1  ceil(lg n!)
#    64                 306                 321             321          296
#   256                1725                1793            1793         1684
#  1024                8942                9217            9217         8770
#  2048               19942               20481           20481        19581

The worst case is not found by chance: the function worst_input builds it by "un-merging" a sorted list, putting the even positions on the left and the odd ones on the right, recursively, so that every merge has to consume both halves all the way to the end. The comparisons on that input coincide with n lg n − n + 1 for n a power of two (321, 1,793, 9,217, 20,481 in the rows of the table). The random case is lower (19,942 for n = 2,048) but not by much: for merge sort the best and the worst case are within a small factor.

The last column is a limit that no comparison sort can beat. A list of n distinct keys can be in n! orders, and each comparison has two outcomes, so at least ⌈log₂(n!)⌉ comparisons are needed in the worst case: 19,581 for n = 2,048. Merge sort uses at most 20,481, that is 4.6% more; this is why it is called asymptotically optimal. The price is space: the auxiliary buffer is as long as the input, so O(n) of extra memory. Quick sort and heap sort, which have no buffer proportional to the input, show the difference in space at equal average time.

Quick sort: the pivot decides everything

Quick sort chooses an element (the pivot), moves the smaller ones to the left and the larger ones to the right, then sorts the two parts. The cost depends on how balanced the two parts are. The listing uses Lomuto's partition with the pivot in the last position, or with a randomly chosen pivot, and instead of recursion it uses an explicit stack: after each partition it pushes the larger part and continues on the smaller one. This rule limits the stack to O(log n) elements, whatever the input.

# Quick sort (Lomuto partition) with an explicit stack: the larger part is pushed, the smaller one is processed
import itertools, random

def quick_sort(a, pivot="last", rng=None):
    a, comparisons, stack, max_stack = a[:], 0, [(0, len(a) - 1)], 1
    while stack:
        lo, hi = stack.pop()
        while lo < hi:
            if pivot == "random":
                r = rng.randint(lo, hi)
                a[r], a[hi] = a[hi], a[r]
            p, store = a[hi][0], lo
            for i in range(lo, hi):
                comparisons += 1
                if a[i][0] < p:
                    a[i], a[store] = a[store], a[i]
                    store += 1
            a[store], a[hi] = a[hi], a[store]
            left, right = (lo, store - 1), (store + 1, hi)
            if left[1] - left[0] < right[1] - right[0]:   # the larger part goes on the stack
                stack.append(right); hi = store - 1
            else:
                stack.append(left); lo = store + 1
            max_stack = max(max_stack, len(stack))
    return a, comparisons, max_stack

for n in range(0, 7):
    for p in itertools.permutations(range(n)):
        assert [k for k, _ in quick_sort([(k, k) for k in p])[0]] == sorted(p)
rng = random.Random(3)
for _ in range(3000):
    a = [(rng.randrange(8), t) for t in range(rng.randrange(0, 60))]
    for mode in ("last", "random"):
        assert [k for k, _ in quick_sort(a, mode, rng)[0]] == sorted(k for k, _ in a)
stable = quick_sort([(0, t) for t in range(2)] + [(1, 9)] + [(0, 5)])[0]
print("quick sort: correct (permutations n<=6, 3000 cases x 2 pivots). Order of the labels on ties:",
      [t for k, t in stable if k == 0], "(not guaranteed stable)")

print(f"{'n':>5} {'sorted, last pivot':>23} {'n(n-1)/2':>9} {'max stack':>9} | {'sorted, random pivot (mean of 10)':>34} {'random, last pivot':>22} {'2(n+1)Hn-4n':>12}")
for n in (256, 512, 1024, 2048):
    srt = [(k, k) for k in range(n)]
    _, cw, depth = quick_sort(srt)
    cr = sum(quick_sort(srt, "random", rng)[1] for _ in range(10)) / 10
    cavg = sum(quick_sort([(k, k) for k in rng.sample(range(n), n)])[1] for _ in range(10)) / 10
    h = sum(1 / i for i in range(1, n + 1))          # harmonic number H_n
    print(f"{n:>5} {cw:>23} {n * (n - 1) // 2:>9} {depth:>9} | {cr:>34.0f} {cavg:>22.0f} {2 * (n + 1) * h - 4 * n:>12.0f}")

# Output:
# quick sort: correct (permutations n<=6, 3000 cases x 2 pivots). Order of the labels on ties: [5, 0, 1] (not guaranteed stable)
#     n      sorted, last pivot  n(n-1)/2 max stack |  sorted, random pivot (mean of 10)     random, last pivot  2(n+1)Hn-4n
#   256                   32640     32640         1 |                               2181                   2203         2124
#   512                  130816    130816         1 |                               4693                   5006         4946
#  1024                  523776    523776         1 |                              11514                  11244        11298
#  2048                 2096128   2096128         1 |                              25108                  26125        25420

On the already sorted input with the pivot in the last position every partition isolates a single element, and the comparisons are exactly n(n − 1)/2 (2,096,128 for n = 2,048, identical to the formula). It is the real worst case, and it happens with data that arrive already almost sorted. With a random pivot the same input costs about 25,000 comparisons, close to the expected value 2(n + 1)Hₙ − 4n of the fourth column (25,420), where Hₙ is the harmonic number: the worst case still exists, but it can no longer be built in advance by looking at the input.

Two details of the output are worth reading. The maximum stack is 1 even in the quadratic case: the rule "continue on the smaller part" keeps the memory small when the time is not, and as in the first article time and space must be counted separately. And the labels of the duplicates come out in the order 5, 0, 1: quick sort is not stable, and the line shows it on a constructed example (it is not a theorem, it is a proof by counterexample).

Heap sort: no quadratic worst case and no extra memory

Heap sort turns the list into a max heap, that is a binary tree stored in the list in which every node is greater than its children, and then extracts the maximum n − 1 times, moving it to the end. The construction proceeds from the bottom, starting from the last node with children, and for this reason it costs O(n) and not O(n log n): the nodes close to the leaves have little to sink. The listing counts separately the comparisons of the construction and the total ones.

# Heap sort: bottom-up construction (heapify), then n-1 extractions, with comparisons counted separately
import itertools, math, random

def heap_sort(a):
    a, comparisons = a[:], 0
    def sift_down(root, end):                # end: excluded index
        nonlocal comparisons
        while 2 * root + 1 < end:
            child = 2 * root + 1
            if child + 1 < end:
                comparisons += 1
                if a[child][0] < a[child + 1][0]:
                    child += 1               # the larger child
            comparisons += 1
            if a[root][0] < a[child][0]:
                a[root], a[child] = a[child], a[root]
                root = child
            else:
                return
    n = len(a)
    for root in range(n // 2 - 1, -1, -1):   # construction: O(n), not O(n log n)
        sift_down(root, n)
    build = comparisons
    for end in range(n - 1, 0, -1):          # the maximum goes to the end, the heap shrinks
        a[0], a[end] = a[end], a[0]
        sift_down(0, end)
    return a, build, comparisons

for n in range(0, 7):
    for p in itertools.permutations(range(n)):
        assert [k for k, _ in heap_sort([(k, k) for k in p])[0]] == sorted(p)
rng = random.Random(4)
for _ in range(3000):
    a = [(rng.randrange(8), t) for t in range(rng.randrange(0, 60))]
    assert [k for k, _ in heap_sort(a)[0]] == sorted(k for k, _ in a)
unst = heap_sort([(0, 0), (0, 1), (0, 2), (0, 3)])[0]
print("heap sort: correct. Labels of four equal keys, in the output:", [t for _, t in unst], "(not stable)")

print(f"{'n':>5} {'construction':>12} {'2n':>6} {'random total':>15} {'sorted':>9} {'reversed':>9} {'2n lg n':>8}")
for n in (64, 256, 1024, 2048):
    _, b, tot = heap_sort([(k, k) for k in rng.sample(range(n), n)])
    _, _, t_sorted = heap_sort([(k, k) for k in range(n)])
    _, _, t_rev = heap_sort([(k, k) for k in range(n, 0, -1)])
    print(f"{n:>5} {b:>12} {2 * n:>6} {tot:>15} {t_sorted:>9} {t_rev:>9} {2 * n * int(math.log2(n)):>8}")

# Output:
# heap sort: correct. Labels of four equal keys, in the output: [1, 2, 3, 0] (not stable)
#     n construction     2n    random total    sorted  reversed  2n lg n
#    64          104    128             573       593       525      768
#   256          460    512            3317      3452      3106     4096
#  1024         1941   2048           17320     18060     16407    20480
#  2048         3839   4096           38733     40204     36973    45056

The construction costs fewer than 2n comparisons in every row (1,941 against 2,048 for n = 1,024), which confirms the linear bound. The total grows as n log n with a constant close to 2 (38,733 for n = 2,048, against 2n lg n = 45,056), because every sinking step makes two comparisons, one between the children and one with the parent. Unlike quick sort the cost is practically identical for sorted, reversed and random inputs, and no memory is needed beyond the list. The trade-off is instability, visible in the four equal keys that came out in the order 1, 2, 3, 0, and poor locality of access, which I do not measure here.

Heaps are also useful outside sorting, to always extract the minimum or the maximum: you will come back to them with priority queues in the article on linear structures and with trees in the article on BSTs and AVLs.

Which one to choose: the table of properties

• Bubble sort: Θ(n²) in the worst case, O(1) space, stable; n − 1 comparisons on the already sorted input only with the early exit.
• Merge sort: n lg n − n + 1 comparisons in the worst case (n a power of two), O(n) space, stable.
• Quick sort: n(n − 1)/2 in the worst case, on average about 2(n + 1)Hₙ − 4n comparisons (the leading term is 1.39 n lg n) with a random pivot, O(log n) stack with the smaller-part rule, not stable.
• Heap sort: O(n log n) in every case, O(1) additional space, not stable, linear construction.

The rows cannot be ranked in a single order: it depends on which property you need. If the duplicates carry information (sorting records by one field, then by another) you need stability, and only merge sort and bubble sort remain. If memory is the constraint, heap sort and quick sort. If you cannot afford a quadratic worst case, neither quick sort with a fixed pivot nor bubble sort. No algorithm wins on all the rows, and that is why a library sort is almost always a hybrid.

Python's sort: it exploits the order that is already there

Python's sorted function uses Timsort, a stable merge sort that looks in the input for already sorted sequences (runs) and merges them. The implementation notes, in Objects/listsort.txt, state that fewer than lg(N!) comparisons are needed on partially ordered data and "as few as N-1". The listing verifies it with a class that counts the calls to __lt__.

# Python's sort: the comparisons are counted with a class that records them in __lt__
import math, random

class Counted:
    calls = 0
    def __init__(self, k): self.k = k
    def __lt__(self, other):
        Counted.calls += 1
        return self.k < other.k

def count(keys):
    Counted.calls = 0
    sorted(Counted(k) for k in keys)
    return Counted.calls

rng = random.Random(5)
print(f"{'n':>5} {'sorted':>9} {'reversed':>9} {'all equal':>13} {'random (mean of 10)':>19} {'ceil(lg n!)':>12}")
for n in (64, 512, 4096):
    asc = count(range(n))
    desc = count(range(n, 0, -1))
    same = count([7] * n)
    avg = sum(count(rng.sample(range(n), n)) for _ in range(10)) / 10
    print(f"{n:>5} {asc:>9} {desc:>9} {same:>13} {avg:>19.0f} {math.ceil(math.log2(math.factorial(n))):>12}")

# Output:
#     n    sorted  reversed     all equal random (mean of 10)  ceil(lg n!)
#    64        63        63            63                 304          296
#   512       511       511           511                3965         3876
#  4096      4095      4095          4095               43988        43251

Sorted, reversed and all-equal cost exactly n − 1 comparisons (4,095 for n = 4,096): a single scan is enough to recognize a run, even a descending one. On the random input the cost is between 1.7% and 2.7% above the limit ⌈lg n!⌉ (43,988 against 43,251), so on data without structure it loses almost nothing compared to the best possible. It is a result you would not get with any of the four algorithms in textbook form: for real data, which almost always have a bit of order, an adaptive hybrid is the sensible library choice. For the exact guarantees on stability and on key see Python's guide.

When the n log n limit does not hold: counting sort

The limit ⌈lg n!⌉ concerns the sorts that compare elements. If the keys are integers in a known range [0, k) you can sort without comparisons: you count how many times each value appears, compute the position where each value starts and place the elements. The listing counts the operations and compares them with 2n + k.

# Counting sort: no comparisons between elements. It works only if the keys are integers in [0, k)
import random

def counting_sort(a, k):
    ops = 0
    count = [0] * k
    for x in a:                              # 1) histogram: n steps
        count[x] += 1; ops += 1
    pos, total = [0] * k, 0
    for v in range(k):                       # 2) starting positions: k steps
        pos[v] = total; total += count[v]; ops += 1
    out = [None] * len(a)
    for x in a:                              # 3) placement in arrival order: n steps, stable
        out[pos[x]] = x; pos[x] += 1; ops += 1
    return out, ops

rng = random.Random(6)
for _ in range(3000):
    k = rng.randrange(1, 12)
    a = [rng.randrange(k) for _ in range(rng.randrange(0, 50))]
    assert counting_sort(a, k)[0] == sorted(a)
print("counting sort: correct on 3000 cases")

print(f"{'n':>7} {'k':>7} {'operations':>11} {'2n + k':>9}")
for n, k in ((1000, 10), (1000, 1000), (1000, 100000), (100000, 10)):
    a = [rng.randrange(k) for _ in range(n)]
    out, ops = counting_sort(a, k)
    assert out == sorted(a)
    print(f"{n:>7} {k:>7} {ops:>11} {2 * n + k:>9}")

# Output:
# counting sort: correct on 3000 cases
#       n       k  operations    2n + k
#    1000      10        2010      2010
#    1000    1000        3000      3000
#    1000  100000      102000    102000
#  100000      10      200010    200010

The operations coincide with 2n + k in every row: two passes over the input and one over the range. With a thousand keys in [0, 1,000) they are 3,000 operations, far fewer than the at least 8,530 comparisons (⌈lg 1000!⌉) that a comparison sort needs; but with the same thousand keys in [0, 100,000) they become 102,000, because the range pays and not the elements. Placement in arrival order makes it stable, which is why it is used as the base step of radix sort. It works only with bounded integer keys: on long strings or real numbers you go back to comparison sorts.

Exercises

1. Remove the early exit from bubble sort and verify that the already sorted case costs n(n − 1)/2 comparisons.
2. Replace < with <= in the merge of merge sort and find a two-element input on which the comparison with sorted fails.
3. Change the quick sort pivot to "median of three" (first, middle, last) and build an input on which it still costs much more than n lg n.
4. Write a stable heap-based sort, adding the arrival index to the key, and count how many more comparisons it costs.

The first exercise closes the reasoning on bubble sort, the second lets you see stability lost for a single character. The third has no single answer: adversarial inputs exist for every fixed pivot-selection rule, and finding one is the quickest way to understand why we randomize.

Frequently asked questions

Which sorting algorithm is the fastest?

There is none: it depends on the input and on the constraints. In the worst case no comparison sort goes below ⌈log₂ n!⌉ comparisons, and merge sort and heap sort come close to it; in practice, for real data, an adaptive hybrid like Timsort is often the best choice.

Why can quick sort become quadratic?

Because if the pivot always falls near an extreme, every partition eliminates a single element. With the pivot in the last position on an already sorted list the comparisons are n(n − 1)/2; with a random pivot that case can no longer be built in advance.

What does it mean for a sort to be stable?

It means that two elements with the same key stay in the order they were in. It is essential when you sort by several fields in succession; Python's sorted guarantees it, quick sort and heap sort do not.

What you can do now

You can write the four classic sorts and prove their correctness on all small inputs, you can build the worst case of each and compare it with the formula, and you can read a comparison count against the ⌈lg n!⌉ limit. You can also choose based on what you really need, stability, memory or a guarantee on the worst case, and you know that a sort without comparisons exists when the keys are bounded integers. In the next article you look for an element instead: linear search, binary search and hash tables, which presupposes precisely a sorted list.

The same series also includes the closing article on dynamic programming.