Thorup's Algorithm, in depth: linear-time undirected shortest paths with integer weights, the component hierarchy, bucketed visiting and a tested implementation

The claim and the model it needs

Precision about the claim comes first. The graph is undirected and connected, with n vertices and m edges, so m >= n - 1. Every weight is a positive integer that fits in a machine word. The model is the word RAM, where arithmetic, shifts and finding the most significant bit of a word cost constant time. Under those assumptions the algorithm runs in O(m).

Change any assumption and the claim goes away. Directed graphs are not covered. Real weights are not covered in the comparison-addition model, where sorting-style lower bounds apply to anything that outputs vertices in order. Zero weights are not covered, but they are easy to handle: contract each zero-weight component into one vertex first. The 2025 result by Duan and co-authors, an O(m log^2/3 n) algorithm for directed graphs with real weights, is a different line of work. It breaks Dijkstra's sorting barrier in the comparison-addition model and won a STOC 2025 best paper award. For the baseline both improve on, see Dijkstra's algorithm, in depth.

The lemma that removes sorting

Everything rests on one observation. Fix a scale 2^i-1 and look at the components formed by edges lighter than that, the components of edges with weight below 2^i-1. Call them level-(i-1) components. Any edge joining two different level-(i-1) components weighs at least 2^i-1.

Now let several of these components be the children of a level-i component C. Put each child c into a bucket by floor(minD(c) / 2^i-1), where minD(c) is the smallest tentative distance of an unvisited vertex in c. Suppose bucket k is the lowest non-empty one. Any path into child c from a sibling must cross an edge of weight at least 2^i-1. It starts from a vertex whose distance is at least k · 2^i-1. So it arrives with a distance of at least (k + 1) · 2^i-1. Inside c, every vertex whose tentative distance is below that threshold cannot be improved by a sibling.

The consequence is the whole trick. You may finish all of c's vertices below (k + 1) · 2^i-1, recursively, before touching any sibling in the same bucket, and in any order between siblings. Vertices in bucket k are not settled in distance order, and they do not need to be. A bucket index is a shift of an integer, and the next bucket is found by adding one. Nothing is ever compared or sorted.

The component hierarchy

The component hierarchy applies the lemma at every scale at once. Leaves are vertices, at level 0. A level-i node is a component of the subgraph of edges with weight below 2ⁱ, and its children are the level-(i-1) components it contains. Levels where nothing merges are skipped, so every internal node has at least two children and the tree has fewer than 2n nodes.

Building it is Kruskal's algorithm grouped by bit length. An edge of weight w joins components at level w.bit_length(), so group the edges by that value, process the groups in increasing order, and union endpoints with a union-find. After each group, every set of old components that merged becomes a new node at that level. Only minimum-spanning-tree edges ever merge anything, which is why Thorup builds the hierarchy from the MST. For the two building blocks see Kruskal's MST algorithm and union-find, in depth.

A tested implementation

The implementation below keeps the algorithm's structure: the hierarchy, the buckets by shifted minimum distance, the monotone bucket index and the recursive visit with an upper limit. It uses ordinary Python sets and dictionaries for the parts that Thorup replaces with specialised structures. It was checked against a heap-based Dijkstra on 3,000 random connected graphs with weights from 1 to 2⁴⁰, and on a 20,000-vertex graph.

INF = float("inf")

def build_hierarchy(n, edges):
    parent = list(range(n))
    def find(x):
        while parent[x] != x:
            parent[x] = parent[parent[x]]
            x = parent[x]
        return x
    level, children, node_of = [0] * n, [[] for _ in range(n)], list(range(n))
    groups = {}
    for u, v, w in edges:
        if w < 1:
            raise ValueError("weights must be positive integers")
        groups.setdefault(w.bit_length(), []).append((u, v))
    for b in sorted(groups):
        olds = {find(x) for e in groups[b] for x in e}
        for u, v in groups[b]:
            ru, rv = find(u), find(v)
            if ru != rv:
                parent[ru] = rv
        merged = {}
        for r in olds:
            merged.setdefault(find(r), set()).add(r)
        for root, parts in merged.items():
            if len(parts) > 1:               # skip levels where nothing merges
                level.append(b)
                children.append([node_of[r] for r in parts])
                node_of[root] = len(level) - 1
    roots = {node_of[find(x)] for x in range(n)}
    if len(roots) != 1:
        raise ValueError("graph must be connected")
    return level, children, roots.pop()

