KM算法模板

KM算法,在此不涉及科学的证明,仅仅是粗略的说明算法的正确性。

题目链接

Uva1411-Ants

题意

第一个参数是代表点的对数。
接下来n行表示蚂蚁的坐标
接下来n行表示apple tree的坐标
要求第i个蚂蚁对应的苹果树的序号,要求任意一种连线不能重合的方案即可

题解

不妨将蚂蚁抽象成白点,苹果树抽象成黑点。黑点和白点构成了二分图的关系。两点间边的权值对应着欧几里得距离。
那么为什么使用KM算法是正确的呢?
假设最佳完美匹配中,存在两条相交的直线\(a_1-b_1\)和\(a_2-b_2\)。那么dist(a1, b1)+dist(a2, b2)一定大于dist(a1, b2)+dist(a2, b1);
而后者是不相交的。而KM算法算的是权值和最大的值.但是我们可以将原理的权值转换为负值来求最大值,再取反从而求得了最小值

AC 代码

说明:这是我在网上找到的一个复杂度为O(\(n^3\))的算法,代码风格不是很好。
去除了原来算法中的slack[]变量。

左右顶标(ex_girl[]和ex_boy[])的和就是最佳完美匹配的权值和
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
#include <iostream>
#include <cstring>
#include <cstdio>
#include <cmath>
using namespace std;
const int MAXN = 110;
const double INF = 1e9;
const double eps = 1e-6;
double love[MAXN][MAXN]; // 记录每个妹子和每个男生的好感度
double ex_girl[MAXN]; // 每个妹子的期望值
double ex_boy[MAXN]; // 每个男生的期望值
bool vis_girl[MAXN]; // 记录每一轮匹配匹配过的女生
bool vis_boy[MAXN]; // 记录每一轮匹配匹配过的男生
int match[MAXN]; // 记录每个男生匹配到的妹子 如果没有则为-1
//double slack[MAXN]; // 记录每个汉子如果能被妹子倾心最少还需要多少期望值
int N;
double get_dis(int x1, int y1, int x2, int y2){
return sqrt((double)(x2-x1)*(x2-x1)+(double)(y2-y1)*(y2-y1));
}
bool dfs(int girl)
{
vis_girl[girl] = true;
for (int boy = 0; boy < N; ++boy) {
if (vis_boy[boy]) continue; // 每一轮匹配 每个男生只尝试一次
double gap = ex_girl[girl] + ex_boy[boy] - love[girl][boy];
if (fabs(gap) <= eps) { // 如果符合要求
vis_boy[boy] = true;
if (match[boy] == -1 || dfs( match[boy] )) { // 找到一个没有匹配的男生 或者该男生的妹子可以找到其他人
match[boy] = girl;
return true;
}
}
}
return false;
}
double KM()
{
memset(match, -1, sizeof(match)); // 初始每个男生都没有匹配的女生
memset(ex_boy, 0, sizeof(ex_boy)); // 初始每个男生的期望值为0
// 每个女生的初始期望值是与她相连的男生最大的好感度
for (int i = 0; i < N; ++i) {
ex_girl[i] = love[i][0];
for (int j = 1; j < N; ++j) {
ex_girl[i] = max(ex_girl[i], love[i][j]);
}
}
// 尝试为每一个女生解决归宿问题
for (int i = 0; i < N; ++i) {
//fill(slack, slack + N, INF); // 因为要取最小值 初始化为无穷大
while (1) {
// 为每个女生解决归宿问题的方法是 :如果找不到就降低期望值,直到找到为止
// 记录每轮匹配中男生女生是否被尝试匹配过
memset(vis_girl, false, sizeof(vis_girl));
memset(vis_boy, false, sizeof(vis_boy));
if (dfs(i)) break; // 找到归宿 退出
// 如果不能找到 就降低期望值
// 最小可降低的期望值
double d = INF;
for(int i=0; i<N; i++)if(vis_girl[i]){
for(int j=0; j<N; j++){
if(!vis_boy[j]){
d = min(d, ex_girl[i]+ex_boy[j]-love[i][j]);
}
}
}
for (int j = 0; j < N; ++j) {
// 所有访问过的女生降低期望值
if (vis_girl[j]) ex_girl[j] -= d;
// 所有访问过的男生增加期望值
if (vis_boy[j]) ex_boy[j] += d;
// 没有访问过的boy 因为girl们的期望值降低,距离得到女生倾心又进了一步!
}
}
}
double ans = 0;
for(int i=0; i<N; i++){
ans+=(ex_boy[i]+ex_girl[i]);
}
return ans;
}
int main()
{
int kase=0;
int pos[MAXN][2];
int pos1[MAXN][2];
while (scanf("%d", &N) == 1){
if(++kase>1)printf("\n");
for(int i=0; i<N; i++){
scanf("%d%d", &pos[i][0], &pos[i][1]);
}
for(int i=0; i<N; i++){
scanf("%d%d", &pos1[i][0], &pos1[i][1]);
}
for (int i = 0; i < N; ++i)
for (int j = 0; j < N; ++j)
love[i][j] = -1.0*get_dis(pos1[i][0], pos1[i][1], pos[j][0], pos[j][1]);
KM();
for(int i=0; i<N; i++){
printf("%d\n", match[i]+1);
}
}
return 0;
}

