【学习笔记】KM 算法

省流：$\mathrm{KM}$ 算法可以认为是 $\text{zkw}$ 费用流放在二分图最大权匹配上了。我们暂且希望它还能有常数上的优势吧。

壹、二分图

这个我实在不想讲。

1. 定义

可以用两种颜色将点染色，使得任意两个有边直接相连的节点的颜色不同，则该图为二分图。

任取一种染色方案（不妨设染成了黑白两色），染成黑色的点的集合是 $X$ 部（或写作 $\mathfrak{X}$ 部），染成白色的点的集合是 $Y$ 部（或写成 $\mathfrak{Y}$ 部）。

Comment. 有人问我为什么要用旧德式字体，我的回答是：我找不出其他任何地方会用到了。

2. 完美匹配

大小为 $|\mathfrak{X}|=|\mathfrak{Y}|$ 的匹配，就是完美的。

3. 完备匹配

指最大匹配等于 $\min (|\mathfrak{X}|,|\mathfrak{Y}|)$ 的情况。完美匹配是完备匹配的特殊情况。

贰、匈牙利算法

用于解决二分图最大匹配问题。

1. 增广路

任意取一个点 $x\in \mathfrak{X}$，你可以试着去跟 $y\in \mathfrak{Y}$ 匹配。如果 $y$ 已经有匹配了呢？就让原本与 $y$ 匹配的 $p\in \mathfrak{X}$ 另寻新欢！也就是说，如果存在另一个匹配，能够维持现有匹配的大小，并且 $y$ 不在匹配中，那就可以加入 $⟨x,y⟩$ 这条边。而想要 $y$ 不在匹配中，就必须让 $p$ 取消这个原有的匹配 $⟨p,y⟩$ 。

如果所有 $p$ 都别无选择呢？那么 $x$ 是否匹配无所谓。毕竟我们需要的只是最大化匹配数量。$x$ 的匹配必然带来别人匹配的丢失，同时 已匹配的 $y$ 的集合不变，所以 对后续匹配无影响。

如果用 $f(x)$ 表示能否让 $x$ 寻到一个匹配（为零则不可，为正数则可以），用 $p(x)$ 表示 $x$ 与谁匹配，则
$$f(x)=\sum _{⟨x,y⟩\in \mathbb{E}}[p(y)=0]+f[p(y)]$$

此时有一个大问题：转移有环！举个例子：原本 $⟨A,B⟩,⟨C,D⟩$ 玩得好好的，突然来个 $E$ 说：$D$ 给我，$C$ 滚蛋。 然后 $C$ 找 $B$ ，说服 $A$ 找别人。$A$ 看到了 $D$ ，结果又让 $C$ 去更换！

考虑 $f(x)$ 的作用：让某个 $x_0$ 获得匹配的机会。既然如此（好像没什么因果关系），我们 规定转移不能出现环。那么，在递归过程中，已经访问了 $S$ 集合中的点，就不能继续访问了。

转移确实无环了。但此时的 $f$ 似乎会算错……么？如果将 $f(x)$ 的计算视为 $\mathrm{dfs}$ 的计算，那么 $f[p(y)]$ 就是递归子状态。一个很重要的特点是，只要一个 $f\ne 0$ 就可以使得 $f(x)\ne 0$，因为它会从所有子节点得到信息。换句话说，我们只需要找到一个 $p(y)=0$ 的递归出口。所以我们实质上只需要找出 $f(x_0)$ 可以到达的状态有哪些；故而状态转移中的环可以直接忽略。

于是暴力 $\mathtt{dfs}$ 求出 $f(x_0)$ 即可，时间复杂度是 $\mathcal{O}(n+m)$ 的。

求出了 $f(x_0)\ne 0$ 却找不到匹配，悲剧啊！当然你可以把 $\mathrm{dfs}$ 树建出来，然后一路往下。事实上，$\mathrm{dfs}$ 的回溯过程（返回值为 $\mathrm{true}$ 的那条链）中即可更改。

总复杂度 $\mathcal{O}(n^2+nm)$ ，但难以达到 $n^3$，因为稠密图的匹配很容易找到。——其实匈牙利的本质就是 $EK$ 求最大流，而网络流的复杂度是玄学。

2. 代码

由于 $f$ 可以保留，在 $\mathrm{dfs}$ 回溯时不清空 $vis$ 数组。但是更改 $x_0$ 就需要清空了！

