Push-Relabel (Goldberg-Tarjan): Preflows, Height Functions & Highest-Label Selection

Preflows, Valid Height Functions & Local Push/Relabel Operations

A preflow $f$ satisfies capacity constraints ($0 \le f(u,v) \le c(u,v)$) and relaxed conservation: for every $u \in V \setminus \{s\}$, the excess flow is non-negative:

A vertex $u \in V \setminus \{s, t\}$ with $e(u) > 0$ is called active (overflowing).

The Two Primitive Operations & Complexity Bounds

• Push($u, v$): Applicable when $u$ is active ($e(u) > 0$), $c_f(u, v) > 0$, and $h(u) = h(v) + 1$ (an admissible downhill edge). Push $\Delta = \min(e(u), c_f(u, v))$ units from $u$ to $v$. A push is saturating if $\Delta = c_f(u, v)$ and nonsaturating otherwise.
• Relabel($u$): Applicable when $u$ is active ($e(u) > 0$) and has no admissible downhill residual edges ($h(u) \le h(v)$ for all $(u, v) \in E_f$). Set $h(u) = 1 + \min_{(u, v) \in E_f} h(v)$.

Since every active node $u$ has a residual path back to $s$ (conservation of initial source push) of length $\le V - 1$ and $h(s) = V$, $h(u)$ can never exceed $2V - 1$. This bounds total relabels to $O(V^2)$, saturating pushes to $O(V E)$, and nonsaturating pushes to $O(V^2 E)$ (or $O(V^2 \sqrt{E})$ under Highest-Label Selection + Gap Heuristic).

Production C++20 Implementation: Highest-Label Push-Relabel with Gap Heuristic

#include <vector>
#include <algorithm>

class PushRelabel {
public:
    struct Arc {
        int to;
        long long cap;
    };

    explicit PushRelabel(int n) : n(n), adj(n), excess(n, 0), h(n, 0), count(2 * n + 1, 0), buckets(2 * n + 1) {}

    void addEdge(int u, int v, long long cap) {
        if (u == v) return;
        adj[u].push_back(static_cast<int>(arcs.size()));
        arcs.push_back({v, cap});
        adj[v].push_back(static_cast<int>(arcs.size()));
        arcs.push_back({u, 0});
    }

    long long maxFlow(int s, int t) {
        h[s] = n;
        excess[t] = 1; // Prevent t from being added to active buckets
        count[0] = n - 1;

        for (int id : adj[s]) {
            push(s, id, arcs[id].cap);
        }

        while (maxHeight >= 0) {
            if (buckets[maxHeight].empty()) {
                --maxHeight;
                continue;
            }
            int u = buckets[maxHeight].back();
            buckets[maxHeight].pop_back();
            discharge(u, s, t);
        }
        return excess[t] - 1;
    }

private:
    int n, maxHeight = 0;
    std::vector<Arc> arcs;
    std::vector<std::vector<int>> adj;
    std::vector<long long> excess;
    std::vector<int> h, count;
    std::vector<std::vector<int>> buckets;

    void push(int u, int id, long long d) {
        if (d == 0) return;
        int v = arcs[id].to;
        if (excess[v] == 0 && h[v] < n) {
            buckets[h[v]].push_back(v);
            maxHeight = std::max(maxHeight, h[v]);
        }
        arcs[id].cap -= d;
        arcs[id ^ 1].cap += d;
        excess[u] -= d;
        excess[v] += d;
    }

    void discharge(int u, int s, int t) {
        while (excess[u] > 0) {
            int minNeighborH = 2 * n;
            for (int id : adj[u]) {
                if (arcs[id].cap == 0) continue;
                int v = arcs[id].to;
                if (h[u] == h[v] + 1) {
                    push(u, id, std::min(excess[u], arcs[id].cap));
                    if (excess[u] == 0) return;
                } else {
                    minNeighborH = std::min(minNeighborH, h[v]);
                }
            }
            // Relabel with Gap Heuristic
            int oldH = h[u];
            --count[oldH];
            if (count[oldH] == 0 && oldH < n) {
                for (int i = 0; i < n; ++i) {
                    if (i != s && i != t && h[i] > oldH && h[i] <= n) {
                        --count[h[i]];
                        h[i] = n + 1;
                    }
                }
            }
            h[u] = minNeighborH + 1;
            ++count[h[u]];
            if (h[u] >= n) return; // Excess will drain back to s in phase 2 or can be ignored for maxFlow value
        }
    }
};