Edmonds-Karp in Depth: Why Shortest Augmenting Paths Bound Max Flow at O(VE²), With a Traced Proof and an Implementation You Can Trust

Why the choice of path matters

Take vertices s, u, v, t with edges s-u, s-v, u-t and v-t of capacity C, and a middle edge u-v of capacity 1. An unlucky depth-first search finds s-u-v-t and pushes 1 unit. The residual graph now has a reverse edge v-u, so the next search can find s-v-u-t and push 1 more, cancelling the middle flow. Alternate like that and you need 2C augmentations to reach the maximum flow of 2C. The running time depends on the numeric value of the capacities, not on the size of the graph, so it is not polynomial in the input length.

Breadth-first search never makes that mistake here: the two-edge paths s-u-t and s-v-t are shorter than any path through the middle edge, so BFS finds them first and finishes in two augmentations. The deeper point is that path length gives the algorithm a potential function that only moves one way. That monotone quantity is what the proof below exploits, and it is also why the bound holds for real-valued capacities.

An implementation built for the proof

Store edges in flat arrays and add each edge together with its reverse, so edge 2k is the forward edge and 2k+1 its reverse. The partner of edge e is then e XOR 1, which makes the residual update two array writes. Each vertex keeps a list of edge ids, so a BFS touches each edge at most once: O(V + E) per search. The code returns as soon as t is discovered, which is safe because BFS discovers vertices in order of distance, so the recorded path to t is already a shortest one.

from collections import deque

class EdmondsKarp:
    def __init__(self, n):
        self.n = n
        self.adj = [[] for _ in range(n)]   # vertex -> list of edge ids
        self.to, self.cap = [], []          # edge 2k is forward, 2k+1 its reverse

    def add_edge(self, u, v, c):
        self.adj[u].append(len(self.to)); self.to.append(v); self.cap.append(c)
        self.adj[v].append(len(self.to)); self.to.append(u); self.cap.append(0)

    def bfs(self, s, t):
        dist, via = [-1] * self.n, [-1] * self.n
        dist[s] = 0
        q = deque([s])
        while q:
            u = q.popleft()
            for e in self.adj[u]:
                v = self.to[e]
                if self.cap[e] > 0 and dist[v] < 0:
                    dist[v], via[v] = dist[u] + 1, e
                    if v == t:
                        return dist, via
                    q.append(v)
        return dist, None

    def max_flow(self, s, t, trace=None):
        flow = 0
        while True:
            dist, via = self.bfs(s, t)
            if via is None:
                self.reachable = [d >= 0 for d in dist]   # source side of a min cut
                return flow
            path, v = [], t
            while v != s:
                path.append(via[v]); v = self.to[via[v] ^ 1]
            path.reverse()
            b = min(self.cap[e] for e in path)
            for e in path:
                self.cap[e] -= b
                self.cap[e ^ 1] += b
            flow += b
            if trace is not None:
                trace(path, b, dist[t], flow, self)

Three details carry more weight than they look. The capacity array is the residual graph, so there is no separate flow array; the flow on original edge 2k is cap[2k+1], the capacity accumulated on its reverse. For an undirected edge, give both directions capacity c in the same pair rather than adding two pairs. And the final BFS, the one that fails to reach t, is not wasted work: the vertices it reached are exactly the source side of a minimum cut.

Worked example: a reverse edge earns its keep

Six vertices and nine edges, added in this order: s-a 9, a-t 7, s-d 8, a-b 1, c-d 5, a-c 6, d-b 6, c-t 3, b-t 2. Running the class above with a trace callback that prints each path, marking reverse edges (odd ids) with an asterisk, produces the following. The table is the program's output, not a hand calculation.

Follow vertex a through the run. At augmentation 2 the path s-a-b-t saturates a-b, so a is at distance 1 and b at distance 2. By augmentation 5 the only route left from s is through d: s-d-b reaches b at distance 2, and the one way onward from b is the reverse residual edge b-a, created when augmentation 2 pushed a unit along a-b. Taking it cancels that unit: the flow that went a-b-t is rerouted, a now sends it through c instead, and b's capacity to t is used by the flow arriving from d. Vertex a's distance has grown from 1 to 3 while the path length jumped from 3 to 5.

After augmentation 5, BFS from s reaches only b and d: s-a is saturated, b-t is saturated, and the reverse edge b-a has been used up. The original edges leaving {s, b, d} are s-a and b-t, capacity 9 + 2 = 11, equal to the flow.

The proof, step by step

Write d_f(v) for the BFS distance from s to v in the residual graph of the current flow f. The whole argument rests on one lemma and one counting trick.

Lemma 1: distances never decrease. After augmenting along a shortest path, d(v) does not decrease for any vertex v. Augmenting changes the residual graph in two ways: it removes the critical edges of the path and adds reverse edges v-u for path edges u-v. Every path edge satisfies d(v) = d(u) + 1, so every new reverse edge goes from a vertex at distance k + 1 back to one at distance k. An edge that steps backwards in distance cannot create a shortcut; a new shorter path would need an edge that jumps forward by two or more levels, and augmentation never creates one. Formally, take the closest vertex whose distance decreased and derive a contradiction.