KM算法

原理及其正确性

算法模拟的详细过程

首先我们定义两个概念:

  • 可行顶标是一个节点函数l,使得对任意的弧(x, y), 有l(x)+l(y)\(\geq\)w(x, y)。
    进一步的说明:l(x)和l(y)函数的值都是人为的自由的设定的, w(x, y)是边的权值,题目中给定。

  • 相等子图是G的生成子图,包含所有的点,但是只包含l(x)+l(y) = w(x, y)的所有的边(x, y)。

结论:如果相等子图有完美匹配, 那么该匹配是原图的最大权匹配

证明上面的结论:

  1. 设\(M^* \)是相等子图的完美匹配,那么\(M^*\)的边权和等于点的权值和(顶标的和);
  2. 任取G的一个完美匹配M,由于M中边满足不等式w(x, y)\(\geq\)l(x)+l(y),M的边权的和不能超过所有的顶标之和
  3. 综合前面两点,完美匹配不能超过顶标之和(第二条的规定保证),并且能够取到等号(第一条的规定保证),那么顶标之和就是权值最大了。

模板代码

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
// LA4043 Ants
// Rujia Liu
#include <cstdio>
#include <cstring>
#include <cmath>
#include <algorithm>
using namespace std;
const int maxn = 100 + 10;
const double INF = 1e9;
int n;
double W[maxn][maxn];
double Lx[maxn], Ly[maxn]; // 顶标
int left[maxn]; // left[i]为右边第i个点的匹配点编号
bool S[maxn], T[maxn]; // S[i]和T[i]为左/右第i个点是否已标记
bool eq(double a, double b) {
return fabs(a-b) < 1e-6;
}
bool match(int i){
S[i] = true;
for(int j = 1; j <= n; j++) if (eq(Lx[i]+Ly[j], W[i][j]) && !T[j]){
T[j] = true;
if (!left[j] || match(left[j])){
left[j] = i;
return true;
}
}
return false;
}
void update(){
double a = INF;
for(int i = 1; i <= n; i++) if(S[i])
for(int j = 1; j <= n; j++) if(!T[j])
a = min(a, Lx[i]+Ly[j] - W[i][j]);
for(int i = 1; i <= n; i++) {
if(S[i]) Lx[i] -= a;
if(T[i]) Ly[i] += a;
}
}
void KM() {
for(int i = 1; i <= n; i++) {
left[i] = Lx[i] = Ly[i] = 0;
Lx[i] = W[i][1];
for(int j = 1; j <= n; j++)
Lx[i] = max(Lx[i], W[i][j]);
}
for(int i = 1; i <= n; i++) {
for(;;) {
for(int j = 1; j <= n; j++) S[j] = T[j] = 0;
if(match(i)) break; else update();
}
}
}
int main(){
int kase = 0;
while(scanf("%d", &n) == 1) {
if(++kase > 1) printf("\n");
int x1[maxn], y1[maxn], x2[maxn], y2[maxn];
for(int i = 1; i <= n; i++) scanf("%d%d", &x1[i], &y1[i]);//ants
for(int i = 1; i <= n; i++) scanf("%d%d", &x2[i], &y2[i]);//apple tree
for(int i = 1; i <= n; i++) // ant colony
for(int j = 1; j <= n; j++) // apple tree
W[i][j] = -sqrt((double)(x1[j]-x2[i])*(x1[j]-x2[i]) + (double)(y1[j]-y2[i])*(y1[j]-y2[i]));
// for(int i=0; i<n; i++)
//// {
//// for(int j=0; j<n; j++){
//// printf("%lf ", W[i][j]);
//// }
//// printf("\n");
//// }
KM(); // 最大权匹配
double ans = 0;
for(int i=1; i<=n; i++) ans+=(Lx[i]+Ly[i]);
for(int i = 1; i <= n; i++) printf("%d\n", left[i]);
}
return 0;
}

