Stoer-Wagner Min-Cut in practice: preparing real graphs, safe contraction rules, recovering the partition, testing against brute force and HCS clustering

The algorithm in one paragraph, and its real cost

Each phase of the algorithm grows a set A from an arbitrary start vertex, always adding the outside vertex with the largest total edge weight into A, the maximum adjacency ordering that resembles Prim's algorithm. When the last vertex t is added, the weight of edges from t to everything else is the cut of the phase, and it is a minimum cut separating t from the vertex s added just before it. The phase records that value and merges s and t into one vertex. After n-1 phases, the smallest recorded value is the global minimum. Cost depends on how you find the most tightly connected vertex:

The stub this page replaces gave the Fibonacci bound for a binary heap; that is wrong, and the difference matters for sparse graphs with millions of edges. Priority queues and heaps covers why the asymptotically better heap usually loses on real hardware.

Preparing a real graph

Real graphs violate the algorithm's assumptions in four routine ways, and each one either crashes your code or silently returns a wrong answer.

• Parallel edges. Two links between the same pair of hosts must be summed into one weight. networkx's stoer_wagner is not implemented for multigraphs and raises; a hand-rolled adjacency dictionary that overwrites instead of adds is silently wrong.
• Self-loops. They never cross a cut. Drop them, or they inflate degrees used by the contraction rules below.
• Negative weights. The correctness argument adds non-negative quantities; with a negative edge the answer can be wrong, and the problem itself becomes much harder. networkx raises. Reject them at load time.
• Disconnected input. The minimum cut of a disconnected graph is 0: one component against the rest. networkx raises on it instead of returning 0, which surprises people running nightly jobs on graphs that are occasionally split. Check connectivity first and short-circuit.

Directed graphs need different algorithms entirely, because a cut's weight depends on direction. If your data is directed but the question is about undirected reliability, symmetrise explicitly and document it.

def normalize(n, edges):
    """Sum parallel edges, drop self-loops, reject negative weights."""
    W = [[0.0] * n for _ in range(n)]
    for u, v, w in edges:
        if w < 0:
            raise ValueError(f"negative weight on ({u}, {v})")
        if u != v:
            W[u][v] += w
            W[v][u] += w
    return W

An implementation that returns the partition

Most library calls give you the cut value, and some give a partition. When you write your own, track which original vertices each merged super-vertex contains; the side of the best cut is then the group of the last-added vertex t in the best phase. The version below uses a linear scan, so it is O(n^3) and very short, which makes it the right reference implementation even if production uses a heap:

def stoer_wagner_dense(W):
    n = len(W)
    if n < 2:
        raise ValueError("need at least two vertices")
    comps = components(W)                  # BFS over W[u][v] > 0
    if len(comps) > 1:
        return 0.0, comps[0]               # disconnected: a free cut
    W = [row[:] for row in W]
    groups = [[v] for v in range(n)]       # originals inside each super-vertex
    active = list(range(n))
    best, best_side = float("inf"), None
    while len(active) > 1:
        conn = {v: 0.0 for v in active}
        order, remaining = [], set(active)
        for _ in range(len(active)):       # maximum adjacency ordering
            u = max(remaining, key=lambda v: conn[v])
            remaining.remove(u)
            order.append(u)
            for v in remaining:
                conn[v] += W[u][v]
        s, t = order[-2], order[-1]
        if conn[t] < best:                 # cut of the phase: t versus the rest
            best, best_side = conn[t], list(groups[t])
        groups[s].extend(groups[t])        # merge t into s
        for v in active:
            W[s][v] += W[t][v]
            W[v][s] = W[s][v]
        W[s][s] = 0.0
        active.remove(t)
    return best, best_side

Ties in the max may be broken arbitrarily without affecting correctness, but they change which minimum cut you get when several exist. With floating-point weights, compare cut values with a tolerance.

Safe contractions before the main loop