Lemma 2: each edge is critical at most V/2 times. Call a residual edge u-v critical when it has the smallest residual capacity on the augmenting path, so the augmentation removes it. When that happens, d(v) = d(u) + 1. For u-v to reappear, some later augmentation must push flow along v-u, which requires d'(u) = d'(v) + 1 at that moment. By Lemma 1, d'(v) ≥ d(v), so d'(u) ≥ d(u) + 2. Each time u-v is critical again, u's distance has grown by at least 2, and a distance is at most V - 2 while u is reachable at all, so u-v is critical at most about V/2 times. The worked example shows the mechanism: a-b went critical with d(a) = 1, and before the reverse edge b-a could be used, d(a) had grown to 3.

Counting. Every augmentation has at least one critical edge, and there are at most 2E residual edges, each critical at most V/2 times. So the number of augmentations is O(VE). Each costs one BFS plus a path walk, O(V + E), which is O(E) when every vertex is reachable. Total: O(VE²). Nothing in the argument mentions capacity values, which is why the bound holds for real capacities. See Big-O notation if the asymptotic bookkeeping is unfamiliar.

What the old version of this page got wrong

The earlier version of this page gave an adjacency-matrix implementation and then quoted O(VE²) for it. The two do not go together. A BFS over a V by V capacity matrix scans a full row for every dequeued vertex, costing O(V²) per search regardless of how many edges exist. With O(VE) augmentations that code runs in O(V³E), which on a sparse graph with E close to V is worse than O(VE²) by a factor of about V. The matrix form is fine for contest-sized dense graphs of a few hundred vertices and wastes memory beyond that: a hundred thousand vertices would need ten billion cells. It also cannot represent parallel edges without merging them. The edge-list form above has the claimed complexity and uses O(V + E) memory.

Special cases worth knowing

Unit capacities and bipartite matching. Every augmentation adds at least one unit, so the number of augmentations is at most the maximum flow value. For bipartite matching the flow is at most min(|L|, |R|) ≤ V, which gives O(VE) directly: better than the general bound, and one reason Edmonds-Karp is a reasonable first matching solver. Dinic does better still on those graphs, O(E sqrt V), matching Hopcroft-Karp.

Integrality. With integer capacities every bottleneck is an integer, so the final flow is integral and decomposes into integral paths. That is what lets you read assignments out of a matching network.

The other rule in the 1972 paper. Edmonds and Karp also analysed augmenting along the path with the largest bottleneck, found with a Dijkstra-style search on a priority queue that maximises the minimum capacity. With integer capacities that variant needs O(E log F*) augmentations, where F* is the maximum flow value, so its bound depends on the capacities; it can beat shortest paths when capacities differ by orders of magnitude. The O(VE) augmentation bound belongs to the shortest-path rule only.

Failure modes and engineering pitfalls

• Overflow. The flow value can reach the sum of the capacities leaving s. In Java or C++ with 32-bit int capacities that sum overflows silently; use 64-bit for flow and for the initial bottleneck sentinel.
• Forgetting reverse edges. Without them the algorithm is a greedy path packer and returns less than the maximum, often only slightly less, so tests with hand-picked answers may still pass. Test on graphs, like the one above, where the optimum requires cancellation.
• Reading the cut from the wrong graph. The source side is the set reachable in the final residual graph. Report the original edges leaving that set; reporting saturated edges is not the same thing, because a saturated edge can have both endpoints on the same side.
• Reusing a solved graph. max_flow consumes capacities. To solve again with a new source, sink or capacity, rebuild or keep a copy of the original capacity array.
• Floating-point capacities. Bottlenecks like 1e-17 keep BFS finding paths that make no real progress. Treat residuals below an epsilon as zero, or scale to integers.

Choosing an algorithm

Dinic runs the same BFS but keeps its level graph and pushes along every shortest path of that length before searching again. If Edmonds-Karp is the bottleneck in your profile, Dinic is the smallest change that fixes it. Measure first: count augmentations and BFS calls on your real inputs, because the worst-case bound is rarely approached and a readable Edmonds-Karp on a graph of a few thousand edges finishes in milliseconds.

What to do next

1. Type in the class above, run it on the nine-edge network, and check that your trace matches the table, including the reverse edge in augmentation 5.
2. Add counters for augmentations and BFS calls, run your real graphs, and compare the counts with VE before deciding the algorithm is too slow.
3. Write a brute-force checker: on small random graphs, compare the returned flow with the capacity of the cut read from reachable, and assert they are equal.
4. Replace any adjacency-matrix implementation on sparse inputs with the paired edge-list form and 64-bit flow values.
5. Model a bipartite assignment problem as a flow network and extract the matching from edges whose reverse holds one unit.
6. When you need more speed, implement Dinic by reusing the BFS and adding a level graph with per-vertex edge pointers.