Weighted matching

The **assignment problem** asks for a perfect matching between n rows and n columns of minimum total cost. The **Hungarian algorithm** solves it by inserting the rows one at a time. Each insertion grows an alternating tree from the new row using only edges whose **reduced cost** c[i][j] − u[i] − v[j] is zero (tight edges). When the tree cannot reach a free column, the algorithm lowers the reduced costs along the tree's boundary by the smallest slack delta, which creates at least one new tight edge, and continues.

The implementation below is the classic O(n^3) version with 1-based potentials u (rows) and v (columns). minv[j] is the smallest reduced cost from the tree to column j, and way[j] remembers which tree column gave that minimum. When a free column is reached, walking back through way[] hands every column along the path to the previous row, which is an augmentation.

Open the animation with the 3 × 3 matrix 4 1 3 / 2 0 5 / 3 2 2. Each cell shows its reduced cost with the original cost underneath. Rows are inserted one at a time (3 dual updates in total here), and the final assignment is row 0 to column 1, row 1 to column 0, row 2 to column 2, costing 1 + 2 + 2 = 5. The last message also reports that the dual variables sum to 5; the reason that matters is explained in the dual-variables section.

1static int[] hungarian(int[][] cost) {                                // n x n costs; returns column assigned to each row
2    int n = cost.length;
3    int[] u = new int[n + 1], v = new int[n + 1], p = new int[n + 1], way = new int[n + 1];   // 1-based; u, v: dual variables
4    for (int i = 1; i <= n; i++) {                                    // insert row i
5        p[0] = i;
6        int j0 = 0;                                                   // column 0 is a dummy that holds the new row
7        int[] minv = new int[n + 1];
8        boolean[] used = new boolean[n + 1];
9        Arrays.fill(minv, Integer.MAX_VALUE);
10        do {
11            used[j0] = true;
12            int i0 = p[j0], delta = Integer.MAX_VALUE, j1 = 0;
13            for (int j = 1; j <= n; j++) {
14                if (used[j]) continue;
15                int cur = cost[i0 - 1][j - 1] - u[i0] - v[j];         // reduced cost (slack) of edge i0 - j
16                if (cur < minv[j]) { minv[j] = cur; way[j] = j0; }
17                if (minv[j] < delta) { delta = minv[j]; j1 = j; }
18            }
19            for (int j = 0; j <= n; j++) {                            // dual update: make the cheapest slack tight
20                if (used[j]) { u[p[j]] += delta; v[j] -= delta; }
21                else minv[j] -= delta;
22            }
23            j0 = j1;
24        } while (p[j0] != 0);                                         // stop at a column nobody owns yet
25        do { int j1 = way[j0]; p[j0] = p[j1]; j0 = j1; } while (j0 != 0);   // augment along the way[] chain
26    }
27    int[] result = new int[n];
28    for (int j = 1; j <= n; j++) result[p[j] - 1] = j - 1;
29    return result;
30}

1static int[] hungarian(int[][] cost) {                                // n x n costs; returns column assigned to each row
2    int n = cost.length;
3    int[] u = new int[n + 1], v = new int[n + 1], p = new int[n + 1], way = new int[n + 1];   // 1-based; u, v: dual variables
4    for (int i = 1; i <= n; i++) {                                    // insert row i
5        p[0] = i;
6        int j0 = 0;                                                   // column 0 is a dummy that holds the new row
7        int[] minv = new int[n + 1];
8        boolean[] used = new boolean[n + 1];
9        Arrays.fill(minv, Integer.MAX_VALUE);
10        do {
11            used[j0] = true;
12            int i0 = p[j0], delta = Integer.MAX_VALUE, j1 = 0;
13            for (int j = 1; j <= n; j++) {
14                if (used[j]) continue;
15                int cur = cost[i0 - 1][j - 1] - u[i0] - v[j];         // reduced cost (slack) of edge i0 - j
16                if (cur < minv[j]) { minv[j] = cur; way[j] = j0; }
17                if (minv[j] < delta) { delta = minv[j]; j1 = j; }
18            }
19            for (int j = 0; j <= n; j++) {                            // dual update: make the cheapest slack tight
20                if (used[j]) { u[p[j]] += delta; v[j] -= delta; }
21                else minv[j] -= delta;
22            }
23            j0 = j1;
24        } while (p[j0] != 0);                                         // stop at a column nobody owns yet
25        do { int j1 = way[j0]; p[j0] = p[j1]; j0 = j1; } while (j0 != 0);   // augment along the way[] chain
26    }
27    int[] result = new int[n];
28    for (int j = 1; j <= n; j++) result[p[j] - 1] = j - 1;
29    return result;
30}

The Hungarian routine wants a **square** matrix of costs to be **minimised**. Real problems rarely arrive that way, but small transformations fix every common variation. The code below pads a rectangular matrix with zero rows or columns and, when the goal is to maximise profit, replaces each value x by best − x so that maximising profit becomes minimising cost.