On large graphs, most of the time goes into phases that merge vertices which obviously belong together. Padberg and Rinaldi described contraction tests that are safe, meaning they never destroy every minimum cut. Two are cheap enough to run first. Keep an upper bound B, initially the smallest weighted degree, since each single vertex is itself a cut. Then contract an edge (u, v) if either holds:

1. w(u, v) >= B. Any cut separating u from v pays at least w(u, v), so it cannot beat the cut you already hold.
2. 2 w(u, v) >= min(deg u, deg v). Let u be the endpoint with the smaller degree, and take any cut separating them with u on side S, where S has other vertices. Moving u across changes the weight by w(u, S without u) minus w(u, other side), which is at most (deg u - w(u, v)) - w(u, v) = deg u - 2 w(u, v), so never positive. The only cut you cannot move to is u alone, and B already covers it.

After each contraction, re-check the new vertex's degree against B. Then run Stoer-Wagner on whatever is left and expand groups to recover the side:

def reduce_graph(W):
    W = [row[:] for row in W]
    n = len(W)
    groups = [[v] for v in range(n)]
    deg = lambda u: sum(W[u])
    bound, side = min((deg(u), [u]) for u in range(n))   # trivial cuts first
    changed = True
    while changed:
        changed = False
        live = [u for u in range(n) if groups[u]]
        for u, v in itertools.combinations(live, 2):
            w = W[u][v]
            if w > 0 and (w >= bound or 2 * w >= min(deg(u), deg(v))):
                contract(W, groups, u, v)               # merge v into u
                if sum(1 for g in groups if g) > 1:
                    d, x = min((deg(y), y) for y in range(n) if groups[y])
                    if d < bound:
                        bound, side = d, list(groups[x])
                changed = True
                break
    live = [u for u in range(n) if groups[u]]
    return bound, side, [[W[a][b] for b in live] for a in live], [groups[a] for a in live]

def min_cut(W):
    bound, side, Wr, groups = reduce_graph(W)
    if len(Wr) >= 2:
        value, s = stoer_wagner_dense(Wr)
        if value < bound:
            bound, side = value, [v for i in s for v in groups[i]]
    return bound, side

Here contract folds row v into row u and empties v, and components is a plain BFS. This quadratic-scan version is for clarity; production code keeps edge lists and a worklist.

Worked example: two sites and a monitoring host

Build the graph in the figure: sites E and W, each a complete graph on four servers with weight 4 per link, cross links E0-W0 (3) and E3-W2 (2), and a monitoring host M attached to E1 with weight 4. Weighted degrees are 15, 16, 12 and 14 for E0 to E3, 15, 12, 14 and 12 for W0 to W3, and 4 for M. Running stoer_wagner_dense returns 4 with side ['M']; brute force over all 255 bipartitions agrees.

Running the reductions first is instructive. The trivial bound is 4, from M. Rule 1 then contracts every intra-site edge, since 4 >= 4, and the M-E1 edge. That leaves two super-vertices joined by weight 3 + 2 = 5, and 5 >= 4 contracts them too. The graph collapses to one vertex and the answer, 4 with side M, is certified without a single Stoer-Wagner phase.

Now the operational lesson. The question an engineer asked was which links, if cut, would partition the sites. The answer they got was to unplug the monitoring host. Global minimum cut is a correct answer to a different question. To get the site split, you must exclude the trivial cut explicitly: contract M into E1 (rule 2 justifies it, since 2 x 4 >= deg M) without letting its degree set the bound, then rerun. The result is 5 with sides E and W. In general, decide whether single-vertex cuts are answers or noise for your problem, and encode that decision.

Testing an implementation

A min-cut routine can return a plausible number that is wrong, and a partition that does not match the number. Test both properties against an oracle on many small random graphs, including the awkward inputs from the preparation section:

def cut_weight(W, side):
    S = set(side)
    return sum(W[u][v] for u in S for v in range(len(W)) if v not in S)

