Kuhn-Munkres 算法

本文接着上一篇博客《匈牙利算法》，前篇讨论的是不带权的二分图求最大匹配，本篇讨论带权的二分图求最大权匹配 ——KM 算法。实际上前者可以看作是后者的特殊情况，即权值只有 0 和 1 的情形。

几个基本事实

• WLOG，通过人为添加权值为 0 的边，任一个二分图 G 都可以被改造成完全图而不影响我们所讨论的问题，因此以下的讨论全都视作完全二分带权图（只考虑非负权值）。

• 如前篇所述，可以通过取反增广路的操作使得匹配边数 + 1

• 最大权匹配一定是完备匹配（如果最大权匹配不完备的话，说明 X 和 Y 至少各有一个点未被匹配，加上这条边之后权和不可能更小）

可行顶标与相等子图

顶标：顶标是一个结点函数 l:V→R

可行顶标：可行顶标是满足下面不等式的顶标（即边权不能超过两端点的顶标和） l(x)+l(y)≥w(x,y),∀x∈X,y∈Y 相等子图：图 G=(V,E) 的（相对于顶标 l 的）相等子图 Gl=(V,El) 包含图 G 的所有顶点，但只包含边权等于两端点顶标和的边。即 El=(x,y):l(x)+l(y)=w(x,y)

Kuhn-Munkres 定理：如果 l 是可行顶标，M 是相等子图 Gl 里的完美匹配，那么 M 是图 G 的最大权匹配。

证明：用 e=(ex,ey) 来记图中的边，ex 和 ey 分别是两端点。任取图 G 的一个完美匹配 M′，容易的到下面的不等式 w(M′)=∑e∈M′w(e)≤∑e∈M′l(ex)+l(ey)=∑v∈Vl(v) 这个不等式说明任何一个完美匹配权和的上界就是所有顶标和 ∑v∈Vl(v)，现在 M 是相等子图 Gl 的完美匹配，从而也是原图 G 的完美匹配，所以有 w(M)=∑e∈Mw(e)=∑v∈Vl(v) 最大权匹配一定是一个完备匹配，易知所有完备匹配的权和上界也是 ∑v∈Vl(v)，所以 M 就是最大权匹配。

• KM 定理将求最大权匹配问题转化成了求完美匹配问题，这是组合优化问题里经典的技巧

• 注意 KM 定理的证明里实际上证明了任意一个匹配 M 和任意一个可行顶标 l 都满足这个不等式

w(M)≤∑v∈Vl(v)

• 和最大流最小割定理有哪些相似之处？

1. 给定初始可行顶标 l 和初始匹配 M
2. 当 M 不是 El 的完美匹配时，重复以下操作
   1. 在相等子图 El 里寻找增广路，然后通过取反操作使匹配数 + 1
   2. 如果找不到增广路，将可行顶标 l 改进为 l′ 使得 El⊂El′，跳转至 1

• 如果 G 有完美匹配，因为每次迭代要么匹配数 + 1 要么相等子图扩大，最终算法一定会停止，而且停止时 M 是 El 的完美匹配。根据 KM 定理，M 就是 G 的最大权匹配。

初始可行顶标和初始匹配

最简单的可行顶标就是取 X 的顶标为邻接边的最大权值，Y 的顶标取 0 ∀y∈Y,l(y)=0∀x∈X,l(x)=maxy∈Y{w(x,y)} 这样显然满足可行顶标的要求 ∀x∈X,y∈Y,w(x)≤l(x)+l(y) 最简单的初始匹配就是 M=ϕ

改进可行顶标

设 l 是一个可行顶标，顶点 u∈V，点集 S⊆V，定义点 u 的邻接集 Nl(u) 和 S 的邻接集 Nl(S)Nl(u)={v:(u,v)∈El}Nl(S)=⋃u∈SNl(u)

引理：设 S⊆X,T=Nl(S)≠Y，令 αl 为

