Hungarian 法
(Calculate/Hungarian.hpp)

• View this file on GitHub
• Last update: 2026-02-15 14:56:42+09:00
• Include: #include "Calculate/Hungarian.hpp"

Outline

以下で表される割当問題を Hungarian 法で解く.

$N \leq M$ とする. $N$ 行 $M$ 列の行列 $A = (A_{i,j})_{\substack{1 \leq i \leq N \ 1 \leq j \leq M}}$ が与えられる.

このとき,

\[\sum_{i=1}^N A_{i, P_i}\]

を最小にする, 以下を満たす長さ $N$ の列を求める.

• $1 \leq P_i \leq M$.
• $P_1, \dots, P_M$ は全て異なる.

このクラスでは, この問題を Hungarian 法で解く.

History

日付	内容
2026/02/15	スライド最大値 実装

Verified with

• verify/yosupo_library_checker/graph/Assignment_Problem.test.cpp

Code

#include <iostream>
#include <vector>
#include <algorithm>
#include <limits>

using namespace std;

template<typename T>
class Hungarian {
    int n, m;
    vector<vector<T>> matrix;
    const T INF = numeric_limits<T>::max();
    vector<T> u, v, minv;
    vector<int> p, way;
    vector<bool> used;
    bool maximize_mode;

    void step(int &j0) {
        used[j0] = true;
        int i0 = p[j0];
        T delta = INF;
        int j1 = -1;

        for (int j = 0; j < m; ++j) {
            if (!used[j]) {
                T cur = matrix[i0][j] - u[i0] - v[j];
                if (cur < minv[j]) {
                    minv[j] = cur;
                    way[j] = j0;
                }
                if (minv[j] < delta) {
                    delta = minv[j];
                    j1 = j;
                }
            }
        }

        for (int j = 0; j <= m; ++j) {
            if (used[j]) {
                u[p[j]] += delta;
                v[j] -= delta;
            } else {
                if (minv[j] < INF) minv[j] -= delta;
            }
        }
        j0 = j1;
    }

    void update_matching(int j0) {
        do {
            int j1 = way[j0];
            p[j0] = p[j1];
            j0 = j1;
        } while (j0 != m);
    }

    void construct_result() {
        total_cost = 0;
        matching.assign(n, -1);
        for (int j = 0; j < m; ++j) {
            int i = p[j];
            if (i != n) {
                matching[i] = j;
                total_cost += matrix[i][j];
            }
        }
        if (maximize_mode) total_cost = -total_cost;
    }

    void augment(int s) {
        p[m] = s;
        int j0 = m;

        fill(minv.begin(), minv.end(), INF);
        fill(used.begin(), used.end(), false);

        do {
            step(j0);
        } while (p[j0] != n);

        update_matching(j0);
    }

    public:
    vector<int> matching;
    T total_cost;

    Hungarian(vector<vector<T>> cost_matrix, bool maximize = false) : matrix(cost_matrix), maximize_mode(maximize) {
        n = matrix.size();
        m = n == 0 ? 0 : matrix[0].size();

        if (maximize_mode) {
            for (int i = 0; i < n; ++i) {
                for (int j = 0; j < m; ++j) {
                    matrix[i][j] = -matrix[i][j];
                }
            }
        }

        // initialize
        // 0-based indexing for internal logic, n is dummy
        u.assign(n + 1, 0);
        v.assign(m + 1, 0);
        p.assign(m + 1, n);
        way.assign(m + 1, -1);
        minv.assign(m + 1, INF);
        used.assign(m + 1, false);

        solve();
    }

    void solve() {
        for (int i = 0; i < n; ++i) augment(i);

        construct_result();
    }
};

#line 1 "Calculate/Hungarian.hpp"
#include <iostream>
#include <vector>
#include <algorithm>
#include <limits>

using namespace std;

template<typename T>
class Hungarian {
    int n, m;
    vector<vector<T>> matrix;
    const T INF = numeric_limits<T>::max();
    vector<T> u, v, minv;
    vector<int> p, way;
    vector<bool> used;
    bool maximize_mode;

    void step(int &j0) {
        used[j0] = true;
        int i0 = p[j0];
        T delta = INF;
        int j1 = -1;

        for (int j = 0; j < m; ++j) {
            if (!used[j]) {
                T cur = matrix[i0][j] - u[i0] - v[j];
                if (cur < minv[j]) {
                    minv[j] = cur;
                    way[j] = j0;
                }
                if (minv[j] < delta) {
                    delta = minv[j];
                    j1 = j;
                }
            }
        }

        for (int j = 0; j <= m; ++j) {
            if (used[j]) {
                u[p[j]] += delta;
                v[j] -= delta;
            } else {
                if (minv[j] < INF) minv[j] -= delta;
            }
        }
        j0 = j1;
    }

    void update_matching(int j0) {
        do {
            int j1 = way[j0];
            p[j0] = p[j1];
            j0 = j1;
        } while (j0 != m);
    }

    void construct_result() {
        total_cost = 0;
        matching.assign(n, -1);
        for (int j = 0; j < m; ++j) {
            int i = p[j];
            if (i != n) {
                matching[i] = j;
                total_cost += matrix[i][j];
            }
        }
        if (maximize_mode) total_cost = -total_cost;
    }

    void augment(int s) {
        p[m] = s;
        int j0 = m;

        fill(minv.begin(), minv.end(), INF);
        fill(used.begin(), used.end(), false);

        do {
            step(j0);
        } while (p[j0] != n);

        update_matching(j0);
    }

    public:
    vector<int> matching;
    T total_cost;

    Hungarian(vector<vector<T>> cost_matrix, bool maximize = false) : matrix(cost_matrix), maximize_mode(maximize) {
        n = matrix.size();
        m = n == 0 ? 0 : matrix[0].size();

        if (maximize_mode) {
            for (int i = 0; i < n; ++i) {
                for (int j = 0; j < m; ++j) {
                    matrix[i][j] = -matrix[i][j];
                }
            }
        }

        // initialize
        // 0-based indexing for internal logic, n is dummy
        u.assign(n + 1, 0);
        v.assign(m + 1, 0);
        p.assign(m + 1, n);
        way.assign(m + 1, -1);
        minv.assign(m + 1, INF);
        used.assign(m + 1, false);

        solve();
    }

    void solve() {
        for (int i = 0; i < n; ++i) augment(i);

        construct_result();
    }
};

Back to top page