[学习笔记] KM算法

前言

这个东西学了我挺久了，我先奉劝各位一定要先搞清楚匈牙利算法。感谢 \(\tt jzm\) 巨佬对我耐心的讲解，因为我太弱了所以卡了很久都不懂。如果你有任何问题请在本篇博客下面留言，我会尽力解答的。

\(\tt KM\) 算法主要用来解决最大权完美匹配，因其稳定的 \(O(n^3)\) 可以吊打玄学费用流，所以出匹配的题时无脑写费用流可能被卡。但是它也只能解决匹配问题，费用流的应用却极其广泛。

我会着重讲解代码，因为我觉得 \(\tt KM\) 的代码才是最难懂的（而且网上的代码很多都有点问题）

思想

神奇的思想，关键是 顶标 这个设计。所谓顶标也就是二分图上的点权，我们把边权转化成点权，设 \(X\) 部的点的顶标为 \(vx[i]\) ，\(Y\) 部的点的顶标为 \(vy[i]\) ，记图上的边为 \(a[i][j]\) ，那么对于任意的 \(i,j\) ，我们要一直保证：

\[vx[i]+vy[j]\geq a[i][j] \]

那么对于 \(a[i][j]=vx[i]+vy[j]\) 的边，贪心地选他是最优的。我们成这样的边为 相等边 ，称相等边连接而成的子图为 相等子图 ，那么当相等子图的最大匹配是完美匹配时，说明我们找到了答案，答案是顶标之和。

\(\tt KM\) 算法和匈牙利算法类似，每次都会进行 增广 操作，如果可以直接用相等边增广的话，自然是极好的。如果不能进行增广操作呢？就说明现在的相等子图不行了，我们要 扩大相等子图 。

怎么扩大相等子图？其实就是修改顶标，顶标是可以按照我们的愿望修改的，但是注意要一直满足 \(vx[i]+vy[j]\geq a[i][j]\) ，设我们现在增广到的集合是 \(S\)（既有 \(X\) 部也有 \(Y\) 部），记 \(i\) 是 \(X\) 部的点，\(j\) 是 \(Y\) 部的点，那么对于 \(i\in S\) ，我们使他的顶标减少 \(d\) ，对于 \(j\in S\) ，我们使他的顶标增加 \(d\) ，考虑这样做的影响：

• 对于 \(x\in S,y\in S\) ，\(vx[x]+vy[y]\) 不变，所以还是满足\(vx[x]+vy[x]\geq a[x][y]\)
• 对于 \(x\in S,y\not\in S\) ，这条边 \(vx[x]+vy[y]\) 会减少 \(d\) ，这条边可能就会成为新的相等边。
• 对于 \(x\not\in S,y\in S\) ，这条边 \(vx[x]+vy[y]\) 会增加 \(d\) ，这条边依然不是相等边，但是我们用不到他。
• 对于 \(x\not\in S,y\not\in S\) ，无影响，不用管。

综上所诉，我们可以知道 \(d\) 应该这样取值：

\[d=\min_{x\in S,y\not\in S} vx[x]+vy[y]-a[x][y] \]

当我们修改顶标之后获得新的相等边之后就可以继续增广了，而且那个限制条件还是一直满足的，挺不错。注意一开始为了满足条件我们应该让 \(vx[i]=\max_j a[i][j]\) ，\(vy[i]=0\)

代码

就算懂了思想之后可能还是难以敲出代码，先给出一些数组的定义：

• \(vx[i],vy[i]\) ，顶标数组，如上所述。
• \(xtoy[i],ytox[i]\) ，表示和 \(X/Y\) 部的点匹配的是哪个点，如果没有匹配点就是 \(0\)
• \(visx[i],visy[i]\) ，表示在当前增广过程中这个点是否被访问到。
• \(pre[i]\) 表示 \(Y\) 部的点想匹配其的 \(X\) 部的点 ，理解这个数组很重要，就是匈牙利算法的思想，\(X\) 部的点 \(a\) 虽然想匹配他，但是碍于这个点已经有了一个匹配的 \(X\) 部的点 \(b\) ，所以我们先记录下这个数组，再去增广 \(b\) ，因为写法是 \(\tt bfs\) 的，所以用这个数组来达到 \(\tt dfs\) 的功能。
• \(sy[i]\) 表示还没有访问的 \(Y\) 部的点的 \(\min_{x_\in S}vx[x]+vy[y]-a[x][y]\) ，这个数组方便算 \(d\) ，而且维护起来也不是很麻烦，在每次 \(X\) 部的点加入 \(S\) 集合时修改一下它就行了。

请先完全理解上面的定义之后再读我带注释的代码，就是这道题【模板】二分图最大权完美匹配

