⚙️AlgorithmAdvanced

Hungarian Algorithm

Key Points

•
The Hungarian algorithm solves the square assignment problem (matching n workers to n jobs) in O( $n^{3}$ ) time using a clever potential (label) function on vertices.
•
It transforms costs with vertex potentials so only zero reduced-cost edges matter, then finds a perfect matching inside this equality graph.
•
When a perfect matching is not yet available, it grows an alternating tree and shifts potentials by the minimum slack to create at least one new zero edge.
•
The method is strongly polynomial, works with any real costs (including negatives), and returns an optimal minimum-cost perfect matching.
•
Maximization problems are handled by converting profits to costs, e.g., c'_{ij} = M - $w_{ij}$ with a large enough M, or by negating costs.
•
Rectangular matrices are supported by padding the smaller side with dummy rows/columns and zero (or large) costs to make the matrix square.
•
Its correctness follows from linear programming duality and complementary slackness between primal matchings and dual vertex potentials.
•
Although min-cost max-flow also solves assignment, Hungarian is typically faster and simpler for dense n×n cases.

Prerequisites

→Bipartite graphs and matchings — The assignment problem is a perfect matching on a bipartite graph, and the algorithm augments a matching.
→Greedy algorithms and invariants — Understanding algorithmic invariants helps to follow why potentials remain feasible and progress is ensured.
→Linear programming duality — Potentials are dual variables; correctness relies on dual feasibility and complementary slackness.
→Shortest path and reduced costs — Label updates mimic reduced-cost adjustments, similar to potentials in shortest-path algorithms for min-cost flow.
→Complexity analysis (Big-O) — To reason about O(n^{3}) behavior and compare with alternatives like min-cost max-flow.
→Matrix operations and indices — The algorithm is implemented over an n×n matrix and is sensitive to 0/1-based indexing.

Detailed Explanation

Tap terms for definitions

01Overview

Hook: Imagine you manage n workers and n jobs. You want to assign exactly one worker to each job so the total cost (or time, or distance) is as small as possible. Brute force tries n! assignments—impossible even for modest n. The Hungarian algorithm cracks this in cubic time. Concept: The Hungarian algorithm (also called Kuhn–Munkres) solves the assignment problem on a complete bipartite graph with edge costs. It maintains a potential (label) for each vertex. By subtracting these labels from the edge costs, we get reduced costs. Only edges with zero reduced cost are considered “free” and form the equality graph. The algorithm repeatedly searches for augmenting paths within this equality graph. If none exists, it shifts the vertex potentials just enough to create new zero-cost edges and continues. Example: With a 3×3 cost matrix representing worker–job costs, the algorithm will iteratively adjust row and column labels so that an entire perfect matching can be chosen among zero-cost edges. The final matching is provably optimal, and the method runs in O(n^{3}) time for an n×n matrix.

02Intuition & Analogies

Think of each worker having a personal coupon (a discount) and each job having a coupon too. The real price of assigning worker i to job j is c_{ij}, but with coupons u_{i} (worker’s) and v_{j} (job’s), the adjusted price is c_{ij} - u_{i} - v_{j}. We never allow adjusted prices to go negative—that would break the rules—so u and v must be chosen carefully. The goal is to make enough adjusted prices exactly zero so that we can pick one zero-price edge per worker and per job (a perfect matching) without conflicts. At any moment, the zero adjusted prices define the equality graph. If you can pick n edges there forming a perfect matching, you’re done: everyone gets a job at zero adjusted price, and the total of the coupons exactly accounts for the true minimal cost. If you cannot, you grow an alternating tree from some unmatched worker along zero edges; when you get stuck, you make a tiny, precise coupon adjustment (the minimum slack) so that at least one new zero edge appears. This is like slowly tuning discounts until new free options open up—and you never make options illegally cheaper than free. This dance between selecting zero edges and tuning coupons is guided by a simple invariant: adjusted prices c_{ij} - u_{i} - v_{j} are always nonnegative. Therefore, the sum of all coupons (the dual objective) is always a valid lower bound on the true cost. When you finally get a perfect matching using only zero edges, that lower bound equals the matching’s cost—so it must be optimal. The algorithm is efficient because each adjustment creates progress, and there are only O(n^{2}) edges and O(n) augmentations to consider.