def brute_force(W):                         # 2^(n-1) - 1 bipartitions
    n = len(W)
    return min(cut_weight(W, [v for v in range(n) if mask >> v & 1])
               for mask in range(1, 2 ** (n - 1)))

rng = random.Random(7)
for _ in range(400):
    n = rng.randint(2, 9)
    edges = [(u, v, rng.choice([0, 1, 2, 3, 5]))
             for u, v in itertools.combinations(range(n), 2) if rng.random() < 0.5]
    edges += [(u, u, 9) for u in range(n) if rng.random() < 0.1]   # self-loops
    edges += edges[: rng.randint(0, 3)]                              # parallel edges
    W = normalize(n, edges)
    value, side = min_cut(W)
    assert 0 < len(side) < n
    assert abs(cut_weight(W, side) - value) < 1e-9     # side matches value
    assert abs(brute_force(W) - value) < 1e-9          # value is minimum

Both the plain and the reduced versions on this page pass this harness. For graphs too big for brute force, a second oracle is n-1 maximum-flow computations from a fixed vertex, using the max-flow min-cut theorem and any solver such as Dinic's algorithm; the smallest s-t cut must equal your answer.

Clustering with Highly Connected Subgraphs

Hartuv and Shamir's Highly Connected Subgraphs (HCS) clustering, published in 2000, uses minimum cut recursively on an unweighted similarity graph. A subgraph on k vertices is highly connected if its edge connectivity exceeds k/2; if not, split it along a minimum cut and recurse on both sides:

def hcs(W, vertices=None):                  # W holds 0/1 adjacency
    if vertices is None:
        vertices = list(range(len(W)))
    if len(vertices) == 1:
        return [vertices]
    sub = [[W[a][b] for b in vertices] for a in vertices]
    value, side = stoer_wagner_dense(sub)
    if value > len(vertices) / 2:
        return [vertices]
    left = [vertices[i] for i in side]
    right = [v for v in vertices if v not in left]
    return hcs(W, left) + hcs(W, right)

On the unweighted version of the example, this returns three clusters: E0 to E3, W0 to W3, and M alone. Each four-server site has connectivity 3, above 4/2. The singleton is typical: min cut peels off weakly attached vertices first, so real uses of HCS add a step that reattaches singletons to their best-connected cluster. If you need balanced parts, global min cut is the wrong objective. Normalised cut, ratio cut or a balanced partitioner such as METIS penalise tiny sides, which min cut rewards.

Choosing a library or algorithm

For edge connectivity of an unweighted graph, the same code applies with every weight 1. For all-pairs minimum cuts, a Gomory-Hu tree from n-1 max-flow runs answers every pair at once, which one global cut cannot. To track connected components while edges are added or removed, look at union-find variants before reaching for repeated cuts.

Failure modes

• Trivial answers. A leaf or low-degree vertex wins, as in the worked example. Decide whether singletons count, and exclude them explicitly if not.
• Library exceptions in batch jobs. A day where the graph is disconnected crashes networkx-based pipelines. Pre-check connectivity.
• Silent edge overwrites. Building adjacency with assignment instead of addition drops parallel capacity and finds cuts that are too cheap.
• Expecting a specific cut. When several cuts tie, tie-breaking decides which you get. Do not build features that assume stability across runs or versions.
• Float noise. Weights from similarity scores produce near-ties; round to a tolerance or compare with one.
• Cubic blow-up. The dense version at 20,000 vertices allocates 400 million matrix cells. Switch to a sparse, heap-based version long before that.

What to do next

1. Normalise every input graph: sum parallel edges, drop self-loops, reject negative weights, short-circuit disconnected graphs.
2. Keep the dense version as your oracle and run the random brute-force harness in CI.
3. Add the two contraction rules and measure the reduction ratio on your real graphs.
4. Write down whether single-vertex cuts are valid answers for your use case.
5. For clustering, prototype HCS and compare it with a balanced partitioner before choosing.
6. Read the proof in the companion article before modifying the phase logic.