Push-Relabel Maximum Flow

The water-and-platforms intuition

Every classic max-flow algorithm before 1988 worked the same way: find a path from source to sink with spare capacity, push as much as it can hold, repeat. Edmonds-Karp finds the shortest such path with BFS; Dinic batches many shortest paths into one blocking flow. They all share one trait — flow moves along complete paths, end to end, one augmentation at a time.

Andrew Goldberg and Robert Tarjan threw that away in 1988. Their idea: forget paths entirely. Picture each vertex as a platform at some height, with reservoirs of water sitting on top. Water can only flow downhill, and only one step at a time, from a platform to a strictly lower neighbor. The source is jacked up to height V and its pipes are flung wide open, so a wall of water — the preflow — slams into its neighbors all at once. From there it's purely local: any platform holding leftover water (excess) pours it onto a lower neighbor. If a platform has excess but every neighbor is at the same height or higher, it relabels — it raises itself just high enough to spill into its lowest neighbor. That's the whole algorithm. Push when you can, relabel when you're stuck.

The miracle is that this purely local, no-global-knowledge process provably converges to a maximum flow. Excess that can never reach the sink eventually climbs back above height V and drains back to the source; everything else ends up at the sink. When no platform has excess, the preflow has become a flow, and it's maximum.

The precise mechanism

On a network with vertices V, edges E, source s and sink t, push-relabel maintains two arrays:

• Height h[v] — an integer label. Initialized to h[s] = V and h[v] = 0 for everyone else.
• Excess e[v] — flow in minus flow out. A vertex is active if e[v] > 0 and it isn't the source or sink.

The two operations obey one invariant — for every residual edge (u,v) with spare capacity, h[u] ≤ h[v] + 1. In words: you can be at most one step above any neighbor you could still push to.