αl=minx∈S,y∉T{l(x)+l(y)−w(x,y)} 令顶标 l′ 为 l′(v)={l(v)−αl,v∈Sl(v)+αl,v∈Tl(v),otherwise 那么 l′ 也是一个可行顶标，并且

1. 若 (x,y)∈El,x∈S,y∈T，那么 (x,y)∈El′
2. 若 (x,y)∈El,x∉S,y∉T，那么 (x,y)∈El′
3. 满足 x∈S,y∉T 的边 (x,y) 可能属于新的相等子图 El′

• 这个引理说明：通过这种方式改进可行顶标，能够使原先在相等子图里的边仍然在相等子图里，有一端在 S 里，另一端不在 T 里的边有可能被纳入新的相等子图

1. 初始化可行顶标 l 和相等子图 El 的初始匹配 M∀y∈Y,l(y)=0∀x∈X,l(x)=maxy∈Y{w(x,y)}M=ϕ

2. 若 M 已经是完美匹配，停止算法；否则，选取未匹配点 u∈X，令 ，S={u}，T=ϕ

3. 若 Nl(S)=T，更新可行顶标 αl=minx∈S,y∉T{l(x)+l(y)−w(x,y)}l′(v)={l(v)−αl,v∈Sl(v)+αl,v∈Tl(v),otherwise

4. 若 Nl(S)≠T，选取 y∈Nl(S)−T

   1. 若 y 是未匹配点，那么从 u 到 y 就找到了一条增广路，取反使匹配数 + 1，跳转至 2
   2. 若 y 是匹配点，设 y 目前与 z 匹配，令 S=S∪{z},T=T∪{y}，跳转至 3

算法正确性

由引理保证了算法在执行过程一致保持着顶标的可行性，算法通过改变可行顶标从而逐渐扩大相等子图。每次找到一个 Y 里的未匹配点，就使得匹配数 + 1。最终算法会在达到完美匹配的时候终止，根据 KM 定理，就得到了原图的最大权匹配。

时间复杂度分析

在具体实现中，通常为每个 y 结点定义一个松弛度存放在数组里，这样可以减少 αl 的计算量 slacky=minx∈S{l(x)+l(y)−w(x,y)} 首先将算法运行的过程分为 |M| 个阶段，使匹配数 + 1 的若干次迭代记作一个阶段，每个阶段都经历了若干次步骤 2、3、4，下面来具体分析每个步骤的时间复杂度。

步骤 2：每个阶段都只会经历一次步骤 2。初始化 S 和 T 显然是 O(1) 的，但此时还要初始化 slack 数组。因为 S 里只有 1 个点，只需要 |S| 次计算就可初始化完毕，所以步骤 2 总共是 O(|V|) 的。

步骤 3：计算 αl=miny∈Tslacky 这一步是 O(|T|) 的，更新可行顶标 l′ 是 O(|T|+|S|) 的，更新完可行顶标之后还要重新更新 slack 数组，根据定义不难推出这里只需要执行 ∀y∉T,slacky=slacky−αl，所以是 O(|Y−T|) 的。要注意的是步骤 3 在每个阶段可能经历多次，但最多只会经历 O(|Y|) 次（因为每经历一次步骤 4 如果需要回到步骤 3 的话，一定使 T 增加一个点）。上面的复杂度全都可以放缩成 O(|V|)，所以步骤 3 总共是 O(|V|2) 的。

步骤 4：第 1 种情况在每个阶段只会经历一次，取反操作是 O(|V|) 的。第 2 种情况同步骤 3 的分析，至多经历 O(|V|) 次，更新集合 S 和 T 只需要 O(1) 的时间。但由于集合 S 多加了一个点，需要更新 slack 数组，∀y∉T,slacky=min{slacky,l(z)+l(y)−w(z,y)}，这需要 O(|Y|) 的时间。步骤 4 总共是 O(|V|2) 的。

最后，因为匹配数 |M| 不可能超过 |V|/2，所以阶段数也是 O(|V|) 的，得出整个算法的时间复杂度是 O(|V|3)

cpp 实现

由于时间紧迫，用松弛度数组优化后的 O(|V|3) 的算法暂时先挖个坑以后再填。这里给出原始的 O(|V|4) 实现。

用邻接矩阵存储二分图，两个数组存储顶标，match 数组用来记录匹配（记录增广路），S 和 T 两个数组来记录算法里的两个集合 S 和 T

int Weight[maxm][maxn];
int Lx[maxm], Ly[maxn]; // 顶标
int match[maxn];    // 记录匹配
bool S[maxm], T[maxn];  // 算法中的两个集合S和T

步骤 1 的初始化可行顶标和初始化匹配写成一个函数

void Init()
{
    // 将X集合的顶标设为最大边权，Y集合的顶标设为0
    for (int i = 1; i <= m; i++)
    {
        Lx[i] = 0;
        for (int j = 1; j <= n; j++)
        {
            match[j] = 0;   // match记录的是Y集合里的点与谁匹配
            Ly[j] = 0;
            Lx[i] = max(Lx[i], Weight[i][j]);
        }
    }
}

步骤 3 的更新顶标就直接用没优化前的原始实现了

void update() // 更新顶标
{
    // 计算a
    int a = 1 << 30;
    for (int i = 1; i <= m; i++)
        if (S[i])
            for (int j = 1; j <= n; j++)
                if (!T[j])
                    a = min(a, Lx[i] + Ly[j] - Weight[i][j]);

// 修改顶标
    for (int i = 1; i <= m; i++)
        if (S[i])
            Lx[i] -= a;
    for (int j = 1; j <= n; j++)
        if (T[j]) 
            Ly[j] += a;
}

步骤 4 可以用递归来实现，写成了 DFS（这里和前篇的匈牙利算法实现非常相似，因为都是在寻找增广路）

bool findPath(int i)
{
    S[i] = true;
    for (int j = 1; j <= n; j++)
    {
        if (Lx[i] + Ly[j] == Weight[i][j] && !T[j]) // 找出在相等子图里又还未被标记的边
        {
            T[j] = true;
            if (!match[j] || findPath(match[j])) // 未被匹配，或者已经匹配又找到增广路
            {
                match[j] = i;
                return true;
            }
        }
    }
    return false;
}

然后就是整个算法的合写

void KM()
{
    Init();

for (int i = 1; i <= m; i++)
    {
        while (true)
        {
            for (int i = 1; i <= m; i++)
                S[i] = 0;
            for (int j = 1; j <= n; j++)
                T[j] = 0;
            if (!findPath(i))
                update();
            else
                break;
        }
    }
}

最后调用即可，完成后查看 match 数组就能查看匹配结果

受益匪浅的一篇博客

另一篇受益匪浅的博客

又一篇受益匪浅的博客

形象描述算法运行过程

线性规划的观点值得一看

一篇介绍的非常详细的博客

本文撰写时的主要参考来源

一篇 Murray State University 的 KM 算法教程

最后放上 Kuhn 和 Munkres 的论文

[1] Munkres, James. "Algorithms for the assignment and transportation problems." Journal of the society for industrial and applied mathematics 5.1 (1957): 32-38.

[2] Kuhn, Harold W. "The Hungarian method for the assignment problem." Naval research logistics quarterly 2.1‐2 (1955): 83-97.