km算法

浏览数：19 / 时间：2015年06月08日

今天花了些时间学了下km算法看了下代码有点大概思路，还是要多做题；

KM算法求最小权二分匹配，模板题，构图很简单，直接把人当作左边的点，房子当作右边的点，
两者之间的曼哈顿距离当作权值即可。第一次搞带权二分匹配的题，就是用KM算法求最小权的时候要加个处，由于KM求的是最大权，
所以在套模板之前把权值都取下相反值最后再把KM算法求出来的最大权值取反即可。

Kuhn－Munkras算法流程：
　　(1)初始化可行顶标的值
　　(2)用匈牙利算法寻找完备匹配
　　(3)若未找到完备匹配则修改可行顶标的值
　　(4)重复(2)(3)直到找到相等子图的完备匹配为止

引用：
KM算法是通过给每个顶点一个标号（叫做顶标）来把求最大权匹配的问题转化为求完备匹配的问题的。
设顶点Xi的顶标为A[i]，顶点Yi的顶标为B [i]，顶点Xi与Yj之间的边权为w[i,j]。在算法执行过程中的任一时刻，对于任一条边(i,j)，A[i]+B[j]>=w[i,j]始终成立。

KM算法的正确性基于以下定理：
　　若由二分图中所有满足A[i]+B[j]=w[i,j]的边(i,j)构成的子图（称做相等子图）有完备匹配，那么这个完备匹配就是二分图的最大权匹配。
　　这个定理是显然的。因为对于二分图的任意一个匹配，如果它包含于相等子图，那么它的边权和等于所有顶点的顶标和；
如果它有的边不包含于相等子图，那么它的边权和小于所有顶点的顶标和。所以相等子图的完备匹配一定是二分图的最大权匹配。
　　初始时为了使A[i]+B[j]>=w[i,j]恒成立，令A[i]为所有与顶点Xi关联的边的最大权，B[j]=0。如果当前的相等子图没有完备匹配，就按下面的方法修改顶标以使扩大相等子图，
直到相等子图具有完备匹配为止。
　　我们求当前相等子图的完备匹配失败了，是因为对于某个X顶点，我们找不到一条从它出发的交错路。
这时我们获得了一棵交错树，它的叶子结点全部是X顶点。现在我们把交错树中X顶点的顶标全都减小某个值d，Y顶点的顶标全都增加同一个值d，那么我们会发现：
两端都在交错树中的边(i,j)，A[i]+B[j]的值没有变化。也就是说，它原来属于相等子图，现在仍属于相等子图。
两端都不在交错树中的边(i,j)，A[i]和B[j]都没有变化。也就是说，它原来属于（或不属于）相等子图，现在仍属于（或不属于）相等子图。
X端不在交错树中，Y端在交错树中的边(i,j)，它的A[i]+B[j]的值有所增大。它原来不属于相等子图，现在仍不属于相等子图。
X端在交错树中，Y端不在交错树中的边(i,j)，它的A[i]+B[j]的值有所减小。也就说，它原来不属于相等子图，现在可能进入了相等子图，因而使相等子图得到了扩大。
　　现在的问题就是求d值了。为了使A[i]+B[j]>=w[i,j]始终成立，且至少有一条边进入相等子图，d应该等于min{A[i]+B[j]-w[i,j]|Xi在交错树中，Yi不在交错树中}。
　　以上就是KM算法的基本思路。但是朴素的实现方法，时间复杂度为O(n4)——需要找O(n)次增广路，每次增广最多需要修改O(n)次顶标，每次修改顶标时由于要枚举边来求d值，
复杂度为O(n2)。实际上KM算法的复杂度是可以做到O(n3)的。我们给每个Y顶点一个“松弛量”函数 slack，每次开始找增广路时初始化为无穷大。
在寻找增广路的过程中，检查边(i,j)时，如果它不在相等子图中，则让slack[j]变成原值与A [i]+B[j]-w[i,j]的较小值。
这样，在修改顶标时，取所有不在交错树中的Y顶点的slack值中的最小值作为d值即可。但还要注意一点：修改顶标后，要把所有的slack值都减去d。
最后把值加起来就ok。

  1 //poj2195 km算法
  2 #include<stdio.h>
  3 #include<string.h>
  4 #include<stdlib.h>
  5 #define INF 99999999
  6 char map[103][103];
  7 int n,m,g[103][103],visr[103],visl[103],pr[103],pl[103],match[103],num,slack[103];
  8 struct node
  9 {
 10     int x;
 11     int y;
 12 }mm[103],hh[103];
 13 int dfs(int u)
 14 {
 15     int i,j,val;
 16     visl[u]=1;
 17     for(i=0;i<num;i++)
 18     {
 19         if(!visr[i])
 20         {
 21             val=pl[u]+pr[i]-g[u][i];
 22             if(val==0)//即满足匹配条件
 23             {
 24                 visr[i]=1;
 25                 if(match[i]==-1||dfs(match[i]))
 26                 {
 27                     match[i]=u;
 28                     return 1;
 29                 }
 30             }
 31             if(val>0&&val<slack[i])
 32                 slack[i]=val;
 33         }
 34     }
 35     return 0;
 36 }
 37 int km()
 38 {
 39     int i,j,res,d;
 40     res=0;
 41     memset(pr,0,sizeof(pr));//右边的值为0
 42     for(i=0;i<num;i++)//左边的值为INF
 43         pl[i]=INF;
 44     memset(match,-1,sizeof(match));//匈牙利算法
 45     for(i=0;i<num;i++)
 46     {
 47 
 48         for(j=0;j<num;j++)//辅助数组slack[]初始为无穷大
 49             slack[j]=INF;
 50 
 51         while(1)//死循环知道有满足的匹配为止
 52         {
 53             memset(visl,0,sizeof(visl));
 54             memset(visr,0,sizeof(visr));
 55             if(dfs(i))break;//匹配成功结束
 56             d=INF;
 57             for(j=0;j<num;j++)
 58                 if(!visr[j]&&d>slack[j])//找到一个改进量dx，dx=min{Li+Lj-wi,j}(i∈S,j不∈T)        
 59                     d=slack[j];    
 60                                         //Li=Li-dx (i∈S)
 61             for(j=0;j<num;j++)            //Li=Li+dx (i∈T)
 62             {                            //重复以上过程，不断的调整L，直到求出完备匹配为止。
 63                 if(visl[j])
 64                     pl[j]-=d;
 65                 if(visr[j])
 66                     pr[j]+=d;
 67             }
 68         }
 69     }
 70     for(j=0;j<num;j++)
 71         res+=g[match[j]][j];
 72     return res;
 73 }
 74 int main()
 75 {
 76     int i,j;
 77     while(scanf("%d%d",&n,&m)!=EOF)
 78     {
 79         if(!n&&!m)break;
 80 
 81         for(i=0;i<n;i++)
 82             scanf("%s",map[i]);
 83 
 84         int numm,numh;
 85         numm=numh=0;
 86         for(i=0;i<n;i++)
 87         {
 88             for(j=0;j<m;j++)
 89             {
 90                 if(map[i][j]==‘m‘)//记录人所在的坐标
 91                 {
 92                     mm[numm].x=i;mm[numm++].y=j;
 93                 }
 94                 if(map[i][j]==‘H‘)//记录房子所在的坐标
 95                 {
 96                     hh[numh].x=i;hh[numh++].y=j;
 97                 }
 98             }
 99         }
100         //建图
101         num=numm;
102         for(i=0;i<numm;i++)
103         {
104             for(j=0;j<numm;j++)
105             {
106                 g[i][j]=-(abs(hh[i].x-mm[j].x)+abs(hh[i].y-mm[j].y));//任何人和任何房子间的距离
107             }
108         }
109         /*for(i=0;i<num;i++)
110         {
111             for(j=0;j<num;j++)
112                 printf("%d ",g[i][j]);
113             printf("\n");
114         }*/
115         int ans=km();
116         printf("%d\n",-ans);
117     }
118 }

 1 //hdu2255 竟然700+ms
 2 #include<stdio.h>
 3 #include<string.h>
 4 #define INF 99999999
 5 #define maxn 303
 6 int map[maxn][maxn],match[maxn],visl[maxn],visr[maxn],pr[maxn],pl[maxn],slack[maxn];
 7 int n,m;
 8 int dfs(int u)
 9 {
10     int i,j,val;
11     visl[u]=1;
12     for(i=1;i<=n;i++)
13     {
14         if(!visr[i])
15         {
16             val=pl[u]+pr[i]-map[u][i];
17             if(val==0)
18             {
19                 visr[i]=1;
20                 if(match[i]==-1||dfs(match[i]))
21                 {
22                     match[i]=u;
23                     return 1;
24                 }
25             }
26             if(val>0&&val<slack[i])
27                 slack[i]=val;
28         }
29     }
30     return 0;
31 }
32 int km()
33 {
34     int i,j,res,d;
35     res=0;
36     memset(pr,0,sizeof(pr));
37     for(i=1;i<=n;i++)
38         pl[i]=INF;
39     memset(match,-1,sizeof(match));
40     for(i=1;i<=n;i++)
41     {
42         for(j=1;j<=n;j++)
43             slack[j]=INF;
44         while(1)
45         {
46             memset(visl,0,sizeof(visl));
47             memset(visr,0,sizeof(visr));
48             if(dfs(i))break;
49             d=INF;
50             for(j=1;j<=n;j++)
51             {
52                 if(!visr[j]&&d>slack[j])
53                     d=slack[j];
54             }
55             for(j=1;j<=n;j++)
56             {
57                 if(visr[j])
58                     pr[j]+=d;
59                 if(visl[j])
60                     pl[j]-=d;
61             }
62         }
63     }
64     for(i=1;i<=n;i++)
65         res+=map[match[i]][i];
66     return res;
67 }
68 int main()
69 {
70     int i,j;
71     while(scanf("%d",&n)!=EOF)
72     {
73         for(i=0;i<=n;i++)
74             for(j=0;j<=n;j++)
75             {
76                 if(i==j)map[i][j]=0;
77                 else map[i][j]=INF;
78             }
79         for(i=1;i<=n;i++)
80         {
81             for(j=1;j<=n;j++)
82             {
83                 int z;
84                 scanf("%d",&z);
85                 map[i][j]=z;
86             }
87         }
88         int ans=km();
89         printf("%d\n",ans);
90     }
91 }