Johnson's Algorithm, in depth: all-pairs shortest paths with negative edges through reweighting

The reweighting idea, and why it is safe

Why not add a constant to every edge? Because paths have different numbers of edges: adding 3 per edge adds 3 to a one-edge path and 9 to a three-edge path, so the shortest path can change. A safe transformation must shift every path between the same two endpoints by the same amount.

Johnson's trick is to give every vertex a potential h(v) and define the reweighted edge as w'(u, v) = w(u, v) + h(u) - h(v). Take any path s = v₀, v₁, ..., v_k = t and sum its new weights. Each inner vertex appears once with a minus sign, as the head of one edge, and once with a plus sign, as the tail of the next, so those terms cancel and the sum telescopes:

w'(path) = sum of (w(v[i-1], v[i]) + h(v[i-1]) - h(v[i]))   for i = 1..k
         = w(path) + h(s) - h(t)

Every s-to-t path shifts by the same amount, h(s) - h(t), so the ordering of s-to-t paths is unchanged and a shortest path under w' is a shortest path under w. The same argument applied to a cycle, where s = t, shows that cycle weights do not change at all. Reweighting cannot remove a negative cycle, which is correct, because a graph with one has no well-defined shortest paths through it.

That works for any h. The remaining question is how to choose h so that every w' is non-negative, which is what Dijkstra needs. Add a new vertex q with a zero-weight edge to every vertex, and let h(v) be the shortest distance from q to v. Shortest distances satisfy the triangle inequality: for every edge (u, v), h(v) <= h(u) + w(u, v), because one way to reach v is to reach u and then take that edge. Rearranged, w(u, v) + h(u) - h(v) >= 0. Every reweighted edge is non-negative, and an edge on some shortest path from q has a reweighted value of exactly zero. Bellman-Ford computes h because it tolerates negative edges, and it detects a negative cycle if one exists. The edges from q make every vertex reachable, so every h(v) is finite and at most zero.

The algorithm and a complete implementation

The full algorithm has four steps:

1. Build G' by adding vertex q and a zero-weight edge from q to every vertex.
2. Run Bellman-Ford from q on G'. If it finds a negative cycle, stop and report it. Otherwise h(v) = dist(q, v).
3. Reweight every original edge: w'(u, v) = w(u, v) + h(u) - h(v). Every w' is at least zero.
4. For each vertex s, run Dijkstra on the reweighted graph to get d'(s, t), and convert back with d(s, t) = d'(s, t) - h(s) + h(t). Pairs that Dijkstra never reaches stay at infinity.

Here is a complete Python implementation. It does not materialise q as a real vertex: initialising every h(v) to 0 is exactly the state after Bellman-Ford's first relaxation of q's zero-weight edges.

import heapq

INF = float("inf")

def johnson(n, edges):
    """n vertices 0..n-1; edges: list of (u, v, w), w may be negative.
    Returns (dist, h) where dist[s][t] is the shortest s->t distance (INF if unreachable).
    Raises ValueError if the graph contains a negative cycle."""
    # Step 1 + 2: Bellman-Ford from a virtual source q with 0-edges to all vertices.
    h = [0] * n                          # dist(q, v) after relaxing q's edges
    for _ in range(n):                   # G' has n + 1 vertices, so n rounds suffice
        changed = False
        for u, v, w in edges:
            if h[u] + w < h[v]:
                h[v] = h[u] + w
                changed = True
        if not changed:
            break
    else:                                # still changing after n rounds
        raise ValueError("negative cycle")

    # Step 3: reweight. Every w2 must be >= 0; assert it, it is a cheap self-check.
    adj = [[] for _ in range(n)]
    for u, v, w in edges:
        w2 = w + h[u] - h[v]
        assert w2 >= 0, (u, v, w2)
        adj[u].append((v, w2))

    # Step 4: Dijkstra from every source, then undo the shift.
    dist = []
    for s in range(n):
        d = [INF] * n
        d[s] = 0
        pq = [(0, s)]
        while pq:
            du, u = heapq.heappop(pq)
            if du > d[u]:
                continue                 # stale heap entry
            for v, w2 in adj[u]:
                if du + w2 < d[v]:
                    d[v] = du + w2
                    heapq.heappush(pq, (d[v], v))
        dist.append([d[t] - h[s] + h[t] if d[t] < INF else INF for t in range(n)])
    return dist, h

The else branch of the for loop runs only if edges were still relaxing after n rounds, which requires a negative cycle. The conversion skips unreachable targets rather than doing arithmetic on infinity.

Worked example: five vertices, two negative edges

Take five vertices and eight directed edges: A→B 4, A→C 2, B→C -3, B→D 2, C→D 3, C→E 1, D→E -2 and E→B 5. There are two negative edges and two cycles, B→C→E→B with weight -3 + 1 + 5 = 3 and B→D→E→B with weight 2 - 2 + 5 = 5. Both are positive, so shortest paths exist.

Potentials. Bellman-Ford from q starts every vertex at 0. B→C lowers C to -3, C→E lowers E to -2, and nothing else improves. So h = (A 0, B 0, C -3, D 0, E -2). No potential is positive, as expected.

Reweighting. Each edge gains h(tail) - h(head):