bool vis[MAXN]; // need to clear
vector<int> G[MAXN]; // graph
int match[MAXN]; // part Y
bool dfs(int x){
	if(vis[x]) return false;
	vis[x] = true; // remove cycle
	for(auto y : G[x]) // ∀ <x,y>
		if(!match[y] || dfs(match[y]))
			return match[y] = x, true;
	return false; // unable to change
}

而我们注意到，所有会访问到的 $x$ 都是 $match(y)$，除了最初的 $x_0$，并且它一定不会被 $match(y)$ 访问到。所以将 $vis$ 数组定义在 $\mathfrak{Y}$ 部上同样可行。效果是相同的，但是代码短了一行，在 $|\mathfrak{Y}|<|\mathfrak{X}|$ 时可以略微减小常数。

bool dfs(int x){
	for(int y : G[x])
		if(!match[y] || (vis[y] == false
		&& (vis[y] = true) && dfs(match[y])))
			return match[y] = x, true;
	return false;
}

还有 $\mathrm{bfs}$ 写法，据说 性能更好，可是比较麻烦——比如没有了函数递归的隐式栈，我们要用 $pre$ 模拟栈来修改。

void getAugment(int y,int x){
	while(y != 0){
		match[y] = x;
		y = prex[x]; x = prey[y];
	}
}
bool bfs(int x){
	int *fro = que, *bac = fro+1;
	for(prex[*que=x]=0; fro!=bac; ++fro)
		for(int y : G[*fro]) if(!vis[y]){
			if(!match[y]) return getAugment(y,*fro), true;
			prey[y] = *fro, prex[match[y]] = y;
			*(bac ++) = match[y], vis[y] = true;
		}
	return false;
}

3. 例题

• 洛谷板题：用上面的 $\mathrm{dfs}$ 代码即可通过。这里只补充主函数。

int main(){
	int X = readint(), Y = readint();
	int m = readint();
	for(int i=1; i<=m; ++i){
		int x = readint();
		G[x].push_back(readint());
	}
	int ans = 0;
	for(int i=1; i<=X; ++i){
		memset(vis+1,0,X);
		if(dfs(i)) ++ ans;
	}
	printf("%d\n",ans);
    return 0;
}

• 劈配：我特意转载的匈牙利做法，为了让大家更熟悉匈牙利。其实是因为我没学懂网络流。

叁、$\mathrm{KM}$ 算法

解决二分图最大权完美匹配。也就是说，保证是完美匹配时，让边权和最大。

1. 思想

把边权变为点权。为什么可以这样？或许是出于 $\mathrm{LP}$ 对偶吧，我猜。

边权完全等价地放在点上，不可能——否则我们将得出，所有的匹配的权值都是一样的，于是最大权匹配就讲完啦！所以我们只能考虑 “最大化” 背景下的转化。等式行不通，考虑不等式。令 $v(x)$ 为点权，那么必须满足
$$\forall ⟨x,y⟩\in \mathbb{E},v(x)+v(y)⩾w(x,y)$$

这是基本性质。因为这一条可以保证 答案不超过点权和。而我们只需要想办法取到点权和。也就是让满足 $v(x)+v(y)=w(x,y)$ 的边构成完美匹配。称这些边为亮边，其余边为暗边。

模拟 $\mathrm{Hungarian}$ 算法的做法，考虑每次让一个 $\mathfrak{X}$ 部的点获得匹配。如果亮边可以找到增广路，那就直接增广。否则，我们要调整点权，这也是算法的核心。实际上就一句话的事儿：对于 $\mathrm{dfs}$ 中访问过的每一个点，如果属于 $\mathfrak{X}$ 部，将其权值减小 $d(d>0)$，如果属于 $\mathfrak{Y}$ 部，将其权值增大 $d$ 。

来看看，这样操作之后会发生什么。不妨设 $\mathrm{dfs}$ 中访问了集合 $S$ 中的点。对于边 $⟨i,j⟩$，不妨设 $i\in \mathfrak{X}$ 而 $j\in \mathfrak{Y}$，分成四类来考虑。

• 如果 $i\in S,j\in S$ 或者 $i\notin S,j\notin S$，显然没有影响（基本性质仍然满足）。
• 如果 $i\in S,j\notin S$，那么这条边是暗边——因为 $\mathrm{dfs}$ 的对象是 $\mathfrak{X}$ 部的点，这条边是亮边则必然被访问。现在 $v(i)$ 减小了，这条边可能变亮，但是也可能变得不满足基本性质。
• 如果 $i\notin S,j\in S$，天知道这条边是亮还是暗。反正 $v(j)$ 增大会让它暗下去，并且不可能导致非法情况。