def thorup_sssp(n, edges, source):
    level, children, root = build_hierarchy(n, edges)
    adj = [[] for _ in range(n)]
    for u, v, w in edges:
        adj[u].append((v, w)); adj[v].append((u, w))
    up = [None] * len(level)
    for p, kids in enumerate(children):
        for k in kids:
            up[k] = p
    D, mind = [INF] * n, [INF] * len(level)
    left = [set(k) for k in children]        # children with unvisited vertices
    bucket, key, ix = [None] * len(level), [None] * len(level), [0] * len(level)

    def lower(x, d):                         # D(x) fell to d: update ancestors
        D[x], a = d, x
        while a is not None and d < mind[a]:
            mind[a] = d
            p = up[a]
            if p is not None and bucket[p] is not None:
                if key[a] is not None:
                    bucket[p][key[a]].discard(a)
                key[a] = d >> (level[p] - 1)
                bucket[p].setdefault(key[a], set()).add(a)
            a = p

    def visit(v, limit):                     # visit v's vertices with D < limit
        if level[v] == 0:                    # a vertex: D[v] is final
            mind[v] = INF
            for x, w in adj[v]:
                if D[v] + w < D[x]:
                    lower(x, D[v] + w)
            return
        s = level[v] - 1                     # children are bucketed by D >> s
        if bucket[v] is None:                # expand on first visit
            bucket[v] = {}
            for c in left[v]:
                if mind[c] < INF:
                    key[c] = mind[c] >> s
                    bucket[v].setdefault(key[c], set()).add(c)
            ix[v] = mind[v] >> s
        while left[v] and (ix[v] << s) < limit:
            b = bucket[v].get(ix[v])
            while b:
                c = b.pop()
                key[c] = None
                visit(c, (ix[v] + 1) << s)
                if mind[c] == INF and not (level[c] and left[c]):
                    left[v].discard(c)       # child finished
                elif mind[c] < INF:
                    key[c] = mind[c] >> s
                    bucket[v].setdefault(key[c], set()).add(c)
            ix[v] += 1                       # buckets only move forward
        mind[v] = min((mind[c] for c in left[v]), default=INF)

    lower(source, 0)
    visit(root, INF)
    return D

Two invariants make it correct. A child is only visited from the parent's current bucket, so the lemma applies. Relaxing an edge can only push a child into the current bucket or a later one, so the bucket index never has to move backwards. When a child returns, its minimum is recomputed and it is re-bucketed in its parent.

Worked example

Take vertices A to F with edges A-B 1, B-C 1, C-D 4, D-E 2, E-F 3, A-F 9 and B-E 6, and source A. The weight-1 edges form X = {A, B, C} at level 1. Edges of weight 2 and 3 have bit length 2, so they form Y = {D, E, F} at level 2. C-D 4 and B-E 6 have bit length 3 and join X and Y into the root Z at level 3. A-F 9 has bit length 4 and merges nothing, so it creates no node.

Z buckets its children by D / 4. X has minD 0 and goes into bucket 0, so Z visits X with limit 4. X's buckets have width 1. It visits A (0), which sets F = 9, so Y moves into Z's bucket 2. Next it visits B (1), which sets C = 2 and E = 7, moving Y to bucket 1. Then it visits C (2), which sets D = 6. X is now finished.