A→C went from 2 to 5 and the path A→B→C from 1 to 4: both gained h(A) - h(C) = 3, so A→B→C is still shorter. That is the telescoping argument on real numbers.

Dijkstra and conversion. Dijkstra from A on the reweighted graph returns 4 for each of B, C, D and E. Converting with d = d' - h(A) + h(t) gives A→B 4, A→C 4 + (-3) = 1, A→D 4, and A→E 4 + (-2) = 2. Repeating from every source gives the full matrix, which matches Bellman-Ford run separately from each vertex:

Column A is infinite below the diagonal because no edge enters A. That is a real answer, not an error, and your output format must be able to represent it.

Negative cycles

Now lower E→B. The cycle B→C→E→B weighs -3 + 1 + w, so it turns negative for any w below 2. At w = 2 the cycle weighs exactly zero and Johnson's algorithm still succeeds, because a zero cycle never makes a path cheaper. At w = 1 Bellman-Ford keeps relaxing after the last round and the function raises. Both behaviours were checked by running the code above.

Raising is right for a library, but callers often want the cycle itself (keep predecessor pointers, walk back n steps from a vertex relaxed in the extra round, then follow predecessors until one repeats) or the set of pairs whose distance is minus infinity. If negative cycles are an expected input, as in currency arbitrage or constraint feasibility, run the dedicated detection described in negative cycle detection first.

Cost, and when to choose something else

The cost is one Bellman-Ford, O(VE), plus V runs of Dijkstra. With a binary heap each run is O(E log V), so the total is O(VE log V). With a Fibonacci heap each run is O(E + V log V) and the total is O(VE + V² log V), the bound usually quoted, though binary heaps usually win in practice on constant factors. How that compares with the alternatives depends on density:

For a road-like graph with V = 10,000 and E = 40,000, Johnson's algorithm does about VE log V, roughly 5 billion basic steps at log V of about 13, against Floyd-Warshall's 10¹². For a complete graph, E is about V² and Johnson's algorithm becomes O(V³ log V), worse than Floyd-Warshall, whose triple loop is also cache-friendly and easy to vectorise. A rule of thumb: below a few hundred vertices, or above about a quarter of all possible edges, start with Floyd-Warshall and measure.

Memory is the other limit: 10,000 vertices means 10⁸ answers, 800 MB as 64-bit values. After the one Bellman-Ford pass each source is independent, so you can compute rows on demand or stream them to disk.

Implementation pitfalls

Most bugs in Johnson implementations come from arithmetic, not from the algorithm.

• Infinity sentinels. In C++ or Java, using INT_MAX for infinity and then computing d' - h(s) + h(t) overflows for unreachable pairs. Test for the sentinel before converting, as the code above does, and use 64-bit integers for sums of many 32-bit weights.
• Floating-point weights. With doubles, w + h(u) - h(v) can come out as -1e-15 on an edge that should be exactly zero, which breaks the strict assertion and can, rarely, make Dijkstra settle a vertex too early. Clamp tiny negatives to zero with a tolerance tied to the weight scale, or scale weights to integers.
• Forgetting to convert back. Reweighted distances are useful internally, but every number you return must be d' - h(s) + h(t). A unit test on a graph with a negative edge catches this at once; a test on non-negative graphs never will, because h is then all zeros.
• Path reconstruction. Predecessor pointers from each Dijkstra run are valid in the original graph, because the set of shortest paths is the same.
• Parallelism. The V Dijkstra runs share only read-only adjacency and h, so they parallelise with no locking. Bellman-Ford is the serial part; an SPFA queue usually computes h faster than the round-based loop, with the same worst case.

Potentials beyond all-pairs

The potentials idea outlives the all-pairs problem. Min-cost flow by successive shortest paths keeps a potential per vertex and reweights residual edges exactly as Johnson does, so that each augmentation can use Dijkstra even though residual graphs contain negative reverse edges. After each Dijkstra run the potentials are updated by adding the new distances, which keeps every reduced cost non-negative for the next round. A* search is the same transformation in disguise: running Dijkstra with weights w(u, v) - H(u) + H(v) for a heuristic H visits vertices in A*'s order, and H is consistent exactly when those reduced weights are non-negative.

If you already hold a feasible potential, for example from yesterday's run, check every reduced weight in one O(E) pass and skip Bellman-Ford if none is negative. For the single-source building block, see Dijkstra's algorithm in depth; for the shorter introduction, Dijkstra's algorithm.

What to do next

1. Implement the function above and test it on the five-vertex example. Check that h is (0, 0, -3, 0, -2), that every reweighted edge is non-negative, and that row A is (0, 4, 1, 4, 2).
2. Add a randomised test: generate small graphs with some negative edges but no negative cycles, and compare Johnson's output with Bellman-Ford from every source.
3. Lower E→B to 2, then 1, and confirm the zero cycle succeeds and the negative one raises. Then extend the error to return the cycle's vertices.
4. Decide your density threshold: time Johnson's algorithm against Floyd-Warshall on your real graph sizes rather than trusting the asymptotic bounds.
5. If you only query some sources, compute h once and run Dijkstra lazily per source, caching rows.
6. If weights are floats, add a tolerance clamp and a test that fails without it.