为了使得第二类边不出问题（仍然满足基本性质：端点权值和不小于边权），取
$$d=\underset{⟨i,j⟩,i\in S\cap \mathfrak{X},j\notin S}{\min }[v(i)+v(j)-w(i,j)]$$

而后 $S$ 必然增大。因为原有的边不消失，新加入的边一定可以用到。$|S|=2n$ 时，当前点必然获得一个匹配。如果我们每次操作都是 $\mathcal{O}(m)$ 暴力找 $d$ 加上 $\mathcal{O}(n+m)$ 求增广路，一个点要操作 $n$ 次使得 $|S|$ 增大到 $2n$，那么总复杂度 $\mathcal{O}[n\cdot n\cdot (n+m)]=\mathcal{O}(n^3+m{n}^2)$ 的。

你可能有这样的疑惑：一定存在二类边吗？答案是肯定的。否则，这说明 $S\cap \mathfrak{X}$ 的邻点都属于 $S$，但是在 $S$ 中仍然找不到 $S\cap \mathfrak{X}$ 的完备匹配——于是原图中也就不存在完美匹配了。

考虑到 $m=\mathcal{O}(n^2)$，整个算法并不是很快。复杂度瓶颈在于找 $d$ 以及求增广路，两个都和 $m$ 相关。我们尝试优化一下，先从求增广路开始吧！

算法流程中，每次 $|S|$ 只增大 $2$ ，也就是 只加入了一条边，但是我们要 从头开始 $\mathrm{dfs}$ 。这难道是合理的吗？而且终极目标都是给 $x_0$ 找到匹配罢了。上面 $\mathrm{Hungarian}$ 算法中也有提到，我们只需要知道从 $x_0$ 出发能走到哪些 $x\in \mathfrak{X}$ 罢了。于是我们考虑 保留 $\mathrm{dfs}$ 树原有形态，也就是说，$x_0$ 能走到的点还是能走到，只看能不能走到更多的点。

避免了多次增广，已经将一个点的 $\mathcal{O}(n)$ 次操作的总复杂度降至 $\mathcal{O}(n+m)$，好耶！

接下来解决求 $d$ 的问题。用 $\mathrm{SY}(x)$，意思是 $\mathrm{Search}\mathfrak{Y}$ 部，存储每个 $\mathfrak{Y}$ 部的点 作为端点的边的最小值。即
$$\mathrm{SY}(y)=\underset{⟨i,y⟩,i\in S}{\min }[v(i)+v(y)-w(i,y)](y\in \mathfrak{Y})$$

为什么不用 $\mathfrak{X}$ 部的点呢？因为它们随着 $S$ 变大，条件 $j\notin S$ 越发难满足，而 $\min$ 是很难删除的。与此相对的是，$\mathrm{SY}(x)$ 是容易维护的，并且第二类边实际上加入的就是 $\mathfrak{Y}$ 部的点，这与 $\mathtt{prim}$ 比较相似。

此时，我们只需要用 $\mathcal{O}(n)$ 扫一遍 $\mathrm{SY}$ 数组，就能找到我们要新加入的点。加入的点直接开始增广。增广的过程中还可以顺便修改 $\mathrm{SY}$ 的值，一举多得啊！

于是两个含有 $m$ 的操作都不见了！时间复杂度就变成 $\mathcal{O}(n^3+nm)$ 的了，基本上是匈牙利算法的复杂度。

2. 代码

洛谷题解区有许多 $\mathtt{bfs}$ 的代码，欲者自取。这里贴一份 $\mathtt{dfs}$ 的代码吧。