Z advances to bucket 1 and visits Y with limit 8. Y's buckets have width 2. D = 6 and E = 7 both fall in bucket 3, so Y may visit E before D, even though E is farther. That is safe. The only edge between them weighs 2, and 6 + 2 = 8 is not below the bucket's end at 8. Y returns with F = 9 left, which is bucket 2 at Z's scale. Z visits Y again and settles F. The final distances are A 0, B 1, C 2, D 6, E 7 and F 9, which match Dijkstra's.

What linear time actually requires

The teaching version is correct, but not linear. Four pieces of it need replacing to reach O(m), and each replacement explains why the theory is hard to deploy.

• The hierarchy in linear time. Thorup builds the MST in linear time on the word RAM, using Fredman and Willard's integer MST algorithm. He then builds the component tree from the MST in linear time, which avoids the inverse-Ackermann factor of a general union-find.
• Minimum distance of unvisited children. The teaching code walks up the tree on every relaxation and rescans children on every return. Thorup maintains these minima with Gabow's split-findmin structure. The linear bound also relies on Fredman and Willard's atomic heaps, word-RAM priority queues with constant-time operations. Atomic heaps only apply when n is astronomically large, which is the main reason exact implementations are rare.
• Bucket arrays. The index only moves forward, so a node's buckets can be an array indexed relative to its starting value. Its range is bounded by the node's diameter, which is at most the sum of its MST edges, divided by the bucket width. Summed over the hierarchy, that is O(n), because an edge's contribution shrinks geometrically at each ancestor level.
• Most-significant-bit in O(1). w.bit_length() is a single instruction on current CPUs, but on the theoretical word RAM it needs a constant-time construction of its own.

In practice

Measured honestly, the hierarchy idea is elegant and slow. On a random graph with 200,000 vertices and 800,000 edges, the teaching version above took 15.8 s in CPython, against 7.2 s for a binary-heap Dijkstra. Published implementation studies point the same way. Asano and Imai (2000) replaced the theoretical structures with practical ones and evaluated the algorithm mainly in a repeated-query setting. There the hierarchy is built once and reused for queries from many sources, because a single query pays the construction cost that Dijkstra never pays. Pettie and Ramachandran later carried the hierarchy idea to real-weighted undirected graphs, where it helps most when the graph is fixed and queried often.

So the practical rule is this. For one query, use Dijkstra with a good heap, or Dial's buckets when weights are small integers. Consider a hierarchy-based method only when the graph is fixed, the queries come from many sources, and profiling shows the priority queue dominating. For heap choices see priority queues and heaps, in depth.

Failure modes

• Directed edges. The lemma needs the sibling edge to work in both directions. On a directed graph the bucket argument fails, and the algorithm returns wrong distances without any error.
• Zero or negative weights. A zero-weight edge has bit length 0 and would join components at level 0, which breaks the level structure. Contract zero-weight components first. Negative weights are out of scope entirely.
• Disconnected input. The hierarchy has no single root. Run the algorithm on the source's component and report the rest as unreachable.
• Float weights. Shifting a float is meaningless. Scale to integers if you know the precision, or use Dijkstra.
• Overflow in fixed-width ports. Distances can reach (n - 1) times the maximum weight. Size the integer type for that, not for one edge.
• Testing only small weights. Bugs in level arithmetic often hide when every weight is 1 or 2. Fuzz with weights spread over many bit lengths, as the tests behind this article did.

What to do next

1. Run the code above against your own Dijkstra on random graphs with weights from 1 to 2⁴⁰, until you trust both.
2. Trace the six-vertex example by hand, then change C-D to 5 and predict which buckets change.
3. Print the hierarchy for one of your real graphs. Its depth and branching show how much weight-scale structure there is to exploit.
4. If you serve many shortest-path queries on a fixed graph, benchmark a reused hierarchy against Dijkstra before you adopt either.
5. Read Thorup's 1999 JACM paper for the linear-time components, then the 2025 sorting-barrier paper to see how the frontier has moved for directed graphs.