The Hungarian Algorithm

The assignment problem, and why brute force fails

You have three taxis and three waiting passengers. Each taxi-to-passenger pairing has a cost — the distance the cab must drive to the pickup. You must assign every taxi to exactly one passenger and every passenger to exactly one taxi, and you want the total driving distance to be as small as possible. That is the assignment problem, and the Hungarian algorithm solves it exactly.

Formally, you are given an n × n cost matrix C, where C[i][j] is the cost of assigning worker i to job j. You must pick exactly one entry in each row and each column — a permutation — so that the sum of the chosen entries is minimized. In graph terms, the workers and jobs form the two sides of a complete bipartite graph, the edge weights are the costs, and you want the minimum-cost perfect matching.

The naive approach is to try every permutation. For a 3×3 matrix that is only 6 candidates, but the factorial explodes fast: a 10×10 matrix has 3.6 million permutations, a 15×15 matrix has 1.3 trillion, and a 20×20 matrix has 2.4 × 10¹⁸ — more than a modern CPU could enumerate in a century. The Hungarian algorithm finds the same provably optimal answer in O(n³): a 20×20 instance is solved in a few thousand operations, instantly.

The key insight: reductions preserve the optimum

The whole method rests on one observation. If you subtract a constant from every entry of a single row, the cost of every possible assignment drops by exactly that constant. Why? Every perfect matching picks exactly one cell from that row, so it pays the subtracted amount exactly once. Subtracting from a column works the same way. These shifts change the total cost of every matching identically, so they never change which matching is cheapest — they only change the numbers.

So we are free to push zeros into the matrix. The algorithm does this in two opening moves:

• Row reduction. Subtract each row's minimum from that row. Now every row contains at least one zero.
• Column reduction. Subtract each column's minimum from that column. Now every column contains at least one zero too.

A zero in the reduced matrix means "this assignment is free, given the discounts we've handed out." If we can pick n zeros, one per row and one per column — a perfect matching made entirely of zeros — that matching costs 0 in the reduced matrix, which is the smallest possible, so it is optimal in the original matrix too. The whole game is to manufacture enough independent zeros that such a matching exists.

The covering step and how it creates new zeros