Padding with zeros means a dummy worker or a dummy job costs nothing, so the real assignment of min(rows, cols) pairs is chosen and the extra rows or columns absorb the leftovers. For a maximisation problem the padded cells have value 0, which becomes cost best; the optimal profit is then n · best − (minimum cost).

Forbidden pairs are modelled with a very large cost, larger than the sum of every real cost. If the optimal total reaches that magnitude, no valid assignment exists. Keep the numbers in int range, or use long, when you add large penalties.

1// pads a rectangular matrix to a square one; maximize = true turns profits into costs (cost = best - profit)
2static int[][] toSquareCost(int[][] value, boolean maximize) {
3    int rows = value.length, cols = value[0].length, n = Math.max(rows, cols), best = 0;
4    for (int[] r : value) for (int x : r) best = Math.max(best, x);
5    int[][] cost = new int[n][n];
6    for (int i = 0; i < n; i++)
7        for (int j = 0; j < n; j++) {
8            int x = i < rows && j < cols ? value[i][j] : 0;           // a missing worker or job contributes nothing
9            cost[i][j] = maximize ? best - x : x;
10        }
11    return cost;
12}

The assignment problem is a linear program: minimise Σ c[i][j] · x[i][j] with every row and column used exactly once. Its **dual** assigns a number u[i] to each row and v[j] to each column and maximises Σ u + Σ v subject to u[i] + v[j] ≤ c[i][j] for every cell. The inequality means the **reduced cost** c[i][j] − u[i] − v[j] is never negative.

The dual gives a cheap **lower bound**: every assignment picks one cell per row and per column, so its cost is at least Σ u + Σ v. If you can find an assignment whose cells are all **tight** (reduced cost 0), its cost equals the dual value, so it is optimal. This is **complementary slackness**, and the Hungarian algorithm is a systematic search for such a pair.

This is why the animation's final message mentions the dual sum: for the 3 × 3 example both the assignment cost and Σ u + Σ v are 5. The check below verifies feasibility and tightness directly, so an assignment together with duals is a **certificate** you can verify in O(n^2) without trusting the algorithm that produced it.

• Feasible duals always give a lower bound on the optimum.
• An assignment that only uses tight cells under feasible duals is optimal.
• A dual update adds delta to the rows of the tree and subtracts it from the tree's columns, which keeps every cell feasible.

1// u[i] + v[j] <= cost[i][j] in every cell: a feasible dual. Its value is a lower bound on every assignment.
2static boolean isFeasibleDual(int[][] cost, int[] u, int[] v) {
3    for (int i = 0; i < cost.length; i++)
4        for (int j = 0; j < cost.length; j++)
5            if (u[i] + v[j] > cost[i][j]) return false;
6    return true;
7}
8static int dualValue(int[] u, int[] v) {
9    int sum = 0;
10    for (int x : u) sum += x;
11    for (int x : v) sum += x;
12    return sum;
13}
14// row i takes column col[i]; optimal when the duals are feasible and every chosen cell is tight (complementary slackness)
15static boolean isOptimalPair(int[][] cost, int[] col, int[] u, int[] v) {
16    if (!isFeasibleDual(cost, u, v)) return false;
17    for (int i = 0; i < col.length; i++) if (u[i] + v[col[i]] != cost[i][col[i]]) return false;
18    return true;
19}

The **equality subgraph** contains only the tight cells, those with u[i] + v[j] = c[i][j]. The Hungarian algorithm keeps a matching inside this subgraph at all times. If the equality subgraph has a **perfect matching**, that matching is an optimal assignment, by the certificate argument of the previous section.

When the equality subgraph has no perfect matching, the alternating search from a free row gets stuck in a tree. The tree's rows and columns form a bottleneck exactly as in Hall's theorem, and the dual update with delta makes at least one new edge tight at the boundary of this tree without breaking any existing tight edge. After at most n such updates per inserted row the search must succeed, which gives the O(n^3) total.

The code below tests the equality-subgraph condition with the Kuhn routine from the previous lessons: it tries to find a perfect matching using tight cells only. Note that when it returns true for duals produced by a finished Hungarian run, the matching it finds may differ from the one the algorithm returned. Both are optimal; ties in the matrix allow several optimal assignments, as the preset with all costs equal shows (every assignment costs 15).

1static boolean tightPerfectMatching(int[][] cost, int[] u, int[] v) { // a perfect matching using only tight cells?
2    int n = cost.length;
3    int[] rowOf = new int[n];
4    Arrays.fill(rowOf, -1);
5    for (int i = 0; i < n; i++) if (!extend(i, cost, u, v, new boolean[n], rowOf)) return false;
6    return true;
7}
8static boolean extend(int i, int[][] cost, int[] u, int[] v, boolean[] seen, int[] rowOf) {
9    for (int j = 0; j < cost.length; j++) {
10        if (seen[j] || u[i] + v[j] != cost[i][j]) continue;           // only tight cells are edges of the equality graph
11        seen[j] = true;
12        if (rowOf[j] == -1 || extend(rowOf[j], cost, u, v, seen, rowOf)) { rowOf[j] = i; return true; }
13    }
14    return false;
15}