m1une's library

This documentation is automatically generated by online-judge-tools/verification-helper

View on GitHub

:heavy_check_mark: Hungarian Algorithm
(optimization/hungarian.hpp)

Overview

hungarian_min(cost) solves the assignment problem for a rectangular cost matrix. It chooses distinct columns for rows, or distinct rows for columns, minimizing the total selected cost.

If the matrix has H rows and W columns, the algorithm creates min(H, W) pairs and completely matches the smaller side. When H <= W, every row is assigned to one column. When H > W, every column is assigned to one row and the remaining rows are left unassigned.

hungarian_max(cost) solves the corresponding maximization problem by negating costs internally. hungarian(cost) is an alias of hungarian_min(cost).

Interface

The input is a rectangular std::vector<std::vector<T>>. All rows must have the same length. T should be a signed numeric type with enough range for subtraction and summation, such as long long or long double.

HungarianResult<T> contains these members:

Member / Method Type / Signature Meaning
cost T Total cost of the returned assignment.
row_to_col std::vector<int> row_to_col[i] is the column assigned to row i, or -1 if row i is unassigned.
col_to_row std::vector<int> col_to_row[j] is the row assigned to column j, or -1 if column j is unassigned.
matching_size int matching_size() const Returns the number of assigned pairs.
matching std::vector<std::pair<int, int>> matching() const Returns (row, column) pairs.

Functions

Function Signature Description Complexity
hungarian_min template <class T> HungarianResult<T> hungarian_min(const std::vector<std::vector<T>>& cost) Returns a minimum-cost assignment. $O(\min(H, W)^2 \max(H, W))$
hungarian_max template <class T> HungarianResult<T> hungarian_max(const std::vector<std::vector<T>>& cost) Returns a maximum-cost assignment. $O(\min(H, W)^2 \max(H, W))$
hungarian template <class T> HungarianResult<T> hungarian(const std::vector<std::vector<T>>& cost) Alias of hungarian_min. $O(\min(H, W)^2 \max(H, W))$

Example

#include "optimization/hungarian.hpp"
#include <iostream>
#include <vector>

int main() {
    std::vector<std::vector<long long>> cost = {
        {4, 1, 3},
        {2, 0, 5},
        {3, 2, 2},
    };

    auto result = m1une::opt::hungarian_min(cost);
    std::cout << result.cost << "\n";  // 5
    for (auto [row, col] : result.matching()) {
        std::cout << row << " " << col << "\n";
    }
}

Required by

Verified with

Code

#ifndef M1UNE_OPTIMIZATION_HUNGARIAN_HPP
#define M1UNE_OPTIMIZATION_HUNGARIAN_HPP 1

#include <algorithm>
#include <cassert>
#include <limits>
#include <utility>
#include <vector>

namespace m1une {
namespace opt {

template <class T>
struct HungarianResult {
    T cost;
    std::vector<int> row_to_col;
    std::vector<int> col_to_row;

    int matching_size() const {
        int result = 0;
        for (int col : row_to_col) {
            if (col != -1) result++;
        }
        return result;
    }

