二分图最佳完美匹配

最佳完美匹配即权值和最大的完美匹配，这里运用到KM算法。
引进两个概念：
1、可行顶标，为一个节点函数l，使任意边（u，v）均满足l(u)+l(v)>=w(u,v);
2、相等子图,为原图的一个生成子图，包含所有点以及满足l(u)+l(v)=w(u,v)的边。
定理：如相等子图G’有完美匹配，则该匹配即为所求。
证明：G’的权值和为所有顶点的顶标和Sum，而此时原图G任意一个完美匹配的权值和均<=Sum。
算法流程：
1.构造一个可行顶标（一般可令lx（u）为与u相连的边的最大权值，ly（v）=0）；
2.搜索寻找增广路，如找到换下一个节点，否则转3；
3.设X为搜索中左边已访问的点，Y为右边已访问的点，令路M(搜索失败)上所有边，在X中的顶标减去d，在Y中的顶标加上d，
d为min{l(x)+l(y)-w(x,y)},x属于X,y不属于Y。下面说明一下为什么这样选取d：
如果比d小不能保证有新边加入，而大于d会顶标变得不可行，破坏性质。

4.调整后l继续执行2；

时间O(n^4)如果加一个小优化，在计算d时用一个数组slack[v]（v不属于Y）来维护与v相连的边min{lx(u)+ly(v)-w(u,v)}可降下一个n。
来道板子POJ3565

为什么上面的A了，下面的WA？（窝调了一天QAQ）

#include <cstdio>
#include <algorithm>
#include <cstring>
#define maxn 1005
#define eps 1e-10
#define INF 0x3f3f3f3f
#include <cmath>
int n,m,match[maxn],head[maxn],vis[maxn],visy[maxn];
using namespace std;
double d,slack[maxn],lx[maxn],ly[maxn],x[maxn],y[maxn],p,q;
struct xx{
    int v,next;
    double q;
}b[maxn*maxn*2];

void add(int u,int v,double q)
{
    b[++m]=(xx){v,head[u],q};
    head[u]=m;
}

int find(int t)
{
    vis[t]=1;
    for (int k=head[t];k!=0;k=b[k].next)
        {
            int v=b[k].v;
            double g=lx[t]+ly[v]-b[k].q;
            if (visy[v]) continue;
            if (abs(g)<=eps)//该边可行 
            {
                visy[v]=1;
                if (!match[v]||find(match[v]))
                    {
                        match[v]=t;
                        return true;
                    }
            }else slack[v]=min(slack[v],g);//t在X中且v不在Y中 
        }
    return false;
}

void km()
{
    memset(match,0,sizeof(match));
    //memset(ly,0,sizeof(ly));
    for (int i=1;i<=n;i++) ly[i]=0;
    for (int i=1;i<=n;i++)
    {
        lx[i]=-INF;
        for (int k=head[i];k!=0;k=b[k].next) 
          lx[i]=max(lx[i],b[k].q);
    }
    for (int i=1;i<=n;i++)
      {
        //memset(slack,INF,sizeof(slack));
        for (int j=1;j<=n;j++) slack[j]=INF;
        for (;;)
         {
            memset(vis,0,sizeof(vis));
            memset(visy,0,sizeof(visy));
            if (find(i)) break;//找到增广路 
            d=INF;             //调整lx和ly 
            for (int j=1;j<=n;j++)
              if (!visy[j]) d=min(d,slack[j]);
            for (int j=1;j<=n;j++){
              if (vis[j]) lx[j]-=d;//在X中 
            }
            for (int j=1;j<=n;j++){
              if (visy[j]) ly[j]+=d;//在Y中 
              else slack[j]-=d;     //不在Y中则维护slack 
            }                  
       }

      }
      return ;
}

int main()
{

    while (scanf("%d",&n)==1)
    {
        m=0;
        memset(head,0,sizeof(head));
        for (int i=1;i<=n;i++) scanf("%lf%lf",&x[i],&y[i]);
        for (int i=1;i<=n;i++){
             scanf("%lf%lf",&p,&q);
             for (int j=1;j<=n;j++){
             double dis=sqrt((x[j]-p)*(x[j]-p)+(y[j]-q)*(y[j]-q));//?????? 
             add(j,i,-dis);//add(i,j,-dis);
            }
        }
        km();
        for (int i=1;i<=n;i++) 
          for (int j=1;j<=n;j++)
            if (match[j]==i)
            {
                printf("%d\n",j);
                break;
            }
    }
    return 0;
}