[图论] KM 算法

· · 题解

KM 算法

二分图最大权完美匹配之所以不容易求解,是因为它不能像最大匹配一样贪心寻找增广路,或者说难以构造一种直接贪心的策略。这时候我们考虑先给原问题确定一个子图,并想方设法保证子图上可以用最大匹配贪心求解,如果这个子图包含了最大权完美匹配,那么答案能直接用最大匹配求解了。

不难发现该算法的难点是保证子图可以贪心这一点。我们用以下方法解决这个问题:

给每个点赋一个顶标值,用 l_p 表示。设 w_{u, v} 表示 u \to v 的边权,我们必须保证对于任意 u, vl_u + l_v \geq w_{u, v}。此时只考虑 l_u + l_v = w_{u, v} 的子图,我们声称如果该子图有完美匹配,那么这一定也是一个最大权完美匹配。

这是因为任意一个完美匹配的权值是 \sum w_{i, match_i},根据 l_u + l_v \geq w_{u, v},这权值肯定是不大于 \sum l_i + l_{match_i},也就是所有点 l 值的总和的。而子图又保证了 l_u + l_v = w_{u, v},这意味着 \sum w_{i, match_i} 已经取到上界。

于是问题转化为如何给顶标数组 l 赋值。我们假设当前已经有一组(不完美)匹配和一组可能不完善的顶标 l,我们继续尝试找到左部的一个未匹配点 now,可能产生两种情况:

我们发现对于后一种情况,存在以下方法修改顶标使得增广路可能增加:

这样 T' 就可能存在一些点与 S 连通了。(我们无需考虑 S'T 会不会断开一些边这种问题,因为它们都不在一棵交错树里,显然是非匹配边)我们重新增广一下即可。

最后的小问题就是如何确定 dt 的值。显然 dt 应该是最小的、使得 ST' 有新边的值,这可以在增广的时候顺便记录一下。

slack_i 表示对于右部点 i,需要加上至少多大的值才能满足其会被扩容进入交错树。

我们用 BFS 写法实现,这样可以保证对于每个点,最多扩容 n 次,且刚好执行完一次匈牙利算法求最大匹配(也就是求最大匹配和扩容交错树同时进行),这样复杂度是 O(n^3)

如果你用 DFS 实现,那么首先你无法做到求解最大匹配和扩容并行,因此复杂度不是并列的,而是 O(n^4),其次你很可能写个错的算法出来,具体原因可以看看其他楼层的题解是怎么说的,这里就不复述一遍了。

代码

实现时注意一点问题:

对于没有边的 (u, v),赋值 -\infin 时这个 \infin 要很大很大,远大于边权几个量级才行。不然可能出现先忍痛走 -\infin 边,然后疯狂连接其他正权边导致求出的答案反而大了的情况。

同理,lx 一开始也要很大很大,不然根本不会被 lx_u + ly_v - mp_{u, v} < slack_v 这个语句更新。

 /* N_Cat x K_crane */
#include <bits/stdc++.h>
#define lowbit(x) ((x) & (-(x)))
using namespace std;

const int N = 510;
const long long INF = 1e10;
int n, m;

int px[N], py[N], rec[N]; long long lx[N], ly[N], mp[N][N], slack[N]; long long sum;
queue < int > q; bool vx[N], vy[N];
inline void get(int pos) // 更新增广路 
{
    int nxt = 0;
    while(pos)
    {
        nxt = px[rec[pos]];
        px[rec[pos]] = pos;
        py[pos] = rec[pos];
        pos = nxt;
    }
}
inline void bfs(int st)
{
    memset(vx, 0, sizeof(vx));
    memset(vy, 0, sizeof(vy));
    for(int i = 1; i <= n; ++i) slack[i] = 3 * INF;
    while(!q.empty()) q.pop(); q.push(st);
    while(true)
    {
        while(!q.empty())
        {
            int pos = q.front(); q.pop(); vx[pos] = true;
            for(int to = 1; to <= n; ++to)
            {
                if(vy[to]) continue; // 之前确定了增广不了就没必要再试一次了 
                if(lx[pos] + ly[to] - mp[pos][to] < slack[to])
                {
                    slack[to] = lx[pos] + ly[to] - mp[pos][to];
                    rec[to] = pos;
                    if(slack[to] == 0) // 在交错树内部
                    {
                        vy[to] = true; 
                        if(!py[to]) {get(to); return ;} // 找到增广路,结束该轮算法
                        else q.push(py[to]); // 对其匹配的左部点继续尝试增广 
                    }
                }
            }
        }
        long long dt = INF;
        for(int i = 1; i <= n; ++i) if(!vy[i]) dt = min(dt, slack[i]);
        for(int i = 1; i <= n; ++i)
        {
            if(vx[i]) lx[i] -= dt; // 在 S 里 
            if(vy[i]) ly[i] += dt; // 在 T 里 
            else slack[i] -= dt; // 在 T' 里 
        }
        for(int i = 1; i <= n; ++i) // 看看调整顶标之后是不是有新的增广方案 
        {
            if(!vy[i] && !slack[i])
            {
                vy[i] = true;
                if(!py[i]) {get(i); return ;}
                else q.push(py[i]);
            }
        }
    }
}

int main()
{
//  freopen("text.in", "r", stdin);
//  freopen("prog.out", "w", stdout);
    ios::sync_with_stdio(false);
    cin.tie(0), cout.tie(0);
    cin >> n >> m;
    for(int i = 1; i <= n; ++i)
        for(int j = 1; j <= n; ++j)
            mp[i][j] = -INF;
    for(int i = 1, x, y, w; i <= m; ++i)
    {
        cin >> x >> y >> w;
        mp[x][y] = max(mp[x][y], 1ll * w);
    }
    for(int i = 1; i <= n; ++i) lx[i] = INF;
    for(int i = 1; i <= n; ++i) bfs(i);
    for(int i = 1; i <= n; ++i) sum += mp[i][px[i]];
    cout << sum << '\n';
    for(int i = 1; i <= n; ++i) cout << py[i] << " "; cout << '\n';
    return 0;
}
/*

*/

其他

KM 算法带给我们的并不只是一个求解二分图的“自动程序”。

正如 Johnson 算法启发我们在网络流中使用势能法优化掉多轮 SPFA 一样,KM 算法带给我们一个结论:

对于这个问题,现在你无需用线性规划的理论去对偶转化它,直接用 KM 算法的组合意义就能证明这个最小值取到以 w 作为邻接矩阵的二分图的最大权值匹配。

学完 KM 之后,你可以尝试一下这道题。