hdu 6346 更快的km模板

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
#include<bits/stdc++.h>
using namespace std;
#define LL long long
const int inf = 0x3f3f3f3f;
const LL INF = 0x3f3f3f3f3f3f3f3f;
const int N = 210;
int val[N][N];
LL lx[N],ly[N];
int linky[N];
LL pre[N];
bool vis[N];
bool visx[N],visy[N];
LL slack[N];
int n;
void bfs(int k) {
LL px, py = 0,yy = 0, d;
memset(pre, 0, sizeof(LL) * (n+2));
memset(slack, inf, sizeof(LL) * (n+2));
linky[py]=k;
do {
px = linky[py],d = INF, vis[py] = 1;
for(int i = 1; i <= n; i++)
if(!vis[i]) {
if(slack[i] > lx[px] + ly[i] - val[px][i])
slack[i] = lx[px] + ly[i] - val[px][i], pre[i] = py;
if(slack[i] < d)
d = slack[i], yy = i;
}
for(int i = 0; i <= n; i++)
if(vis[i])
lx[linky[i]] -= d, ly[i] += d;
else
slack[i] -= d;
py = yy;
}
while(linky[py]);
while(py)
linky[py] = linky[pre[py]], py=pre[py];
}
void KM() {
memset(lx, 0, sizeof(int)*(n+2));
memset(ly, 0, sizeof(int)*(n+2));
memset(linky, 0, sizeof(int)*(n+2));
for(int i = 1; i <= n; i++)
memset(vis, 0, sizeof(bool)*(n+2)), bfs(i);
}
int main() {
int T;
scanf("%d", &T);
for(int _i = 1; _i <= T; _i++) {
scanf("%d", &n);
for(int i = 1; i <= n; i++) {
for(int j = 1; j <= n; j++) {
scanf("%d", &val[i][j]);
val[i][j] = -val[i][j];
}
}
KM();
LL ans = 0;
for(int i = 1; i <= n; ++i)
ans += lx[i] + ly[i];
printf("Case #%d: %I64d\n", _i, -ans);
}
return 0;
}

适用的场景范围

题目中的二分图的两个集合中元素的个数要相同,即要首先要存在完美匹配。且两点间一定要有刻量的边权关系。(比如好感度,若两者间没有好感,那么边权赋值为0,不能是无法连接的状态。)

未解决的问题

KM算法还可以用费用来实现,这样的思想还不会。
适应的场景范围需要进一步的刷题求证。

文章目录
  1. 1. 题目链接
  2. 2. 题意
  3. 3. 题解
  4. 4. AC 代码
  5. 5. KM算法
    1. 5.1. 原理及其正确性
    2. 5.2. 模板代码
  6. 6. hdu 6346 更快的km模板
  7. 7. 适用的场景范围
  8. 8. 未解决的问题
{{ live2d() }}