const int INFTY = (1<<30)-1;
const int NO_SUCH_THING = -INFTY;
const int MAXN = 502;
int g[MAXN][MAXN]; ///< if no edge exists, be @a NO_SUCH_THING
int match[MAXN], prex[MAXN], prey[MAXN];
void getAugment(int y,int x){
	while(y) match[y] = x, swap(y,prex[x]), x = prey[y];
}
int vx[MAXN], vy[MAXN], dis[MAXN];
bool vis[MAXN]; // defined on X part
int dfs(int x,const int &n){
	vis[x] = true; int d;
	for(int i=1,o; i<=n; ++i){
		if(g[x][i] == NO_SUCH_THING) continue;
		if(match[i] && vis[match[i]]) continue;
		if((d = vx[x]+vy[i]-g[x][i]) == 0){
			prey[i] = x; if(!match[i]) return i;
			if(!!(o = dfs(match[i],n))) return o;
		}
		else if(dis[i] > d) dis[i] = d, prey[i] = x;
	}
	return 0; // no augment path found
}
int main(){
	int n = readint();
	rep(i,1,n) rep(j,1,n)
		g[i][j] = NO_SUCH_THING;
	for(int m=readint(),a,b; m; --m){
		a = readint(), b = readint();
		g[a][b] = readint();
	}
	rep(i,1,n) rep(j,1,n)
		if(g[i][j] != NO_SUCH_THING)
			vx[i] = max(vx[i],g[i][j]);
	for(int x=1,o; x<=n; ++x){
		memset(vis+1,false,n);
		fill(dis+1,dis+n+1,INFTY);
		for(o=dfs(x,n); o==0; ){
			int_ slack = INFTY;
			rep(i,1,n) if(!vis[match[i]])
				getMin(slack,dis[i]);
			rep(i,1,n) if(vis[i]) vx[i] -= slack;
			rep(i,1,n) // Y part
				if(vis[match[i]]) vy[i] += slack;
				else if(dis[i] != INFTY) dis[i] -= slack;
			for(int i=1; i<=n&&o==0; ++i)
				if(!vis[match[i]] && !dis[i]){
					if(!match[i]) o = i;
					else o = dfs(match[i],n);
				}
		}
		getAugment(o,prey[o]);
	}
}

3. 注解

边权转为点权，还是不等式，这太糟糕了……它到底想干什么？仔细想想 $v(x)+v(y)⩾w(x,y)$，有点像松弛？难道说它是某种 最短路？也就是说——它其实就是最小费用流吗？

从最小费用流的角度考虑，边权是 $-w(x,y)$，令 $v(y)$ 变为原来的相反数，由上式有 $[-w(x,y)]+v(x)-v(y)⩾0$ 。这就是最小费用流的 $\mathtt{dijkstra}$ 方法中的 “势函数” 嘛！

然而 $\mathtt{dijkstra}$ 方法中 $h(x)$ 就是残留网络上从 $S$ 出发的最短路，这里则应该是到汇点 $T$ 的距离。$\mathrm{KM}$ 相当于假定 $dis(x)=v(x)$，然后慢慢调整。用狗狗的话说：“好蠢呐！”

另外一个困扰了我很久的问题：顶标的范围有多大？值得注意的是，每次进行 $v(x)$ 减 $d$ 同时 $v(y)$ 加 $d$ 的操作时，$\sum v(x)+\sum v(y)$ 都减小了 $d$，因为 $x$ 的数量比 $y$ 多一个。显然 $\sum d=\mathcal{O}(nw)$，其中 $w$ 是边权的值域。那么 $v(x)$ 每次都减小 $d$，就可能达到 $nw$ 的修改量！

而这个上界是可以达到的——只需要让一条增广路径经过 $(n-1)$ 条边权最大的边。也就是 $\mathfrak{X}$ 部的 $i$ 号点向 $\mathfrak{Y}$ 部的 $(i+1)$ 号点连边，权值为 $\max v$，它们会最初作为匹配边；然后 $\mathfrak{X}$ 部的 $i$ 号点向 $\mathfrak{Y}$ 部的 $i$ 号点连边，权值为 $\min v$，为了达到完美匹配，不得不使用这些边。显然此时到 $T$ 的最短路长度是 $\mathcal{O}(nw)$，于是顶标也是这个范围了！

另外说一句：当 $w\ll \max v$ 时，顶标的瓶颈是 $\mathcal{O}(\max v)$ 的初始顶标。解决方案是将每条边的权值都减去 $\min v$ 。

4. 推广

这个 $\mathrm{KM}$ 算法的限制太紧了——必须是完美匹配，这在很多情况下都难以实现。来看看这两个拓展：

• 求原图的最大权值匹配（所有匹配中权值最大的一个）。
• 求原图的最大权值最大匹配（最大匹配中权值最大的一个）。

对于第一个，不妨设 $|\mathfrak{X}|⩾|\mathfrak{Y}|$，把 $\mathfrak{Y}$ 部补充一些虚点，使得 $|\mathfrak{X}|=|\mathfrak{Y}|$ 。然后用权值为 $0$ 的边将图补充成完全图。于是我们就有完美匹配啦！由于权值是 $0$，肯定是最后被选择的，仅用来补全完美匹配；负数边权可以直接扔掉嘛。

对于第二个，操作方式更神奇。与上面相同的补边，但是边权要设置为 $-\infty$ 才行！只有这样，我们才会心甘情愿地选出最大匹配——总的边数是定值，存在的边构成的匹配数越大，$-\infty$ 即不存在的边数量越少。似乎 可以在上面的代码上做微调，鉴于我从未实现过，不妄加评论。