Push(u, v). Applicable when u is active, the residual edge (u,v) has capacity, and h[u] = h[v] + 1. Move δ = min(e[u], residual(u,v)) units across. If δ equals the full residual capacity it's a saturating push (the edge fills up); otherwise it's non-saturating (u's excess hits zero first).

Relabel(u). Applicable when u is active but no admissible push exists — every neighbor reachable in the residual graph is at height ≥ h[u]. Set h[u] = 1 + min{ h[v] : (u,v) has residual capacity }. This raises u by at least one so a push becomes possible next time.

The cost is dominated by counting these operations. There are at most 2V − 1 distinct heights any vertex reaches, so relabels total O(V²). Saturating pushes total O(VE). The expensive bucket is non-saturating pushes: O(V²E) for the generic algorithm. Choosing which active vertex to process next is what tames this — see the comparison below.

When to reach for push-relabel

• Dense graphs where E ≈ V². Push-relabel's local work dominates augmenting-path methods that re-scan the whole network per path.
• Hard adversarial instances — the genrmf, washington, and acyclic-dense generators from the DIMACS max-flow challenge, where augmenting-path counts explode.
• You can afford the heuristics. The gap and global-relabel heuristics are what make push-relabel the fastest practical solver. Bare push-relabel without them is often slower than Dinic.
• Parallel or distributed settings. Because pushes are local, push-relabel parallelizes far more naturally than augmenting-path methods — the asynchronous variant underlies GPU and grid-cut max-flow used in computer-vision image segmentation.

If your graph is unit-capacity or you're doing bipartite matching, prefer Dinic's algorithm at O(E√V). If you just need correct-and-simple, Edmonds-Karp is hard to get wrong.

Push-relabel vs the augmenting-path family

	Push-relabel (highest-label)	Push-relabel (FIFO)	Dinic	Edmonds-Karp	Ford-Fulkerson (BFS-free)
Time complexity	O(V²√E)	O(V³)	O(V²E)	O(VE²)	O(E · max-flow)
Core idea	Local push of excess	Local push of excess	Level graph + blocking flow	Shortest augmenting path	Any augmenting path
Maintains a valid flow throughout?	No — a preflow	No — a preflow	Yes	Yes	Yes
Per-node state	height + excess	height + excess	level only	none	none
Best graph class	Dense / hard instances	General	Unit-capacity, bipartite	Teaching, small graphs	Integer caps, small flow
Heuristics needed for speed	Gap + global relabel	Gap + global relabel	none (current-arc)	none	none
Parallelizable	Yes (asynchronous)	Partly	Hard	Hard	Hard
Terminates on irrational caps?	Yes	Yes	Yes	Yes	Not guaranteed

The headline is the worst-case bound: push-relabel's O(V²√E) beats Edmonds-Karp's O(VE²) outright, and matches or beats Dinic's O(V²E) on dense graphs. But the bound undersells it — with heuristics, push-relabel is the empirical champion in most benchmark suites.

What the numbers actually say

• Heights are bounded by 2V − 1. A vertex with excess that can't reach the sink must route back to the source, which means climbing above height V; it can never exceed 2V − 1. That single fact bounds total relabel work at O(V²) heap-decrement operations.
• The gap heuristic is nearly free and routinely 2–10× faster. Goldberg & Cherkassky's 1997 experimental study showed that combining the gap heuristic with periodic global relabeling cut runtimes by an order of magnitude on the hardest DIMACS families — without changing the asymptotic bound.
• Highest-label vs FIFO is a genuine asymptotic gap. On a sparse graph with V = 10,000 and E = 50,000, the worst-case bounds are V²√E ≈ 2.2 × 10¹⁰ for highest-label against V³ = 10¹² for FIFO — about a 45× difference, and it only widens as the graph grows. In practice both run far below their worst case, but the highest-label rule wins consistently because draining the tallest reservoirs first keeps non-saturating pushes from piling up.
• Memory is O(V + E). You store one residual-capacity int per directed edge, plus height and excess per vertex, plus a current-arc pointer per vertex. No path stack, no level-graph copy.

JavaScript implementation

This is the FIFO variant (process active vertices in queue order), the easiest to read. The graph uses a flat edge list where edge i and its reverse i^1 share capacity bookkeeping.

class MaxFlow {
  constructor(n) {
    this.n = n;
    this.head = Array(n).fill(-1);   // adjacency list head per vertex
    this.to = []; this.cap = []; this.next = [];
  }
  addEdge(u, v, c) {                 // directed cap c, reverse cap 0
    this.to.push(v); this.cap.push(c); this.next.push(this.head[u]); this.head[u] = this.to.length - 1;
    this.to.push(u); this.cap.push(0); this.next.push(this.head[v]); this.head[v] = this.to.length - 1;
  }

  pushRelabel(s, t) {
    const n = this.n, { to, cap, next, head } = this;
    const h = Array(n).fill(0);      // heights
    const e = Array(n).fill(0);      // excess
    const cur = head.slice();        // current-arc pointer per vertex
    h[s] = n;

    // Initial preflow: saturate every edge out of the source.
    for (let i = head[s]; i !== -1; i = next[i]) {
      const c = cap[i];
      if (c > 0) { e[to[i]] += c; e[s] -= c; cap[i] -= c; cap[i ^ 1] += c; }
    }

    const q = [];
    const inq = Array(n).fill(false);
    for (let v = 0; v < n; v++)
      if (v !== s && v !== t && e[v] > 0) { q.push(v); inq[v] = true; }

    while (q.length) {
      const u = q.shift(); inq[u] = false;
      // Discharge u: push then relabel until its excess is gone.
      while (e[u] > 0) {
        if (cur[u] === -1) {                 // no admissible arc — relabel
          let minH = Infinity;
          for (let i = head[u]; i !== -1; i = next[i])
            if (cap[i] > 0) minH = Math.min(minH, h[to[i]] + 1);
          h[u] = minH;
          cur[u] = head[u];
          continue;
        }
        const i = cur[u], v = to[i];
        if (cap[i] > 0 && h[u] === h[v] + 1) {  // admissible: push
          const d = Math.min(e[u], cap[i]);
          cap[i] -= d; cap[i ^ 1] += d; e[u] -= d; e[v] += d;
          if (v !== s && v !== t && !inq[v]) { q.push(v); inq[v] = true; }
        } else {
          cur[u] = next[i];                    // advance current arc
        }
      }
    }
    return e[t];                                // value of the max flow
  }
}

// Example: 6-node network, max flow s=0 → t=5
const g = new MaxFlow(6);
[[0,1,16],[0,2,13],[1,2,10],[2,1,4],[1,3,12],
 [3,2,9],[2,4,14],[4,3,7],[3,5,20],[4,5,4]].forEach(([u,v,c]) => g.addEdge(u,v,c));
console.log(g.pushRelabel(0, 5));   // 23

Two details that bite people. First, the cur[] "current-arc" pointers are not an optimization you can skip — without them, each discharge re-scans all of a vertex's edges and the running time balloons. Second, after a relabel you reset cur[u] to the list head, because the height change may have opened earlier arcs that were previously inadmissible.

Python implementation with the gap heuristic

This version uses the highest-label selection rule plus the gap heuristic — close to what a real solver ships. count[h] tracks how many vertices sit at each height so a gap can be detected in O(1).

from collections import defaultdict

class MaxFlow:
    def __init__(self, n):
        self.n = n
        self.graph = [[] for _ in range(n)]   # each entry: [to, cap, rev_index]

    def add_edge(self, u, v, c):
        self.graph[u].append([v, c, len(self.graph[v])])
        self.graph[v].append([u, 0, len(self.graph[u]) - 1])

    def push_relabel(self, s, t):
        n = self.graph
        N = self.n
        height = [0] * N
        excess = [0] * N
        count = [0] * (2 * N + 1)     # how many vertices at each height
        cur = [0] * N                 # current-arc pointer

        height[s] = N
        count[0] = N - 1
        count[N] = 1

        # Initial preflow out of the source.
        for edge in n[s]:
            v, c, rev = edge[0], edge[1], edge[2]
            if c > 0:
                edge[1] -= c
                n[v][rev][1] += c
                excess[v] += c
                excess[s] -= c

        # Active vertices, bucketed by height for highest-label selection.
        buckets = defaultdict(list)
        max_active = -1                   # -1 means "no active vertex yet"
        for v in range(N):
            if v != s and v != t and excess[v] > 0:
                buckets[height[v]].append(v)
                max_active = max(max_active, height[v])

        def relabel(u):
            old = height[u]
            min_h = 2 * N
            for v, c, _ in n[u]:
                if c > 0:
                    min_h = min(min_h, height[v] + 1)
            count[old] -= 1
            height[u] = min_h
            count[min_h] += 1
            # Gap heuristic: a now-empty height level isolates everything above it.
            if count[old] == 0 and old < N:
                for w in range(N):
                    if old < height[w] < N:
                        count[height[w]] -= 1
                        height[w] = N + 1
                        count[N + 1] += 1

        while max_active >= 0:            # active vertices can sit at height 0 too
            if not buckets[max_active]:
                max_active -= 1
                continue
            u = buckets[max_active].pop()
            if excess[u] == 0 or u in (s, t):
                continue
            while excess[u] > 0:
                if cur[u] == len(n[u]):
                    relabel(u)
                    cur[u] = 0
                    if height[u] > max_active:
                        max_active = height[u]
                    break
                v, c, rev = n[u][cur[u]]
                if c > 0 and height[u] == height[v] + 1:
                    d = min(excess[u], c)
                    n[u][cur[u]][1] -= d
                    n[v][rev][1] += d
                    excess[u] -= d
                    excess[v] += d
                    if v != s and v != t and excess[v] == d:   # just became active
                        buckets[height[v]].append(v)
                        max_active = max(max_active, height[v])
                else:
                    cur[u] += 1
            if excess[u] > 0:                 # still active after relabel
                buckets[height[u]].append(u)
                max_active = max(max_active, height[u])
        return excess[t]

g = MaxFlow(6)
for u, v, c in [(0,1,16),(0,2,13),(1,2,10),(2,1,4),(1,3,12),
                (3,2,9),(2,4,14),(4,3,7),(3,5,20),(4,5,4)]:
    g.add_edge(u, v, c)
print(g.push_relabel(0, 5))   # 23

The gap-heuristic block is the difference between a textbook toy and a competitive solver. When count[old] drops to zero, no vertex remains at height old, so anything taller than old (but below N) is provably cut off from the sink — we lift it straight to N + 1 rather than relabeling it one painful unit at a time.

Variants worth knowing

FIFO push-relabel. Process active vertices in queue order, discharging each fully before re-enqueueing newly-activated neighbors. Clean to implement and provably O(V³).

Highest-label push-relabel. Always discharge the active vertex with the greatest height. This is the O(V²√E) variant and the usual choice for a serious solver. Intuition: draining the highest reservoirs first avoids shoving water up and down repeatedly.

Excess-scaling push-relabel. Only push when excess exceeds a scaling threshold Δ that halves each phase. Ahuja & Orlin's variant achieves O(VE + V²log U) where U is the maximum capacity — strong when capacities are large integers.

The gap and global-relabel heuristics. Not separate algorithms but bolt-ons. Global relabeling periodically recomputes every height as the exact BFS distance to the sink (in the residual graph), correcting the lazy unit relabels. Gap detection lifts isolated levels instantly. Together they're responsible for push-relabel's real-world dominance.

The BK / Boykov-Kolmogorov algorithm. Not push-relabel, but its main rival on the specific grid graphs of computer-vision image segmentation, where it often outruns push-relabel despite a worse worst-case bound. Worth knowing because that's the one domain where the textbook winner sometimes loses.

Common bugs and edge cases

• Forgetting to skip the source and sink as active vertices. Both can carry nonzero "excess" by definition (the source is deeply negative, the sink accumulates the answer), but neither should ever be discharged. Enqueueing the sink is the classic infinite loop.
• Resetting the height instead of decrementing the count. With the gap heuristic you must update count[] on every height change, or the gap detector fires spuriously and corrupts the flow.
• Omitting current-arc pointers. Re-scanning every edge on every discharge silently turns a fast algorithm into a slow one — it still returns the right answer, so the bug hides until a large input times out.
• Initializing h[s] to something other than V. The source must start at exactly V so the preflow can't immediately push back; a wrong value violates the height invariant and the proof of correctness collapses.
• Reading the answer from the wrong place. The max-flow value is excess[t] at termination — not the sum of source out-edges, because some of that preflow flows back when excess can't reach the sink.
• Assuming the preflow is ever a valid flow mid-run. It isn't until the very end. Don't read partial answers or terminate early on a "looks done" check — wait until no active vertex remains.