After reduction, we test whether a zero-only perfect matching exists. The classic test (König's theorem) is: cover all the zeros with the minimum number of horizontal and vertical lines. If the minimum number of covering lines equals n, a perfect matching of zeros exists and we are done. If it is fewer than n, there are not yet enough independent zeros, and we must create more without destroying the optimum.

To create more zeros, look at every entry not covered by any line and find the smallest one, call it δ. Then:

• Subtract δ from every uncovered entry — this drives at least one new entry to zero.
• Add δ to every entry covered by two lines (a row line and a column line crossing) — this keeps those entries non-negative and consistent.
• Entries covered by exactly one line are unchanged.

This is exactly equivalent to a row-and-column reduction in disguise, so the optimum is preserved. Repeat the cover-and-adjust loop. Each pass strictly increases the number of covering lines required, so after at most n iterations the line count hits n and an optimal assignment of zeros is read off the matrix.

The fast O(n³) implementation reframes this with potentials: a value u[i] for each row and v[j] for each column, where the reduced cost is C[i][j] − u[i] − v[j] ≥ 0. Each of the n phases adds one worker to the matching by finding an augmenting path of tight (zero reduced-cost) edges, maintaining a slack array so each phase costs O(n²). That is the duality view: the potentials are the optimal solution to the dual linear program, and complementary slackness proves optimality.

When to reach for the Hungarian algorithm

• One-to-one assignment with weights. Workers to tasks, surgeons to operating rooms, runners to relay legs — anytime each thing on the left must pair with exactly one thing on the right and you care about total cost.
• Multi-object tracking. The standard data-association step in radar and computer-vision trackers (SORT, the classic Kalman-filter tracker) is a Hungarian solve: match this frame's detections to last frame's tracks by predicted distance.
• Small-to-medium dense problems. O(n³) is ideal up to a few thousand on each side. A 1000×1000 instance is a fraction of a second.
• When you need the proven optimum, not a heuristic. Unlike greedy matching or simulated annealing, the Hungarian algorithm returns the exact minimum every time.

Reach elsewhere when the structure differs: for huge sparse bipartite graphs, min-cost max-flow with successive shortest paths can be faster; for unweighted matching, Hopcroft–Karp is the right tool; for non-bipartite (general) graphs, you need Edmonds' blossom algorithm.

Hungarian vs other matching and assignment methods

	Hungarian	Brute force	Greedy	Hopcroft–Karp	Min-cost max-flow	Auction algorithm
Problem solved	Min-cost perfect matching	Min-cost perfect matching	Approximate matching	Max cardinality (unweighted)	Min-cost flow (general)	Min-cost assignment
Time complexity	O(n³)	O(n · n!)	O(n² log n)	O(E√V)	O(V²E) typical	O(n³) (ε-scaling)
Returns optimum?	Yes	Yes	No	Yes (cardinality)	Yes	Yes
Handles weights?	Yes	Yes	Yes	No	Yes	Yes
Sparse-friendly?	No (dense matrix)	No	Yes	Yes	Yes	Yes
Best when	Dense n ≤ few thousand	n ≤ 10, teaching	Need a fast estimate	No weights	Sparse / capacitated	Parallel / large n

The headline contrast is with brute force: both return the exact optimum, but the Hungarian algorithm trades factorial time for cubic time. Against min-cost max-flow, the Hungarian method is the specialized, lower-overhead solver for the specific case of a square dense cost matrix; min-cost flow is the generalization that also handles capacities and unequal side sizes.

What the numbers actually say

• Brute force is hopeless past n ≈ 12. A 12×12 matrix is 479 million permutations; at 100 million scored permutations per second that is ~5 seconds. A 13×13 is over a minute. A 15×15 is several hours, and a 16×16 stretches into days. The Hungarian algorithm solves all of these in microseconds.
• O(n³) is concrete. A 1000×1000 dense instance is on the order of 10⁹ operations — well under a second in C++. A 100×100 instance is ~10⁶ operations, effectively instant.
• The history step-down: O(n⁴) → O(n³). Kuhn and Munkres' original analysis was O(n⁴). The slack-and-potential formulation, widely attributed to refinements in the 1970s, brought it to O(n³) — the version used everywhere today, including SciPy's linear_sum_assignment.
• It is a polynomial-time solution to an LP-hard-looking problem. The assignment problem is an integer program, yet its constraint matrix is totally unimodular, so its LP relaxation has integer optimal vertices — which is the deep reason a combinatorial O(n³) method can reach the exact integer optimum.

JavaScript implementation

This is the O(n³) potential-based version (the Jonker–Volgenant / Kuhn–Munkres formulation). It assumes a square cost matrix and returns, for each column, the row assigned to it.

// Minimum-cost perfect matching of an n×n cost matrix.
// Returns assignment[j] = row matched to column j, and the total cost.
function hungarian(cost) {
  const n = cost.length;
  const INF = Infinity;
  // 1-indexed potentials and helper arrays (index 0 is a virtual start).
  const u = new Array(n + 1).fill(0);   // row potentials
  const v = new Array(n + 1).fill(0);   // column potentials
  const p = new Array(n + 1).fill(0);   // p[j] = row currently matched to column j
  const way = new Array(n + 1).fill(0); // back-pointers for augmenting path

  for (let i = 1; i <= n; i++) {
    p[0] = i;                 // start a phase for worker i via virtual column 0
    let j0 = 0;
    const minv = new Array(n + 1).fill(INF);
    const used = new Array(n + 1).fill(false);

    do {                      // grow the shortest augmenting path
      used[j0] = true;
      const i0 = p[j0];
      let delta = INF, j1 = -1;
      for (let j = 1; j <= n; j++) {
        if (used[j]) continue;
        const cur = cost[i0 - 1][j - 1] - u[i0] - v[j];
        if (cur < minv[j]) { minv[j] = cur; way[j] = j0; }
        if (minv[j] < delta) { delta = minv[j]; j1 = j; }
      }
      for (let j = 0; j <= n; j++) {       // re-weight by the slack delta
        if (used[j]) { u[p[j]] += delta; v[j] -= delta; }
        else { minv[j] -= delta; }
      }
      j0 = j1;
    } while (p[j0] !== 0);    // stop at a free column

    do {                      // flip edges along the augmenting path
      const j1 = way[j0];
      p[j0] = p[j1];
      j0 = j1;
    } while (j0);
  }

  const assignment = new Array(n);   // assignment[j] = matched row (0-indexed)
  let total = 0;
  for (let j = 1; j <= n; j++) {
    assignment[j - 1] = p[j] - 1;
    total += cost[p[j] - 1][j - 1];
  }
  return { assignment, total };
}

// Example: 3 taxis × 3 passengers, costs in km
const C = [
  [4, 1, 3],
  [2, 0, 5],
  [3, 2, 2],
];
console.log(hungarian(C));
// { assignment: [ 1, 0, 2 ], total: 5 }
// assignment[col] = row, so col0→row1, col1→row0, col2→row2 (total cost 5)

The structure mirrors the math directly: the inner do/while grows a shortest augmenting path of tight edges using the minv slack array, the re-weighting loop shifts potentials by delta so a new edge becomes tight, and the final flip reverses the matching along the path. n phases, each O(n²), gives the O(n³) bound.

Python implementation (and the famous one-liner)

In practice almost nobody hand-writes the Hungarian algorithm in production Python — SciPy ships a battle-tested O(n³) implementation:

import numpy as np
from scipy.optimize import linear_sum_assignment

C = np.array([[4, 1, 3],
              [2, 0, 5],
              [3, 2, 2]])

rows, cols = linear_sum_assignment(C)   # the Hungarian algorithm
print(list(zip(rows, cols)))            # [(0, 1), (1, 0), (2, 2)]
print(C[rows, cols].sum())              # 5  → minimum total cost

# Maximization? Negate the matrix, or subtract from the max.
profit = np.array([[7, 9, 2],
                   [6, 4, 3],
                   [5, 1, 8]])
r, c = linear_sum_assignment(profit, maximize=True)
print(profit[r, c].sum())               # 23  → maximum total profit

For learning, here is a compact self-contained version using the same potential method as the JavaScript above:

def hungarian(cost):
    n = len(cost)
    INF = float('inf')
    u = [0] * (n + 1)          # row potentials
    v = [0] * (n + 1)          # column potentials
    p = [0] * (n + 1)          # p[j] = row matched to column j
    way = [0] * (n + 1)

    for i in range(1, n + 1):
        p[0] = i
        j0 = 0
        minv = [INF] * (n + 1)
        used = [False] * (n + 1)
        while True:
            used[j0] = True
            i0, delta, j1 = p[j0], INF, -1
            for j in range(1, n + 1):
                if used[j]:
                    continue
                cur = cost[i0 - 1][j - 1] - u[i0] - v[j]
                if cur < minv[j]:
                    minv[j], way[j] = cur, j0
                if minv[j] < delta:
                    delta, j1 = minv[j], j
            for j in range(n + 1):
                if used[j]:
                    u[p[j]] += delta
                    v[j]   -= delta
                else:
                    minv[j] -= delta
            j0 = j1
            if p[j0] == 0:
                break
        while j0:
            j1 = way[j0]
            p[j0] = p[j1]
            j0 = j1

    assignment = [p[j] - 1 for j in range(1, n + 1)]   # column -> row
    total = sum(cost[p[j] - 1][j - 1] for j in range(1, n + 1))
    return assignment, total

print(hungarian([[4, 1, 3], [2, 0, 5], [3, 2, 2]]))   # ([1, 0, 2], 5)

Variants worth knowing

Jonker–Volgenant (LAPJV). A 1987 refinement that adds a clever column-reduction and shortest-augmenting-path initialization. Same O(n³) worst case but dramatically faster constants on real data; it's what many computer-vision libraries actually use.

Auction algorithm. Bertsekas' market-based method where "workers" bid for "jobs" and prices rise until a stable assignment emerges. With ε-scaling it is O(n³) and parallelizes far better than the Hungarian method, making it the choice for very large or distributed assignment.

Rectangular / unbalanced assignment. When the sides are unequal (say m workers, n jobs, m ≠ n), pad to a square with dummy rows/columns of zero cost, or use a variant of LAPJV that handles rectangular matrices directly without the padding overhead.

Min-cost max-flow generalization. Model workers and jobs as a flow network with a source, a sink, unit capacities, and edge costs equal to the assignment costs. Min-cost max-flow then solves assignment as a special case — and also handles capacities, multiple jobs per worker, and sparse graphs the dense Hungarian method can't.

Edmonds' blossom for general graphs. The Hungarian algorithm is strictly for bipartite graphs. If left and right are not distinct sets — e.g., pairing players for a tournament where anyone can play anyone — you need the blossom algorithm, which contracts odd cycles ("blossoms") to extend augmenting-path matching to general graphs.

Common bugs and edge cases

• Forgetting to make the matrix square. The algorithm needs a perfect matching to exist. With unequal sides you must pad with dummy rows or columns first; skipping this is the most common source of wrong answers.
• Maximizing by reusing the minimization code unchanged. You must negate the costs (or subtract from the matrix max) before minimizing. Feeding raw profits into a minimizer finds the worst assignment.
• Off-by-one in the 1-indexed implementation. The potential-based version uses a virtual column 0 and 1-indexed arrays; mixing those with 0-indexed matrix access (cost[i0-1][j-1]) is where most hand-ported bugs hide.
• Assuming the optimal assignment is unique. Ties produce multiple optimal matchings with the same total cost. If you need a deterministic result, break ties consistently (e.g., prefer lower column index).
• Using floats with exact-equality zero tests. The "is this a zero?" and "is the reduced cost tight?" checks assume exact arithmetic. With floating-point costs, accumulated error can mislabel a near-zero slack; use a tolerance or integer-scaled costs.
• Reaching for it on a huge sparse graph. A dense O(n³) matrix solve on a graph that is 99% non-edges wastes enormous time and memory. Use min-cost max-flow or the auction algorithm on sparse instances instead.