03Formal Definition

We are given a cost matrix C = [

c_{ij}

] for i,j

\in

\{1,

\dots

,n\}. The assignment problem is to minimize

\sum_{i = 1}^{n}

\sum_{j = 1}^{n}

c_{ij}

x_{ij}

subject to

\sum_{j = 1}^{n}

x_{ij} = 1

for all i,

\sum_{i = 1}^{n}

x_{ij} = 1

for all j, and

x_{ij}

\in

\{0,1\}. Its linear relaxation (

x_{ij}

\geq

0) is integral by total unimodularity, so solving the LP suffices. The dual LP is to maximize

\sum_{i = 1}^{n}

u_{i}

\sum_{j = 1}^{n}

v_{j}

subject to

u_{i}

v_{j}

\leq

c_{ij}

for all i,j. Define reduced costs

r_{ij} = c_{ij} - u_{i} - v_{j}

\geq

0. The equality graph contains edges (i,j) with

r_{ij} = 0

. Complementary slackness requires that if

x_{ij} > 0

, then

r_{ij} = 0

, i.e., the final matching must lie entirely on equality edges. The Hungarian algorithm maintains feasibility of the dual (u,v) and a matching on equality edges. If the matching is not perfect, it builds an alternating tree rooted at some unmatched row i using zero edges. When expansion gets stuck, it computes a slack

δ

min_{i \in S, j \in / T}

r_{ij}

and updates labels as

u_{i}

\leftarrow

u_{i}

δ

for i

\in

S and

v_{j}

\leftarrow

v_{j}

δ

for j

\in

