帶花樹算法
先放上大神的blog，個人認為沒辦法比這位dalao解釋的更清楚。
在北京冬令營的時候，yby提到了“帶花樹開花”算法來解非二分圖的最大匹配。
於是，我打算看看這是個什么玩意。其實之前，我已經對這個算法了解了個大概，但是。。。真的不敢去寫。
有一個叫Galil Zvi的人（應該叫計算機科學家），寫了篇論文：
Efficient Algorithms for Finding Maximal Matching in Graphs
（如果你在網上搜不到，可以：
http://builtinclz.abcz8.com/art/2012/Galil%20Zvi.pdf）
這篇論文真神啊，它解決了4個問題：
（一般圖+二分圖）的（最大匹配+最大權匹配）問題。
算法的思想、故事，請自己看論文吧。
這個論文告訴了我們很多有趣的東西，例如：
用Dinic實現的二分圖匹配的時間復雜度其實是O(M*N^0.5)，這也許能夠解釋為什么一般網絡流算法比Hungry要快了。
另外，帶花樹算法的正確性的證明比較困難；而其時間復雜度是可以做到O(M*N^0.5)的，不過要詳細實現，那么就快能到“ACM最長論文獎”了。
我寫了一個實例代碼：
http://builtinclz.abcz8.com/art/2012/ural1099.cpp
沒錯，這是用來解決URAL 1099 Work Schedule那題的。時間復雜度是O（N＾３）
簡述一下“帶花樹”算法吧：
它的核心思想還是找增廣路。假設已經匹配好了一堆點，我們從一個沒有匹配的節點s開始，使用BFS生成搜索樹。每當發現一個節點u，如果u還沒有被匹配，那么就可以進行一次成功的增廣；否則，我們就把節點u和它的配偶v一同接到樹上，之后把ｖ丟進隊列繼續搜索。我們給每個在搜索樹上的點一個類型：S或者T。當ｕ把它的配偶ｖ扔進隊列的時候，我們把ｕ標記為T型，ｖ標記為S型。於是，搜索樹的樣子是這樣的：
ｓ
／
＼
ａ
ｂ
｜ ｜
ｃ ｄ
／
＼
／
＼
ｅ
ｆ
ｕ ｊ
｜ ｜ ｜ ｜
ｉ ｊ ｖ ｋ
其中，黑色豎線相連的兩個點是已經匹配好的，藍色斜線表示兩個點之間有邊，但是沒有配對。T型的用紅色，S型的用黑色。
這里有個小問題：一個S型點ｄ在某一步擴展的時候發現了點ｕ，如果ｕ已經在搜索樹上了（即，出現了環），怎么辦？
我們規定，如果ｕ的類型是T型，就無視這次發現；（這意味着我們找到了一個長度為偶數的環，直接無視）
ｓ
／
＼
ａ
ｂ
｜ ｜
ｃ ｄ 如果連出來的邊是指向T型點的，就無視這個邊。
／
＼
／
＼
ｅ
ｆ
＜
－
ｇ
｜ ｜ ｜
ｉ ｊ ｋ
否則，我們找到了一個長度為奇數的環，就要進行一次“縮花”的操作！所謂縮花操作，就是把這個環縮成一個點。
ｓ
／
＼
ａ
ｂ
｜ ｜
ｃ ｄ
／
＼
／
＼
ｅ
ｆ
｜
ｇ
｜ ｜
｜ ｜
ｉ ｕ
＜－＋ ｋ
這個圖縮花之后變成了５個點（一個大點，或者叫一朵花，加原來的４個點）：
縮點完成之后，還要把原來環里面的T型點統統變成S型點，之后扔到隊列里去。
＋－－－－－－－－－－－－－＋
｜ ｜
｜ ｓ ｜
｜
／
＼ ｜
｜
ａ
ｂ
｜
｜ ｜ ｜ ｜ 現在是一個點了！還是一個S點。
｜ ｃ ｄ ｜
｜
／
＼
／ ＼ ｜
＋－｜－－
ｆ－－ｕ
－－－｜－－－＋
｜ ｜
｜ ｜
｜ ｜ ｜ ｜
｜ ｜ ｜ ｜
｜ ＋－－－－－－－－－－－－－＋ ｜
｜ ｜
ｅ ｇ
｜ ｜
ｉ ｋ
為什么能縮成一個點呢？我們看一個長度為奇數的環（例如上圖中的ｓ－ｂ－ｄ－ｕ－ｆ－ｃ－ａ－），如果我們能夠給它中的任意一個點找一個出度（配偶），那么環中的其他點正好可以配成對，這說明，每個點的出度都是等效的。例如，假設我們能夠給圖中的點ｄ另找一個配偶（例如ｄ＇好了），那么，環中的剩下６個點正好能配成３對，一個不多，一個不少（算上ｄ和ｄ＇一共４對剛剛好）。
ａ－ｓ－ｂ－ｄ
－ｄ＇ ａ ｓ－ｂ ｄ－ｄ＇
＼ ｜ ＝＞ ＼
ｃ
－ｆ－ｕ ｃ ｆ－ｕ
這就是我們縮點的思想來源。有一個勞苦功高的計算機科學家證明了：縮點之前和縮點之后的圖是否有增廣路的情況是相同的。
縮起來的點又叫一朵花（ｂｌｏｓｓｏｍ）．
注意到，組成一朵花的里面可能嵌套着更小的花。
當我們最終找到一條增廣路的時候，要把嵌套着的花層層展開，還原出原圖上的增廣路出來。
嗯，現在你對實現這個算法有任何想法嗎？
天啊，還要縮點……寫死誰。。。。。。
我一開始也是這么想的。
我看了一眼網上某個大牛的程序，之后結合自己的想法，很努力地寫出了一個能AC的版本。
實現的要點有什么呢？
首先，我們不“顯式”地表示花。我們記錄一個Nｅｘｔ數組，表示最終增廣的時候的路徑上的后繼。同時，我們維護一個並查集，表示每個點現在在以哪個點為根的花里（一個花被縮進另一朵花之后就不算花了）。還要記錄每個點的標記。
主程序是一段BFS。對於每個由ｘ發展出來的點ｙ，分４種情況討論：
１。ｘｙ是配偶（不說夫妻，這是非二分圖。。。）或者ｘｙ現在是一個點（在一朵花里）：直接無視
２。ｙ是T型點：直接無視
３。ｙ目前單身：太好了，進行增廣！
４。ｙ是一個S型點：縮點！縮點！
縮點的時候要進行的工作：
１。找ｘ和ｙ的LCA（的根）ｐ。找LCA可以用各種方法。。。直接朴素也行。
２。在Nｅｘｔ數組中把ｘ和ｙ接起來（表示它們形成環了！）
３。從ｘ、ｙ分別走到ｐ，修改並查集使得它們都變成一家人，同時沿路把Nｅｘｔ數組接起來。
Nｅｘｔ數組很奇妙。每時每刻，它實際形成了若干個掛在一起的雙向鏈表來表示一朵花內部的走法。
－－－－
／ ＼
＜－－＋
｜ ｜ ｜
｜ ｜
＜－－＋
ｖ
ｖ
－－－－－－－－－－
／ ＼
＋－ －－＋
｜ ｜
｜ ｜
＋－－－－
＞ｓ
＜－－－－－－＋
有權圖的最大匹配怎么做？
看論文吧。。。用類似KM的方法，不過，是給每個花再來一個權值。真的很復雜。。。
有一個人寫了代碼，好像是GPL許可證。。。你最好想辦法搜到它的網站來看看版權的問題；總之，我先貼出來：
一開始總覺得這東西么。。。沒有什么用。。。因為到了考場上肯定有很多特殊的性質（比如沒有奇環），或者數據規模很小（暴力枚舉所有點貪心）。這樣可以用更少的時間得到一個比較讓人滿意的分數。但是就有這么倆集訓隊爺非得在論文里討論一般圖的最大匹配。。。這就使省選出這類題目的可能性無限變大了。考慮到今年出了圓方樹，所以還是去學習了一下。。。
筆者認為理解了上面dalao的講解之后，對於帶花樹算法的實質內容以及實現形式都應該有了一個比較深的認識。但是代碼確實不是很好打以至於筆者也是各種搜羅才寫出了一份看上去比較優美的。然后加上了一些筆者自己的見解，也就是注釋，希望能夠幫到各位理解代碼吧。。
#include <iostream> #include <cstdlib> #include <cmath> #include <string> #include <cstring> #include <cstdio> #include <algorithm> #include <queue> #include <set> #include <map> #define re register #define max(a,b) ((a)>(b)?(a):(b)) #define min(a,b) ((a)<(b)?(a):(b)) #define MAXN 5001 #define MAXM 1000001 using namespace std; typedef long long ll; typedef unsigned long long ull; #define ms(arr) memset(arr, 0, sizeof(arr)) const int inf = 0x3f3f3f3f; int f[MAXN],head[MAXN],next[MAXN],num,match[MAXN],t,vis[MAXN],mark[MAXN],n; //f是並查集，match是匹配數組，mark表示S/T點。 int x[MAXN],y[MAXN],L; struct po { int from,nxt,to; }edge[MAXM]; queue<int> q; inline int read() { int x=0,c=1; char ch=' '; while((ch>'9'||ch<'0')&&ch!='-')ch=getchar(); while(ch=='-') c*=-1,ch=getchar(); while(ch<='9'&&ch>='0')x=x*10+ch-'0',ch=getchar(); return x*c; } int find(int x) { return f[x]==x?f[x]:f[x]=find(f[x]); } inline void add_edge(int from,int to) { edge[++num].nxt=head[from]; edge[num].to=to; head[from]=num; } inline void add(int from,int to) { add_edge(from,to); add_edge(to,from); } inline int lca(int x,int y)//朴素找lca，雖然筆者自己是死記硬背的。。貌似只會樹剖找lca。。 { t++; while(1){ if(x){ x=find(x);//點要對應到相應的花上。 if(vis[x]==t) return x; vis[x]=t; if(match[x]) x=next[match[x]]; else x=0; } swap(x,y); } } inline void flower(int a,int p) { while(a!=p){ int b=match[a],c=next[b]; if(find(c)!=p) next[c]=b; //next數組是用來標記花中的路徑，結合match就實現了雙向鏈表的功能。 //我們可以知道一個奇環里的所有點只要向外面連邊就有可能與外面匹配，所以所有奇環中的點都
要變成S型點扔到隊列中。又因為環內部肯定是匹配飽和的，所以這一堆點最多只能匹配出來一
個點，所以每次匹配到一個點就可以跳出循環了。 if(mark[b]==2) q.push(b),mark[b]=1; if(mark[c]==2) q.push(c),mark[c]=1; f[find(a)]=find(b);f[find(b)]=find(c); a=c; } } inline void work(int s) { for(re int i=1;i<=n;i++) next[i]=mark[i]=vis[i]=0,f[i]=i;//每個階段都需要重新標記 while(!q.empty()) q.pop(); mark[s]=1;q.push(s); while(!q.empty()){ if(match[s]) return; int x=q.front();q.pop();//隊列中的點均為S型 for(re int i=head[x];i;i=edge[i].nxt){ int y=edge[i].to; if(match[x]==y) continue;//x與y是配偶，忽略掉 if(find(x)==find(y)) continue;//x與y在一朵花里，忽略掉 if(mark[y]==2) continue;//y是T形點，忽略掉。 if(mark[y]==1) {//y是s型點，需要縮點 int r=lca(x,y); if(find(x)!=r) next[x]=y;//x和r不在同一個花朵， next標記花朵內的路徑 if(find(y)!=r) next[y]=x;//同上 flower(x,r);flower(y,r);//將r--x--y--r縮成一個點 } else if(!match[y]){//y沒有匹配過，可以增廣 next[y]=x; for(re int u=y;u;){//增廣操作 int v=next[u],w=match[v]; match[v]=u;match[u]=v;u=w; } break; } else {//y點已經匹配過，將其配偶作為S型點拉進來，y點本身為T型點，繼續搜索。 next[y]=x; mark[match[y]]=1; q.push(match[y]); mark[y]=2; } } } } inline void out(){ for(re int i=1;i<=n;i++) if(!match[i]){ printf("NO\n"); return; } printf("YES\n"); } int main() { //freopen("date.in","r",stdin); while(cin>>n){ memset(match,0,sizeof(match)); memset(head,0,sizeof(head));num=0;t=0; for(re int i=1;i<=n;i++) x[i]=read(),y[i]=read(); L=read(); for(re int i=1;i<=n;i++) for(re int j=i+1;j<=n;j++){ if(abs(x[i]-x[j])+abs(y[i]-y[j])<=L){ add(i,j); } } for(re int i=1;i<=n;i++) if(!match[i]) work(i); out(); } return 0; }
然后還有一個就是帶權帶花樹的代碼以及實現。。筆者找到了一份比價簡單明白的代碼，雖然自己還沒有研究，不過給大家放上來，希望對大家有所幫助。
#include<bits/stdc++.h> using namespace std; //from vfleaking //自己進行一些進行一些小修改 #define INF INT_MAX #define MAXN 400 struct edge{ int u,v,w; edge(){} edge(int u,int v,int w):u(u),v(v),w(w){} }; int n,n_x; edge g[MAXN*2+1][MAXN*2+1]; int lab[MAXN*2+1]; int match[MAXN*2+1],slack[MAXN*2+1],st[MAXN*2+1],pa[MAXN*2+1]; int flower_from[MAXN*2+1][MAXN+1],S[MAXN*2+1],vis[MAXN*2+1]; vector<int> flower[MAXN*2+1]; queue<int> q; inline int e_delta(const edge &e){ // does not work inside blossoms return lab[e.u]+lab[e.v]-g[e.u][e.v].w*2; } inline void update_slack(int u,int x){ if(!slack[x]||e_delta(g[u][x])<e_delta(g[slack[x]][x]))slack[x]=u; } inline void set_slack(int x){ slack[x]=0; for(int u=1;u<=n;++u) if(g[u][x].w>0&&st[u]!=x&&S[st[u]]==0)update_slack(u,x); } void q_push(int x){ if(x<=n)q.push(x); else for(size_t i=0;i<flower[x].size();i++)q_push(flower[x][i]); } inline void set_st(int x,int b){ st[x]=b; if(x>n)for(size_t i=0;i<flower[x].size();++i) set_st(flower[x][i],b); } inline int get_pr(int b,int xr){ int pr=find(flower[b].begin(),flower[b].end(),xr)-flower[b].begin(); if(pr%2==1){//檢查他在前一層圖是奇點還是偶點 reverse(flower[b].begin()+1,flower[b].end()); return (int)flower[b].size()-pr; }else return pr; } inline void set_match(int u,int v){ match[u]=g[u][v].v; if(u>n){ edge e=g[u][v]; int xr=flower_from[u][e.u],pr=get_pr(u,xr); for(int i=0;i<pr;++i)set_match(flower[u][i],flower[u][i^1]); set_match(xr,v); rotate(flower[u].begin(),flower[u].begin()+pr,flower[u].end()); } } inline void augment(int u,int v){ for(;;){ int xnv=st[match[u]]; set_match(u,v); if(!xnv)return; set_match(xnv,st[pa[xnv]]); u=st[pa[xnv]],v=xnv; } } inline int get_lca(int u,int v){ static int t=0; for(++t;u||v;swap(u,v)){ if(u==0)continue; if(vis[u]==t)return u; vis[u]=t;//這種方法可以不用清空v陣列 u=st[match[u]]; if(u)u=st[pa[u]]; } return 0; } inline void add_blossom(int u,int lca,int v){ int b=n+1; while(b<=n_x&&st[b])++b; if(b>n_x)++n_x; lab[b]=0,S[b]=0; match[b]=match[lca]; flower[b].clear(); flower[b].push_back(lca); for(int x=u,y;x!=lca;x=st[pa[y]]) flower[b].push_back(x),flower[b].push_back(y=st[match[x]]),q_push(y); reverse(flower[b].begin()+1,flower[b].end()); for(int x=v,y;x!=lca;x=st[pa[y]]) flower[b].push_back(x),flower[b].push_back(y=st[match[x]]),q_push(y); set_st(b,b); for(int x=1;x<=n_x;++x)g[b][x].w=g[x][b].w=0; for(int x=1;x<=n;++x)flower_from[b][x]=0; for(size_t i=0;i<flower[b].size();++i){ int xs=flower[b][i]; for(int x=1;x<=n_x;++x) if(g[b][x].w==0||e_delta(g[xs][x])<e_delta(g[b][x])) g[b][x]=g[xs][x],g[x][b]=g[x][xs]; for(int x=1;x<=n;++x) if(flower_from[xs][x])flower_from[b][x]=xs; } set_slack(b); } inline void expand_blossom(int b){ // S[b] == 1 for(size_t i=0;i<flower[b].size();++i) set_st(flower[b][i],flower[b][i]); int xr=flower_from[b][g[b][pa[b]].u],pr=get_pr(b,xr); for(int i=0;i<pr;i+=2){ int xs=flower[b][i],xns=flower[b][i+1]; pa[xs]=g[xns][xs].u; S[xs]=1,S[xns]=0; slack[xs]=0,set_slack(xns); q_push(xns); } S[xr]=1,pa[xr]=pa[b]; for(size_t i=pr+1;i<flower[b].size();++i){ int xs=flower[b][i]; S[xs]=-1,set_slack(xs); } st[b]=0; } inline bool on_found_edge(const edge &e){ int u=st[e.u],v=st[e.v]; if(S[v]==-1){ pa[v]=e.u,S[v]=1; int nu=st[match[v]]; slack[v]=slack[nu]=0; S[nu]=0,q_push(nu); }else if(S[v]==0){ int lca=get_lca(u,v); if(!lca)return augment(u,v),augment(v,u),true; else add_blossom(u,lca,v); } return false; } inline bool matching(){ memset(S+1,-1,sizeof(int)*n_x); memset(slack+1,0,sizeof(int)*n_x); q=queue<int>(); for(int x=1;x<=n_x;++x) if(st[x]==x&&!match[x])pa[x]=0,S[x]=0,q_push(x); if(q.empty())return false; for(;;){ while(q.size()){ int u=q.front();q.pop(); if(S[st[u]]==1)continue; for(int v=1;v<=n;++v) if(g[u][v].w>0&&st[u]!=st[v]){ if(e_delta(g[u][v])==0){ if(on_found_edge(g[u][v]))return true; }else update_slack(u,st[v]); } } int d=INF; for(int b=n+1;b<=n_x;++b) if(st[b]==b&&S[b]==1)d=min(d,lab[b]/2); for(int x=1;x<=n_x;++x) if(st[x]==x&&slack[x]){ if(S[x]==-1)d=min(d,e_delta(g[slack[x]][x])); else if(S[x]==0)d=min(d,e_delta(g[slack[x]][x])/2); } for(int u=1;u<=n;++u){ if(S[st[u]]==0){ if(lab[u]<=d)return 0; lab[u]-=d; }else if(S[st[u]]==1)lab[u]+=d; } for(int b=n+1;b<=n_x;++b) if(st[b]==b){ if(S[st[b]]==0)lab[b]+=d*2; else if(S[st[b]]==1)lab[b]-=d*2; } q=queue<int>(); for(int x=1;x<=n_x;++x) if(st[x]==x&&slack[x]&&st[slack[x]]!=x&&e_delta(g[slack[x]][x])==0) if(on_found_edge(g[slack[x]][x]))return true; for(int b=n+1;b<=n_x;++b) if(st[b]==b&&S[b]==1&&lab[b]==0)expand_blossom(b); } return false; } inline pair<long long,int> weight_blossom(){ memset(match+1,0,sizeof(int)*n); n_x=n; int n_matches=0; long long tot_weight=0; for(int u=0;u<=n;++u)st[u]=u,flower[u].clear(); int w_max=0; for(int u=1;u<=n;++u) for(int v=1;v<=n;++v){ flower_from[u][v]=(u==v?u:0); w_max=max(w_max,g[u][v].w); } for(int u=1;u<=n;++u)lab[u]=w_max; while(matching())++n_matches; for(int u=1;u<=n;++u) if(match[u]&&match[u]<u) tot_weight+=g[u][match[u]].w; return make_pair(tot_weight,n_matches); } inline void init_weight_graph(){ for(int u=1;u<=n;++u) for(int v=1;v<=n;++v) g[u][v]=edge(u,v,0); } int main(){ int m; scanf("%d%d",&n,&m); init_weight_graph(); for(int i=0;i<m;++i){ int u,v,w; scanf("%d%d%d",&u,&v,&w); g[u][v].w=g[v][u].w=w; } printf("%lld\n",weight_blossom().first); for(int u=1;u<=n;++u)printf("%d ",match[u]);puts(""); return 0; }