    std::vector<std::pair<int, int>> matching() const {
        std::vector<std::pair<int, int>> result;
        for (int row = 0; row < int(row_to_col.size()); row++) {
            if (row_to_col[row] != -1) result.push_back({row, row_to_col[row]});
        }
        return result;
    }
};

namespace detail {

template <class T>
T assignment_cost(const std::vector<std::vector<T>>& cost, const std::vector<int>& row_to_col) {
    T result = T();
    for (int row = 0; row < int(row_to_col.size()); row++) {
        if (row_to_col[row] != -1) result += cost[row][row_to_col[row]];
    }
    return result;
}

}  // namespace detail

template <class T>
HungarianResult<T> hungarian_min(const std::vector<std::vector<T>>& cost) {
    int row_count = int(cost.size());
    int col_count = row_count == 0 ? 0 : int(cost[0].size());
    for (const auto& row : cost) assert(int(row.size()) == col_count);

    HungarianResult<T> result;
    result.cost = T();
    result.row_to_col.assign(row_count, -1);
    result.col_to_row.assign(col_count, -1);
    if (row_count == 0 || col_count == 0) return result;

    bool transposed = row_count > col_count;
    int n = transposed ? col_count : row_count;
    int m = transposed ? row_count : col_count;
    T inf = std::numeric_limits<T>::max() / T(4);

    std::vector<T> u(n + 1, T()), v(m + 1, T()), minv(m + 1);
    std::vector<int> p(m + 1, 0), way(m + 1, 0);

    auto value = [&](int i, int j) -> T {
        return transposed ? cost[j][i] : cost[i][j];
    };

    for (int i = 1; i <= n; i++) {
        p[0] = i;
        int j0 = 0;
        std::fill(minv.begin(), minv.end(), inf);
        std::vector<char> used(m + 1, false);

        do {
            used[j0] = true;
            int i0 = p[j0];
            int j1 = 0;
            T delta = inf;

            for (int j = 1; j <= m; j++) {
                if (used[j]) continue;
                T cur = value(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 (int j = 0; j <= m; j++) {
                if (used[j]) {
                    u[p[j]] += delta;
                    v[j] -= delta;
                } else {
                    minv[j] -= delta;
                }
            }
            j0 = j1;
        } while (p[j0] != 0);

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

    for (int j = 1; j <= m; j++) {
        if (p[j] == 0) continue;
        int i = p[j] - 1;
        int matched = j - 1;
        if (transposed) {
            int row = matched;
            int col = i;
            result.row_to_col[row] = col;
            result.col_to_row[col] = row;
        } else {
            int row = i;
            int col = matched;
            result.row_to_col[row] = col;
            result.col_to_row[col] = row;
        }
    }
    result.cost = detail::assignment_cost(cost, result.row_to_col);
    return result;
}

template <class T>
HungarianResult<T> hungarian_max(const std::vector<std::vector<T>>& cost) {
    std::vector<std::vector<T>> negated = cost;
    for (auto& row : negated) {
        for (auto& x : row) x = -x;
    }
    auto result = hungarian_min(negated);
    result.cost = detail::assignment_cost(cost, result.row_to_col);
    return result;
}

template <class T>
HungarianResult<T> hungarian(const std::vector<std::vector<T>>& cost) {
    return hungarian_min(cost);
}

}  // namespace opt
}  // namespace m1une

#endif  // M1UNE_OPTIMIZATION_HUNGARIAN_HPP
#line 1 "optimization/hungarian.hpp"



#include <algorithm>
#include <cassert>
#include <limits>
#include <utility>
#include <vector>

namespace m1une {
namespace opt {

template <class T>
struct HungarianResult {
    T cost;
    std::vector<int> row_to_col;
    std::vector<int> col_to_row;

    int matching_size() const {
        int result = 0;
        for (int col : row_to_col) {
            if (col != -1) result++;
        }
        return result;
    }

    std::vector<std::pair<int, int>> matching() const {
        std::vector<std::pair<int, int>> result;
        for (int row = 0; row < int(row_to_col.size()); row++) {
            if (row_to_col[row] != -1) result.push_back({row, row_to_col[row]});
        }
        return result;
    }
};

namespace detail {

template <class T>
T assignment_cost(const std::vector<std::vector<T>>& cost, const std::vector<int>& row_to_col) {
    T result = T();
    for (int row = 0; row < int(row_to_col.size()); row++) {
        if (row_to_col[row] != -1) result += cost[row][row_to_col[row]];
    }
    return result;
}

}  // namespace detail

template <class T>
HungarianResult<T> hungarian_min(const std::vector<std::vector<T>>& cost) {
    int row_count = int(cost.size());
    int col_count = row_count == 0 ? 0 : int(cost[0].size());
    for (const auto& row : cost) assert(int(row.size()) == col_count);

    HungarianResult<T> result;
    result.cost = T();
    result.row_to_col.assign(row_count, -1);
    result.col_to_row.assign(col_count, -1);
    if (row_count == 0 || col_count == 0) return result;

    bool transposed = row_count > col_count;
    int n = transposed ? col_count : row_count;
    int m = transposed ? row_count : col_count;
    T inf = std::numeric_limits<T>::max() / T(4);

    std::vector<T> u(n + 1, T()), v(m + 1, T()), minv(m + 1);
    std::vector<int> p(m + 1, 0), way(m + 1, 0);

    auto value = [&](int i, int j) -> T {
        return transposed ? cost[j][i] : cost[i][j];
    };

    for (int i = 1; i <= n; i++) {
        p[0] = i;
        int j0 = 0;
        std::fill(minv.begin(), minv.end(), inf);
        std::vector<char> used(m + 1, false);

        do {
            used[j0] = true;
            int i0 = p[j0];
            int j1 = 0;
            T delta = inf;

            for (int j = 1; j <= m; j++) {
                if (used[j]) continue;
                T cur = value(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 (int j = 0; j <= m; j++) {
                if (used[j]) {
                    u[p[j]] += delta;
                    v[j] -= delta;
                } else {
                    minv[j] -= delta;
                }
            }
            j0 = j1;
        } while (p[j0] != 0);

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

    for (int j = 1; j <= m; j++) {
        if (p[j] == 0) continue;
        int i = p[j] - 1;
        int matched = j - 1;
        if (transposed) {
            int row = matched;
            int col = i;
            result.row_to_col[row] = col;
            result.col_to_row[col] = row;
        } else {
            int row = i;
            int col = matched;
            result.row_to_col[row] = col;
            result.col_to_row[col] = row;
        }
    }
    result.cost = detail::assignment_cost(cost, result.row_to_col);
    return result;
}

template <class T>
HungarianResult<T> hungarian_max(const std::vector<std::vector<T>>& cost) {
    std::vector<std::vector<T>> negated = cost;
    for (auto& row : negated) {
        for (auto& x : row) x = -x;
    }
    auto result = hungarian_min(negated);
    result.cost = detail::assignment_cost(cost, result.row_to_col);
    return result;
}

template <class T>
HungarianResult<T> hungarian(const std::vector<std::vector<T>>& cost) {
    return hungarian_min(cost);
}

}  // namespace opt
}  // namespace m1une
Back to top page