T, leaving others unchanged. This preserves dual feasibility, decreases some slacks to zero, and enables new equality edges. Repeating this process yields a perfect matching in O(

n^{3}

04When to Use

Use the Hungarian algorithm when you must assign exactly one item from each of n sources to n targets with minimum total cost, and the graph is dense or represented as a matrix. Typical cases include task assignment to workers, matching students to projects with additive costs, matching predicted objects to detections in tracking (with costs as distances), and aligning features in computer vision or NLP. Prefer Hungarian when: (1) you have a dense n×n cost matrix; (2) edge costs can be arbitrary reals (including negatives); (3) you need deterministic, exact optimality in strongly polynomial O(n^{3}). It is also straightforward to handle rectangular matrices by padding and to convert maximization to minimization. Consider min-cost max-flow when: (1) the bipartite graph is very sparse (far fewer than n^{2} edges); (2) you need to incorporate capacities other than 1, additional side constraints, or multiple units of flow; (3) you need to solve many related problems with small incremental changes and can reuse a flow network. For extremely large n with special structure (e.g., metric costs), specialized or approximate methods (auction algorithm, network simplex, or cost-scaling) might be preferable.

⚠️Common Mistakes

• Using Hungarian directly on a rectangular matrix without padding. Always pad the smaller dimension with dummy rows/columns whose costs are zero (or a neutral constant). The algorithm then finds a perfect matching in the square matrix. • Forgetting to convert maximization to minimization. For profits w_{ij}, set c'{ij} = M - w{ij} where M \ge \max w_{ij}, or use c'{ij} = -w{ij} if your implementation tolerates negatives. The assignment from c' minimizes cost and thus maximizes the original profit. • Overflow due to using an insufficiently large INF or narrow integer type. Keep costs in 64-bit integers (long long) and use INF \approx 9e18/4 for safety when initializing slacks. • Mishandling indices and initialization in the core O(n^{3}) loop (arrays p, way, minv, used). Every augmentation must reset minv and used; p[0] must be the current row; and the backtracking loop must correctly rewire the matching. • Assuming row/column reductions are required as a separate preprocessing step. The standard O(n^{3}) implementation implicitly maintains dual feasibility and makes such reductions on the fly; adding ad-hoc reductions can complicate correctness. • Mixing up 1-based vs 0-based. Many reference codes are 1-based for convenience; porting to 0-based without care often introduces off-by-one errors. • Treating missing edges by putting an extremely large cost but then choosing too small a constant, allowing the algorithm to pick them. If an edge is forbidden, ensure its cost is truly prohibitive (or remove it in a sparse solver like min-cost flow).

Key Formulas

Primal LP (assignment)

min i = 1 \sum n j = 1 \sum n c_{ij} x_{ij} s.t. j = 1 \sum n x_{ij} = 1, i = 1 \sum n x_{ij} = 1, x_{ij} \geq 0

Explanation: This is the linear relaxation of the assignment problem. Because the constraint matrix is totally unimodular, any optimal solution is integral, so solving this LP yields an exact assignment.

Dual LP (potentials)

max i = 1 \sum n u_{i} + j = 1 \sum n v_{j} s.t. u_{i} + v_{j} \leq c_{ij} \forall i, j

Explanation: The dual introduces vertex potentials u and v bounded by edge costs. It provides a lower bound on the primal cost; equality at optimum certifies optimality.

Reduced cost

r_{ij} = c_{ij} - u_{i} - v_{j} \geq 0

Explanation: Reduced costs quantify how far an edge is from being tight (eligible). Zero reduced cost edges form the equality graph from which we build the matching.

Complementary slackness condition

x_{ij} > 0 \Rightarrow r_{ij} = 0

Explanation: Any edge used in the optimal matching must be tight under the optimal potentials. This is the key correctness condition exploited by the algorithm.

Minimum slack for label update

δ = i \in S, j \in / T min r_{ij}

Explanation: When the search gets stuck, we increase u on rows in S and decrease v on columns in T by $δ$ . This preserves feasibility and creates at least one new zero edge.

Potential update rule

u_{i}^{'} = u_{i} + δ (i \in S), v_{j}^{'} = v_{j} - δ (j \in T)

Explanation: This update leaves $r_{ij}$ unchanged on edges inside S $\times$ T, increases tightness for frontier edges, and maintains $r_{ij}$ $\geq$ 0 everywhere.

Time complexity

T (n) = O (n^{3})

Explanation: For an n×n dense matrix, the Hungarian algorithm performs O(n) augmentations, each costing O( $n^{2}$ ) work for scanning/updating slacks, leading to cubic time.

Maximization to minimization

c_{ij}^{'} = M - w_{ij}

Explanation: To maximize profits $w_{ij}$ , choose M $\geq$ $max$ $w_{ij}$ and minimize c' instead. The assignment minimizing c' maximizes the original total profit.

Padding a rectangular matrix

C_{p a d} = [C 0 00]

Explanation: Add dummy rows/columns with zero (or neutral) costs to make the matrix square without changing the true optimal assignment among real items.

Strong duality at optimum

i = 1 \sum n c_{i, π (i)} = i = 1 \sum n (u_{i} + v_{π (i)})

Explanation: When a perfect matching lies on equality edges, its cost equals the sum of potentials, proving optimality by duality.

Complexity Analysis

For an n×n dense cost matrix, the classic Hungarian implementation runs in O(

n^{3}

) time. The algorithm performs exactly n augmentations, each adding one more matched row. During one augmentation, it maintains arrays of slacks (minimum reduced cost to each unmatched column) and a visited mask. Each step decreases some slacks and updates potentials by the global minimum slack. The dominant work per augmentation is scanning all O(n) columns multiple times while growing the alternating tree, plus O(n) updates to potentials and slacks; aggregated, this yields O(

n^{2}

) work per augmentation and O(

n^{3}

) total. Space complexity depends on representation. The dense version stores the n×n matrix (O(

n^{2}

)) plus O(n) auxiliary arrays (potentials u and v, match arrays p and way, minv and used). If you stream the matrix or compute costs on the fly, you can reduce storage to O(n) beyond the matrix itself, but most implementations keep O(

n^{2}

) for simplicity. For rectangular n×m inputs (assuming

n \leq m

after padding or transposition), the running time is O(

n^{2}

m), which becomes O(

n^{3}

) when

m = n

. In contrast, a min-cost max-flow solution on a sparse graph may be advantageous when the number of allowed edges E ≪

n^{2}

; its complexity depends on the chosen shortest-path routine and potentials (commonly O(F E log V) for F units of flow), but for dense matrices Hungarian is typically faster and simpler. Numerically, costs can be any 64-bit integers (including negatives). Potential updates use addition/subtraction and comparisons, so numerical stability is not an issue; just select a sufficiently large INF to avoid overflow when initializing slacks.

Code Examples

Hungarian algorithm (dense O(n^3)) for minimum-cost assignment

1 #include <iostream>
2 #include <vector>
3 #include <limits>
4 #include <iomanip>
5 using namespace std;
6 
7 // Hungarian algorithm for rectangular matrices (n <= m):
8 // Returns pair<min_cost, assignment>, where assignment[i] is the column matched to row i (1-based indices internally, but we output 0-based).
9 // Time: O(n*m*(n+m)) -> O(n^3) for square matrices. Space: O(n+m).
10 pair<long long, vector<int>> hungarian(const vector<vector<long long>>& a) {
11     int n = (int)a.size();
12     int m = (int)a[0].size();
13     // We assume n <= m. If not, transpose before calling.
14     const long long INF = numeric_limits<long long>::max() / 4;
15 
16     vector<long long> u(n + 1, 0), v(m + 1, 0); // potentials for rows and columns
17     vector<int> p(m + 1, 0), way(m + 1, 0);     // matching and parent pointers
18 
19     for (int i = 1; i <= n; ++i) {
20         p[0] = i;                               // column 0 is matched with row i (temporary)
21         int j0 = 0;                             // current column on the augmenting path
22         vector<long long> minv(m + 1, INF);     // best reduced cost to reach column j
23         vector<char> used(m + 1, false);        // whether column j is in the alternating tree
24         do {
25             used[j0] = true;
26             int i0 = p[j0];                     // current row
27             long long delta = INF;
28             int j1 = 0;
29             // Try to relax all columns not yet in the tree
30             for (int j = 1; j <= m; ++j) if (!used[j]) {
31                 // Reduced cost of edge (i0, j)
32                 long long cur = a[i0 - 1][j - 1] - u[i0] - v[j];
33                 if (cur < minv[j]) {            // better way to reach column j
34                     minv[j] = cur;
35                     way[j] = j0;                // remember predecessor column
36                 }
37                 if (minv[j] < delta) {          // keep track of minimum slack
38                     delta = minv[j];
39                     j1 = j;
40                 }
41             }
42             // Update potentials: shift by minimal slack so at least one new zero-cost edge appears
43             for (int j = 0; j <= m; ++j) {
44                 if (used[j]) { u[p[j]] += delta; v[j] -= delta; }
45                 else { minv[j] -= delta; }
46             }
47             j0 = j1;                            // advance to the column with new zero slack
48         } while (p[j0] != 0);                   // until we reach an unmatched column
49 
50         // Augment along the path, rewiring the matching
51         do {
52             int j1 = way[j0];
53             p[j0] = p[j1];
54             j0 = j1;
55         } while (j0 != 0);
56     }
57 
58     vector<int> assignment(n, -1);
59     for (int j = 1; j <= m; ++j) if (p[j] != 0 && p[j] <= n) assignment[p[j] - 1] = j - 1;
60 
61     long long min_cost = 0;
62     for (int i = 0; i < n; ++i) {
63         if (assignment[i] >= 0) min_cost += a[i][assignment[i]];
64     }
65     return {min_cost, assignment};
66 }
67 
68 int main() {
69     ios::sync_with_stdio(false);
70     cin.tie(nullptr);
71 
72     int n; // square case for demo; rectangular also works if n <= m
73     cin >> n;
74     vector<vector<long long>> cost(n, vector<long long>(n));
75     for (int i = 0; i < n; ++i) for (int j = 0; j < n; ++j) cin >> cost[i][j];
76 
77     // If rows > columns, transpose the matrix before calling (not shown here).
78     auto [min_cost, match] = hungarian(cost);
79 
80     cout << "Minimum cost: " << min_cost << "\nAssignment (row -> col):\n";
81     for (int i = 0; i < n; ++i) cout << i << " -> " << match[i] << "\n";
82     return 0;
83 }
84

This is the classic dense Hungarian implementation using 1-based arrays for potentials u, v and helper arrays p, way. For each row i, we build an alternating tree over columns and maintain minv as the current slack to each column. We repeatedly pick the column with minimum slack, shift potentials by that slack to create new zero edges, and continue until we reach an unmatched column, at which point we augment. After n augmentations, we recover the assignment and its minimum total cost.

Time: O(n^3) for n×n; O(n^2 m) for n×m with n <= mSpace: O(n + m) auxiliary plus the input matrix O(nm)

Maximization and rectangular matrices via padding and cost transform

1 #include <iostream>
2 #include <vector>
3 #include <limits>
4 #include <algorithm>
5 using namespace std;
6 
7 pair<long long, vector<int>> hungarian(const vector<vector<long long>>& a); // forward decl (reuse from previous example)
8 
9 // Utility: pad a rectangular profit matrix to square and transform to costs for maximization
10 pair<long long, vector<int>> max_assignment(vector<vector<long long>> profit) {
11     int n = (int)profit.size();
12     int m = (int)profit[0].size();
13 
14     // Pad to square by adding dummy rows/cols with 0 profit
15     int N = max(n, m);
16     vector<vector<long long>> P(N, vector<long long>(N, 0));
17     for (int i = 0; i < n; ++i)
18         for (int j = 0; j < m; ++j)
19             P[i][j] = profit[i][j];
20 
21     // Convert maximization to minimization: c' = M - w
22     long long M = 0;
23     for (int i = 0; i < N; ++i)
24         for (int j = 0; j < N; ++j)
25             M = max(M, P[i][j]);
26 
27     vector<vector<long long>> cost(N, vector<long long>(N));
28     for (int i = 0; i < N; ++i)
29         for (int j = 0; j < N; ++j)
30             cost[i][j] = M - P[i][j];
31 
32     auto res = hungarian(cost);
33 
34     // Recover maximum profit from the assignment (ignore dummy rows/cols)
35     long long max_profit = 0;
36     vector<int> assign(n, -1);
37     for (int i = 0; i < n; ++i) {
38         int j = res.second[i];
39         if (j < m) {
40             assign[i] = j;
41             max_profit += profit[i][j];
42         }
43     }
44     return {max_profit, assign};
45 }
46 
47 // Minimal in-file Hungarian implementation for completeness
48 pair<long long, vector<int>> hungarian(const vector<vector<long long>>& a) {
49     int n = (int)a.size();
50     int m = (int)a[0].size();
51     const long long INF = numeric_limits<long long>::max() / 4;
52     vector<long long> u(n + 1), v(m + 1);
53     vector<int> p(m + 1), way(m + 1);
54     for (int i = 1; i <= n; ++i) {
55         p[0] = i; int j0 = 0; vector<long long> minv(m + 1, INF); vector<char> used(m + 1, false);
56         do {
57             used[j0] = true; int i0 = p[j0]; long long delta = INF; int j1 = 0;
58             for (int j = 1; j <= m; ++j) if (!used[j]) {
59                 long long cur = a[i0 - 1][j - 1] - u[i0] - v[j];
60                 if (cur < minv[j]) { minv[j] = cur; way[j] = j0; }
61                 if (minv[j] < delta) { delta = minv[j]; j1 = j; }
62             }
63             for (int j = 0; j <= m; ++j) if (used[j]) { u[p[j]] += delta; v[j] -= delta; } else { minv[j] -= delta; }
64             j0 = j1;
65         } while (p[j0] != 0);
66         do { int j1 = way[j0]; p[j0] = p[j1]; j0 = j1; } while (j0);
67     }
68     vector<int> assignment(n, -1);
69     for (int j = 1; j <= m; ++j) if (p[j] && p[j] <= n) assignment[p[j] - 1] = j - 1;
70     long long min_cost = 0; for (int i = 0; i < n; ++i) if (assignment[i] >= 0) min_cost += a[i][assignment[i]];
71     return {min_cost, assignment};
72 }
73 
74 int main() {
75     ios::sync_with_stdio(false);
76     cin.tie(nullptr);
77 
78     int n, m; cin >> n >> m; // rectangular profits
79     vector<vector<long long>> profit(n, vector<long long>(m));
80     for (int i = 0; i < n; ++i) for (int j = 0; j < m; ++j) cin >> profit[i][j];
81 
82     auto [best, assign] = max_assignment(profit);
83     cout << "Maximum profit: " << best << "\nAssignment (row -> col):\n";
84     for (int i = 0; i < n; ++i) cout << i << " -> " << assign[i] << "\n";
85     return 0;
86 }
87

We show how to solve a maximization and rectangular assignment. First, pad the profit matrix to square by adding dummy rows/columns with zero profit. Then transform to a minimization by setting c' = M - w, where M is the maximum profit in the padded matrix. Running Hungarian on c' yields an assignment that maximizes total profit in the original problem.

Time: O(N^3) after padding to N×N where N = max(n, m)Space: O(N^2) for the padded matrix plus O(N) auxiliary

Assignment via Min-Cost Max-Flow (alternative approach)

1 #include <bits/stdc++.h>
2 using namespace std;
3 
4 struct MinCostMaxFlow {
5     struct Edge { int v, cap, rev; long long cost; };
6     int N; vector<vector<Edge>> G; vector<long long> dist, pot; vector<int> pv_v, pv_e;
7     MinCostMaxFlow(int n): N(n), G(n), dist(n), pot(n), pv_v(n), pv_e(n) {}
8     void addEdge(int u, int v, int cap, long long cost) {
9         Edge a{v, cap, (int)G[v].size(), cost};
10         Edge b{u, 0,   (int)G[u].size(), -cost};
11         G[u].push_back(a); G[v].push_back(b);
12     }
13     pair<int,long long> minCostMaxFlow(int s, int t, int maxf = INT_MAX) {
14         int flow = 0; long long cost = 0; fill(pot.begin(), pot.end(), 0);
15         // Optionally, run Bellman-Ford to initialize potentials for negative edges (skipped if costs are nonnegative)
16         while (flow < maxf) {
17             fill(dist.begin(), dist.end(), (long long)4e18); dist[s] = 0;
18             priority_queue<pair<long long,int>, vector<pair<long long,int>>, greater<pair<long long,int>>> pq; pq.push({0, s});
19             fill(pv_v.begin(), pv_v.end(), -1); fill(pv_e.begin(), pv_e.end(), -1);
20             while (!pq.empty()) {
21                 auto [d, u] = pq.top(); pq.pop(); if (d != dist[u]) continue;
22                 for (int i = 0; i < (int)G[u].size(); ++i) {
23                     auto &e = G[u][i]; if (e.cap <= 0) continue;
24                     long long nd = d + e.cost + pot[u] - pot[e.v];
25                     if (nd < dist[e.v]) { dist[e.v] = nd; pv_v[e.v] = u; pv_e[e.v] = i; pq.push({nd, e.v}); }
26                 }
27             }
28             if (pv_v[t] == -1) break; // no augmenting path
29             for (int v = 0; v < N; ++v) if (dist[v] < (long long)4e18) pot[v] += dist[v];
30             int addf = maxf - flow; for (int v = t; v != s; v = pv_v[v]) addf = min(addf, G[pv_v[v]][pv_e[v]].cap);
31             flow += addf; cost += 1LL * addf * pot[t];
32             for (int v = t; v != s; v = pv_v[v]) { auto &e = G[pv_v[v]][pv_e[v]]; e.cap -= addf; G[v][e.rev].cap += addf; }
33         }
34         return {flow, cost};
35     }
36 };
37 
38 int main(){
39     ios::sync_with_stdio(false); cin.tie(nullptr);
40     int n; cin >> n; // square assignment
41     vector<vector<long long>> c(n, vector<long long>(n));
42     for(int i=0;i<n;++i) for(int j=0;j<n;++j) cin>>c[i][j];
43 
44     int S = 2*n, T = 2*n+1; MinCostMaxFlow mcmf(2*n+2);
45     for(int i=0;i<n;++i) mcmf.addEdge(S, i, 1, 0);
46     for(int j=0;j<n;++j) mcmf.addEdge(n+j, T, 1, 0);
47     for(int i=0;i<n;++i) for(int j=0;j<n;++j) mcmf.addEdge(i, n+j, 1, c[i][j]);
48 
49     auto [flow, cost] = mcmf.minCostMaxFlow(S, T, n);
50     cout << "Minimum cost: " << cost << "\nFlow sent: " << flow << "\n";
51 
52     // Recover assignment from residual graph (edges with cap==0 on i->n+j used once)
53     vector<int> assign(n, -1);
54     for(int i=0;i<n;++i){
55         for(auto &e : mcmf.G[i]) if(e.v>=n && e.v<2*n && e.cap==0) assign[i] = e.v - n;
56     }
57     for(int i=0;i<n;++i) cout << i << " -> " << assign[i] << "\n";
58     return 0;
59 }
60

This builds a standard flow network for assignment: source to each left node with capacity 1, left-to-right edges with capacity 1 and cost c_{ij}, and each right node to sink with capacity 1. We run successive shortest augmenting paths with reduced costs (potentials) to achieve min-cost max-flow of value n. This is an alternative to Hungarian, often better for sparse graphs or when adding constraints.

Time: Approximately O(n * E log V) with Dijkstra and potentials; for dense assignment E = n^2 so roughly O(n^3 log n)Space: O(V + E) ~ O(n^2) for dense graphs

1	#include <iostream>
2	#include <vector>
3	#include <limits>
4	#include <iomanip>
5	using namespace std;
6
7	// Hungarian algorithm for rectangular matrices (n <= m):
8	// Returns pair<min_cost, assignment>, where assignment[i] is the column matched to row i (1-based indices internally, but we output 0-based).
9	// Time: O(nm(n+m)) -> O(n^3) for square matrices. Space: O(n+m).
10	pair<long long, vector<int>> hungarian(const vector<vector<long long>>& a) {
11	int n = (int)a.size();
12	int m = (int)a[0].size();
13	// We assume n <= m. If not, transpose before calling.
14	const long long INF = numeric_limits<long long>::max() / 4;
15
16	vector<long long> u(n + 1, 0), v(m + 1, 0); // potentials for rows and columns
17	vector<int> p(m + 1, 0), way(m + 1, 0); // matching and parent pointers
18
19	for (int i = 1; i <= n; ++i) {
20	p[0] = i; // column 0 is matched with row i (temporary)
21	int j0 = 0; // current column on the augmenting path
22	vector<long long> minv(m + 1, INF); // best reduced cost to reach column j
23	vector<char> used(m + 1, false); // whether column j is in the alternating tree
24	do {
25	used[j0] = true;
26	int i0 = p[j0]; // current row
27	long long delta = INF;
28	int j1 = 0;
29	// Try to relax all columns not yet in the tree
30	for (int j = 1; j <= m; ++j) if (!used[j]) {
31	// Reduced cost of edge (i0, j)
32	long long cur = a[i0 - 1][j - 1] - u[i0] - v[j];
33	if (cur < minv[j]) { // better way to reach column j
34	minv[j] = cur;
35	way[j] = j0; // remember predecessor column
36	}
37	if (minv[j] < delta) { // keep track of minimum slack
38	delta = minv[j];
39	j1 = j;
40	}
41	}
42	// Update potentials: shift by minimal slack so at least one new zero-cost edge appears
43	for (int j = 0; j <= m; ++j) {
44	if (used[j]) { u[p[j]] += delta; v[j] -= delta; }
45	else { minv[j] -= delta; }
46	}
47	j0 = j1; // advance to the column with new zero slack
48	} while (p[j0] != 0); // until we reach an unmatched column
49
50	// Augment along the path, rewiring the matching
51	do {
52	int j1 = way[j0];
53	p[j0] = p[j1];
54	j0 = j1;
55	} while (j0 != 0);
56	}
57
58	vector<int> assignment(n, -1);
59	for (int j = 1; j <= m; ++j) if (p[j] != 0 && p[j] <= n) assignment[p[j] - 1] = j - 1;
60
61	long long min_cost = 0;
62	for (int i = 0; i < n; ++i) {
63	if (assignment[i] >= 0) min_cost += a[i][assignment[i]];
64	}
65	return {min_cost, assignment};
66	}
67
68	int main() {
69	ios::sync_with_stdio(false);
70	cin.tie(nullptr);
71
72	int n; // square case for demo; rectangular also works if n <= m
73	cin >> n;
74	vector<vector<long long>> cost(n, vector<long long>(n));
75	for (int i = 0; i < n; ++i) for (int j = 0; j < n; ++j) cin >> cost[i][j];
76
77	// If rows > columns, transpose the matrix before calling (not shown here).
78	auto [min_cost, match] = hungarian(cost);
79
80	cout << "Minimum cost: " << min_cost << "\nAssignment (row -> col):\n";
81	for (int i = 0; i < n; ++i) cout << i << " -> " << match[i] << "\n";
82	return 0;
83	}
84

1	#include <bits/stdc++.h>
2	using namespace std;
3
4	struct MinCostMaxFlow {
5	struct Edge { int v, cap, rev; long long cost; };
6	int N; vector<vector<Edge>> G; vector<long long> dist, pot; vector<int> pv_v, pv_e;
7	MinCostMaxFlow(int n): N(n), G(n), dist(n), pot(n), pv_v(n), pv_e(n) {}
8	void addEdge(int u, int v, int cap, long long cost) {
9	Edge a{v, cap, (int)G[v].size(), cost};
10	Edge b{u, 0, (int)G[u].size(), -cost};
11	G[u].push_back(a); G[v].push_back(b);
12	}
13	pair<int,long long> minCostMaxFlow(int s, int t, int maxf = INT_MAX) {
14	int flow = 0; long long cost = 0; fill(pot.begin(), pot.end(), 0);
15	// Optionally, run Bellman-Ford to initialize potentials for negative edges (skipped if costs are nonnegative)
16	while (flow < maxf) {
17	fill(dist.begin(), dist.end(), (long long)4e18); dist[s] = 0;
18	priority_queue<pair<long long,int>, vector<pair<long long,int>>, greater<pair<long long,int>>> pq; pq.push({0, s});
19	fill(pv_v.begin(), pv_v.end(), -1); fill(pv_e.begin(), pv_e.end(), -1);
20	while (!pq.empty()) {
21	auto [d, u] = pq.top(); pq.pop(); if (d != dist[u]) continue;
22	for (int i = 0; i < (int)G[u].size(); ++i) {
23	auto &e = G[u][i]; if (e.cap <= 0) continue;
24	long long nd = d + e.cost + pot[u] - pot[e.v];
25	if (nd < dist[e.v]) { dist[e.v] = nd; pv_v[e.v] = u; pv_e[e.v] = i; pq.push({nd, e.v}); }
26	}
27	}
28	if (pv_v[t] == -1) break; // no augmenting path
29	for (int v = 0; v < N; ++v) if (dist[v] < (long long)4e18) pot[v] += dist[v];
30	int addf = maxf - flow; for (int v = t; v != s; v = pv_v[v]) addf = min(addf, G[pv_v[v]][pv_e[v]].cap);
31	flow += addf; cost += 1LL * addf * pot[t];
32	for (int v = t; v != s; v = pv_v[v]) { auto &e = G[pv_v[v]][pv_e[v]]; e.cap -= addf; G[v][e.rev].cap += addf; }
33	}
34	return {flow, cost};
35	}
36	};
37
38	int main(){
39	ios::sync_with_stdio(false); cin.tie(nullptr);
40	int n; cin >> n; // square assignment
41	vector<vector<long long>> c(n, vector<long long>(n));
42	for(int i=0;i<n;++i) for(int j=0;j<n;++j) cin>>c[i][j];
43
44	int S = 2n, T = 2n+1; MinCostMaxFlow mcmf(2*n+2);
45	for(int i=0;i<n;++i) mcmf.addEdge(S, i, 1, 0);
46	for(int j=0;j<n;++j) mcmf.addEdge(n+j, T, 1, 0);
47	for(int i=0;i<n;++i) for(int j=0;j<n;++j) mcmf.addEdge(i, n+j, 1, c[i][j]);
48
49	auto [flow, cost] = mcmf.minCostMaxFlow(S, T, n);
50	cout << "Minimum cost: " << cost << "\nFlow sent: " << flow << "\n";
51
52	// Recover assignment from residual graph (edges with cap==0 on i->n+j used once)
53	vector<int> assign(n, -1);
54	for(int i=0;i<n;++i){
55	for(auto &e : mcmf.G[i]) if(e.v>=n && e.v<2*n && e.cap==0) assign[i] = e.v - n;
56	}
57	for(int i=0;i<n;++i) cout << i << " -> " << assign[i] << "\n";
58	return 0;
59	}
60