如果你读懂了代码，那么不难知道复杂度是 \(O(n^3+nm)\) ，因为均摊下来花费 \(O(n)\) 的时间可以使一条边变成相等边，最差情况当所有边变成相等边时一定可以得到完美匹配，\(O(n^3)\) 就是每次枚举增广哪个点，二分图的每个点又要扩展一次，扩展的消耗时 \(O(n)\) ，所以总共时 \(O(n^3)\) 的。

#include <cstdio>
#include <cstring>
#include <iostream>
#include <queue>
using namespace std;
#define int long long
const int inf = 1e18;
const int M = 505;
int read()
{
	int x=0,f=1;char c;
	while((c=getchar())<'0' || c>'9') {if(c=='-') f=-1;}
	while(c>='0' && c<='9') {x=(x<<3)+(x<<1)+(c^48);c=getchar();}
	return x*f;
}
int n,m,vx[M],vy[M],xtoy[M],ytox[M];
int visx[M],visy[M],pre[M],sy[M],a[M][M];
void aug(int y)
//也就是y这个点是空闲的，可以把以前记录的匹配修改了，可以增加一个匹配 
{
	for(int nxt,x;y;y=nxt)//我觉得这个函数很像"链式反应" 
	{
		nxt=xtoy[x=pre[y]];
		//nxt表示x的原配y，由于x已经找到了新的y所以原配就被解放出来了 
		ytox[y]=x;xtoy[x]=y;//修改匹配关系 
	}
}
void jzm(int x)//每次都只解决x的问题，就和匈牙利很类似啊 
{
	queue<int> q;
	memset(visx,0,sizeof visx);//一定要清空 
	memset(visy,0,sizeof visy);//一定要情况 
	for(int i=1;i<=n;i++) sy[i]=inf;//一开始都是最大值，注意一下定义哦 
	q.push(x);
	while(1)//直到找到了x的匹配再离开 
	{
		while(!q.empty())
		{
			int u=q.front();q.pop();
			visx[u]=1;
			for(int i=1;i<=n;i++)
			{
				if(visy[i] || a[u][i]==-inf) continue;
				//如果已经在S中或者边没有，不能访问这个y 
				int d=vx[u]+vy[i]-a[u][i];//算d 
				if(d==0)//直接就是相等边了 
				{
					visy[i]=1;pre[i]=u;//这个点就被增广过了，点u是想匹配i的，所以修改pre 
					if(!ytox[i]) {aug(i);return ;}//这个点空闲，可以直接达到目的，做完之后就走了 
					else q.push(ytox[i]);//否则继续从它的原配开始找 
				}
				else if(sy[i]>d)//不是相等边，所以可以修改sy 
					sy[i]=d,pre[i]=u;//如果这个相等边被激活了，那么u就想匹配i，多打了标记也没关系 
			}
		}
		int del=inf;//求那个min(....)，就是修改量 
		for(int i=1;i<=n;i++)
			if(!visy[i]) del=min(del,sy[i]);//只能在没有访问的y里面找 
		for(int i=1;i<=n;i++)//这里就是按照定义修改顶标 
		{
			if(visx[i]) vx[i]-=del; 
			if(visy[i]) vy[i]+=del;
			else if(sy[i]!=inf) sy[i]-=del;
			//如果sy有意义，那么他一定被哪个x访问过，那么由于vx的减少，他就一定会减小 
		}
		for(int i=1;i<=n;i++)
			if(!visy[i] && !sy[i])//再修改顶标之后成为了相等边 
			{
				visy[i]=1;//和上面类似的操作，扩展 
				if(!ytox[i]) {aug(i);return ;}
				else q.push(ytox[i]);
			}
	}
}
int KM()
{
	int ans=0;
	for(int i=1;i<=n;i++)
		for(int j=1;j<=n;j++)
			vx[i]=max(vx[i],a[i][j]);//初始的vx 
	for(int i=1;i<=n;i++)//每次就增广一个点 
		jzm(i);
	for(int i=1;i<=n;i++)
		ans+=vx[i]+vy[i];//最后的答案一定是这个 
	return ans;
}
signed main()
{
	n=read();m=read();
	for(int i=1;i<=n;i++)
		for(int j=1;j<=n;j++)
			a[i][j]=-inf;//注意要初始化 
	for(int i=1;i<=m;i++)
	{
		int u=read(),v=read(),c=read();
		a[u][v]=max(a[u][v],c);//选最大的边 
	}
	printf("%lld\n",KM());
	for(int i=1;i<=n;i++)
		printf("%lld ",ytox[i]);
}

应用

注意我再 前言 中说它可以解决匹配问题，所以不只能解决最大权完美匹配问题哦，可以看看 oneindark 巨佬博客的最后一部分，然后如果我做到了更多关于 \(\tt KM\) 的题是会补充到这里的：