Gusfield演算法學習

luo_shen發表於2023-04-29

演算法詳解

等價流樹正如其名,樹上兩點間的路徑上的邊權最小值為圖上兩點間的最小割。

Gusfield演算法就是建等價流樹的一種演算法。設當前正在處理的集合為 \(S(|S|\ge 2)\),從 \(S\) 中任選兩個點 \(x,y\),求出 \(x,y\) 間的最小割也就是最大流 \(flow\),此時在最小割樹中加入一條從 \(x\)\(y\),邊權為 \(flow\)。在原圖中去掉所有的滿流邊後,\(S\) 會變成兩部分,一部分是在殘量網路上從起點能遍歷到的點,另一部分是遍歷不到的點,分治下去求解即可,很容易可以發現,最多求解 \(n-1\) 次,所以得到的必然是棵樹。

證明這個做法建出來的樹必然滿足樹上路徑邊權最小值等於原圖上兩點間的最小割。

有一種不嚴謹的證明,因為遞迴出來的兩邊都是合法的等價流樹,假設這兩棵等價流樹中的邊權全部大於 \((s,t)\) 這條邊,那麼顯然割掉 \((s,t)\) 代表的原圖中的邊集會優於割掉其它邊代表的原圖中的邊集。即若存在 \(u,v\) 分別再 \(s,t\) 側,\((s,t)\) 這條邊的邊權必然為 \(u,v\) 的最小割。至於更詳細的證明可以研讀2016年王文濤國集論文

例題

1.【模板】最小割樹(Gomory-Hu Tree)

其實我並不覺得這是GHT的模板,它只用到了大小,並沒有用到方案等其它資訊。直接建出等價流樹,每次詢問在等價流樹上找路徑邊權最小值即可,倍增預處理皆可。

點選檢視程式碼
#include<bits/stdc++.h>
#define ull unsigned long long
#define ll long long
#define pdi pair<double,int>
#define pii pair<int,int>
#define pb push_back
#define mp make_pair
#define eps 1e-9
using namespace std;
namespace IO{
    template<typename T>
    inline void read(T &x){
        x=0;
        int f=1;
        char ch=getchar();
        while(ch>'9'||ch<'0'){
            if(ch=='-'){
                f=-1;
            }
            ch=getchar();
        }
        while(ch>='0'&&ch<='9'){
            x=x*10+(ch-'0');
            ch=getchar();
        }
        x=(f==1?x:-x);
    }
    template<typename T>
    inline void write(T x){
        if(x<0){
            putchar('-');
            x=-x;
        }
        if(x>=10){
            write(x/10);
        }
        putchar(x%10+'0');
    }
    template<typename T>
    inline void write_endl(T x){
        write(x);
        putchar('\n');
    }
    template<typename T>
    inline void write_space(T x){
        write(x);
        putchar(' ');
    }
}
using namespace IO;
const int N=510,M=6010,Lg=11,inf=1e9;
int n,m,head[N],hd[N],tot1=0,tot=1,id[N];
struct edge{
    int v,w,nxt;
}e[M],G[M];
void add(int u,int v,int w){
    e[++tot].v=v;
    e[tot].nxt=head[u];
    e[tot].w=w;
    head[u]=tot;
}
void add1(int u,int v,int w){
    G[++tot1].v=v;
    G[tot1].w=w;
    G[tot1].nxt=hd[u];
    hd[u]=tot1;
}
namespace GHT{
    int dep[N],cur[N],tmp1[N],tmp2[N];
    bool bfs(int S,int T){
        for(int i=1;i<=n;i++){
            dep[i]=0;
        }
        queue<int>q;
        q.push(S);
        dep[S]=1;
        while(!q.empty()){
            int u=q.front();
            q.pop();
            for(int i=head[u];i;i=e[i].nxt){
                int v=e[i].v,w=e[i].w;
                if(w&&!dep[v]){
                    dep[v]=dep[u]+1;
                    q.push(v);
                    if(v==T){
                        return 1;
                    }
                }
            }
        }
        return 0;
    }
    int dfs(int u,int flow,int T){
        if(u==T){
            return flow;
        }
        int s=0;
        for(int i=cur[u];i;i=e[i].nxt){
            cur[u]=i;
            int v=e[i].v,w=e[i].w;
            if(w&&dep[v]==dep[u]+1){
                int res=dfs(v,min(flow,w),T);
                e[i].w-=res;
                e[i^1].w+=res;
                flow-=res;
                s+=res;
            }
            if(!flow){
                break;
            }
        }
        if(!s){
            dep[u]=0;
        }
        return s;
    }
    void init(){
        for(int i=0;i<=tot;i+=2){
            e[i].w+=e[i^1].w;
            e[i^1].w=0;
        }
    }
    int dinic(int S,int T){
        init();
        int ans=0;
        while(bfs(S,T)){
            memcpy(cur,head,sizeof(head));
            ans+=dfs(S,inf,T);
        }
        return ans;
    }
    void build(int l,int r){
        if(l>=r){
            return;
        }
        int x=id[l],y=id[l+1];
        int flow=dinic(x,y);
        add1(x,y,flow);
        add1(y,x,flow);
        int cnt1=0,cnt2=0;
        for(int i=l;i<=r;i++){
            if(dep[id[i]]){
                tmp1[++cnt1]=id[i];
            }
            else{
                tmp2[++cnt2]=id[i];
            }
        }
        for(int i=1;i<=cnt1;i++){
            id[l+i-1]=tmp1[i];
        }
        for(int i=1;i<=cnt2;i++){
            id[cnt1+i+l-1]=tmp2[i];
        }
        build(l,l+cnt1-1);
        build(l+cnt1,r);
    }
}
int f[N][Lg+5],mn[N][Lg+5],dep[N];
void make_tree(int u,int father){
    f[u][0]=father;
    dep[u]=dep[father]+1;
    for(int i=1;i<=Lg;i++){
        f[u][i]=f[f[u][i-1]][i-1];
        mn[u][i]=min(mn[u][i-1],mn[f[u][i-1]][i-1]);
    }
    for(int i=hd[u];i;i=G[i].nxt){
        int v=G[i].v,w=G[i].w;
        if(v==father){
            continue;
        }
        mn[v][0]=w;
        make_tree(v,u);
    }
}
int lca(int u,int v){
    if(dep[u]<dep[v]){
        swap(u,v);
    }
    int mini=inf;
    for(int i=Lg;i>=0;i--){
        if(dep[f[u][i]]>=dep[v]){
            mini=min(mini,mn[u][i]);
            u=f[u][i];
        }
    }
    if(u==v){
        return mini;
    }
    for(int i=Lg;i>=0;i--){
        if(f[u][i]!=f[v][i]){
            mini=min(mini,min(mn[u][i],mn[v][i]));
            u=f[u][i],v=f[v][i];
        }
    }
    return min(mini,min(mn[u][0],mn[v][0]));
}
signed main(){
    #ifndef ONLINE_JUDGE
        freopen("1.in","r",stdin);
        freopen("1.out","w",stdout);
    #endif
    memset(mn,0x3f,sizeof(mn));
    read(n),read(m);
    for(int i=1,u,v,w;i<=m;i++){
        read(u),read(v),read(w);
        add(u,v,w);
        add(v,u,0);
        add(v,u,w);
        add(u,v,0);
    }
    for(int i=1;i<=n;i++){
        id[i]=i;
    }
    GHT::build(1,n);
    make_tree(1,0);
    int q;
    read(q);
    while(q--){
        int u,v;
        read(u),read(v);
        write_endl(lca(u,v));
    }
    return 0;
}

2.[ZJOI2011]最小割

也是個板題,建出等價流樹,如果不連通,就每個連通塊之間建一條長度為 \(0\) 的邊,最後詢問時在樹上以每個點為起點掃一遍即可。

點選檢視程式碼
#include<bits/stdc++.h>
#define ull unsigned long long
#define int long long
#define pdi pair<double,int>
#define pii pair<int,int>
#define pb push_back
#define mp make_pair
#define eps 1e-9
using namespace std;
namespace IO{
    template<typename T>
    inline void read(T &x){
        x=0;
        int f=1;
        char ch=getchar();
        while(ch>'9'||ch<'0'){
            if(ch=='-'){
                f=-1;
            }
            ch=getchar();
        }
        while(ch>='0'&&ch<='9'){
            x=x*10+(ch-'0');
            ch=getchar();
        }
        x=(f==1?x:-x);
    }
    template<typename T>
    inline void write(T x){
        if(x<0){
            putchar('-');
            x=-x;
        }
        if(x>=10){
            write(x/10);
        }
        putchar(x%10+'0');
    }
    template<typename T>
    inline void write_endl(T x){
        write(x);
        putchar('\n');
    }
    template<typename T>
    inline void write_space(T x){
        write(x);
        putchar(' ');
    }
}
using namespace IO;
const int N=510,M=60010,inf=1e12;
int n,m,head[N],hd[N],tot1=0,tot=1,id[N],fa[N];
struct edge{
    int v,w,nxt;
}e[M],G[M];
int getfa(int x){
    if(fa[x]!=x){
        fa[x]=getfa(fa[x]);
    }
    return fa[x];
}
void merge(int u,int v){
    u=getfa(u),v=getfa(v);
    if(u!=v){
        fa[v]=u;
    }
}
void add(int u,int v,int w){
    e[++tot].v=v;
    e[tot].nxt=head[u];
    e[tot].w=w;
    head[u]=tot;
}
void add1(int u,int v,int w){
    G[++tot1].v=v;
    G[tot1].w=w;
    G[tot1].nxt=hd[u];
    hd[u]=tot1;
    merge(u,v);
}
namespace GHT{
    int dep[N],cur[N],tmp1[N],tmp2[N];
    bool bfs(int S,int T){
        for(int i=1;i<=n;i++){
            dep[i]=0;
        }
        queue<int>q;
        q.push(S);
        dep[S]=1;
        while(!q.empty()){
            int u=q.front();
            q.pop();
            for(int i=head[u];i;i=e[i].nxt){
                int v=e[i].v,w=e[i].w;
                if(w&&!dep[v]){
                    dep[v]=dep[u]+1;
                    q.push(v);
                    if(v==T){
                        return 1;
                    }
                }
            }
        }
        return 0;
    }
    int dfs(int u,int flow,int T){
        if(u==T){
            return flow;
        }
        int s=0;
        for(int i=cur[u];i;i=e[i].nxt){
            cur[u]=i;
            int v=e[i].v,w=e[i].w;
            if(w&&dep[v]==dep[u]+1){
                int res=dfs(v,min(flow,w),T);
                e[i].w-=res;
                e[i^1].w+=res;
                flow-=res;
                s+=res;
            }
            if(!flow){
                break;
            }
        }
        if(!s){
            dep[u]=0;
        }
        return s;
    }
    void init(){
        for(int i=0;i<=tot;i+=2){
            e[i].w+=e[i^1].w;
            e[i^1].w=0;
        }
    }
    int dinic(int S,int T){
        init();
        int ans=0;
        while(bfs(S,T)){
            memcpy(cur,head,sizeof(head));
            ans+=dfs(S,inf,T);
        }
        return ans;
    }
    void build(int l,int r){
        if(l>=r){
            return;
        }
        int x=id[l],y=id[l+1];
        int flow=dinic(x,y);
        add1(x,y,flow);
        add1(y,x,flow);
        int cnt1=0,cnt2=0;
        for(int i=l;i<=r;i++){
            if(dep[id[i]]){
                tmp1[++cnt1]=id[i];
            }
            else{
                tmp2[++cnt2]=id[i];
            }
        }
        for(int i=1;i<=cnt1;i++){
            id[l+i-1]=tmp1[i];
        }
        for(int i=1;i<=cnt2;i++){
            id[cnt1+i+l-1]=tmp2[i];
        }
        build(l,l+cnt1-1);
        build(l+cnt1,r);
    }
}
void clear(){
    for(int i=1;i<=n;i++){
        head[i]=hd[i]=0;
        id[i]=i;
        fa[i]=i;
    }
    tot=tot1=1;
}
int ans,x;
void dfs(int u,int father,int mn){
    if(mn<=x){
        ans++;
    }
    for(int i=hd[u];i;i=G[i].nxt){
        int v=G[i].v,w=G[i].w;
        if(v==father){
            continue;
        }
        dfs(v,u,min(mn,w));
    }
}
void solve(){
    read(n),read(m);
    clear();
    for(int i=1,u,v,w;i<=m;i++){
        read(u),read(v),read(w);
        add(u,v,w);
        add(v,u,0);
        add(v,u,w);
        add(u,v,0);
    }
    GHT::build(1,n);
    for(int i=2;i<=n;i++){
        if(getfa(i)!=getfa(1)){
            add1(i,1,0);
            merge(1,i);
        }
    }
    int q;
    read(q);
    while(q--){
        read(x);
        ans=0;
        for(int i=1;i<=n;i++){
            dfs(i,0,inf);
        }
        write_endl(ans/2);
    }

}
signed main(){
    #ifndef ONLINE_JUDGE
        freopen("1.in","r",stdin);
        freopen("1.out","w",stdout);
    #endif
    int t;
    read(t);
    while(t--){
        solve();
        puts("");
    }
    return 0;
}

3.[CQOI2016]不同的最小割

又一個板題,建出等價流樹,根據性質可得,兩點間的最小割為兩點在等價流樹上的路徑上的邊權最小值,所以統計有等價流樹多少種不同邊權即可。

點選檢視程式碼
#include<bits/stdc++.h>
#define ull unsigned long long
#define ll long long
#define pdi pair<double,int>
#define pii pair<int,int>
#define pb push_back
#define mp make_pair
#define eps 1e-9
using namespace std;
namespace IO{
    template<typename T>
    inline void read(T &x){
        x=0;
        int f=1;
        char ch=getchar();
        while(ch>'9'||ch<'0'){
            if(ch=='-'){
                f=-1;
            }
            ch=getchar();
        }
        while(ch>='0'&&ch<='9'){
            x=x*10+(ch-'0');
            ch=getchar();
        }
        x=(f==1?x:-x);
    }
    template<typename T>
    inline void write(T x){
        if(x<0){
            putchar('-');
            x=-x;
        }
        if(x>=10){
            write(x/10);
        }
        putchar(x%10+'0');
    }
    template<typename T>
    inline void write_endl(T x){
        write(x);
        putchar('\n');
    }
    template<typename T>
    inline void write_space(T x){
        write(x);
        putchar(' ');
    }
}
using namespace IO;
const int N=1010,M=600010,Lg=11,inf=1e9;
int n,m,head[N],hd[N],tot1=0,tot=1,id[N];
struct edge{
    int v,w,nxt;
}e[M];
void add(int u,int v,int w){
    e[++tot].v=v;
    e[tot].nxt=head[u];
    e[tot].w=w;
    head[u]=tot;
}
set<int>s;
namespace GHT{
    int dep[N],cur[N],tmp1[N],tmp2[N];
    bool bfs(int S,int T){
        for(int i=1;i<=n;i++){
            dep[i]=0;
        }
        queue<int>q;
        q.push(S);
        dep[S]=1;
        while(!q.empty()){
            int u=q.front();
            q.pop();
            for(int i=head[u];i;i=e[i].nxt){
                int v=e[i].v,w=e[i].w;
                if(w&&!dep[v]){
                    dep[v]=dep[u]+1;
                    q.push(v);
                    if(v==T){
                        return 1;
                    }
                }
            }
        }
        return 0;
    }
    int dfs(int u,int flow,int T){
        if(u==T){
            return flow;
        }
        int s=0;
        for(int i=cur[u];i;i=e[i].nxt){
            cur[u]=i;
            int v=e[i].v,w=e[i].w;
            if(w&&dep[v]==dep[u]+1){
                int res=dfs(v,min(flow,w),T);
                e[i].w-=res;
                e[i^1].w+=res;
                flow-=res;
                s+=res;
            }
            if(!flow){
                break;
            }
        }
        if(!s){
            dep[u]=0;
        }
        return s;
    }
    void init(){
        for(int i=0;i<=tot;i+=2){
            e[i].w+=e[i^1].w;
            e[i^1].w=0;
        }
    }
    int dinic(int S,int T){
        init();
        int ans=0;
        while(bfs(S,T)){
            memcpy(cur,head,sizeof(head));
            ans+=dfs(S,inf,T);
        }
        return ans;
    }
    void build(int l,int r){
        if(l>=r){
            return;
        }
        int x=id[l],y=id[l+1];
        int flow=dinic(x,y);
        s.insert(flow);
        int cnt1=0,cnt2=0;
        for(int i=l;i<=r;i++){
            if(dep[id[i]]){
                tmp1[++cnt1]=id[i];
            }
            else{
                tmp2[++cnt2]=id[i];
            }
        }
        for(int i=1;i<=cnt1;i++){
            id[l+i-1]=tmp1[i];
        }
        for(int i=1;i<=cnt2;i++){
            id[cnt1+i+l-1]=tmp2[i];
        }
        build(l,l+cnt1-1);
        build(l+cnt1,r);
    }
}
signed main(){
    #ifndef ONLINE_JUDGE
        freopen("1.in","r",stdin);
        freopen("1.out","w",stdout);
    #endif
    read(n),read(m);
    for(int i=1,u,v,w;i<=m;i++){
        read(u),read(v),read(w);
        add(u,v,w);
        add(v,u,0);
        add(v,u,w);
        add(u,v,0);
    }
    for(int i=1;i<=n;i++){
        id[i]=i;
    }
    GHT::build(1,n);
    write_endl(s.size());
    return 0;
}

4.cf343e

先想最暴力的做法,是不是列舉排列,然後跑最小割,求最小割的和的最大值。然後我們發現這樣要把任意兩個點之間的最小割都求出來,於是建出等價流樹。
此時我們再看一下一個排列表示什麼,表示的是給每個點標一個dfs序,然後按照dfs序在最小割樹上搜一遍,因為等價流樹的性質,所以最小割的和等於dfs序相鄰的兩點的路徑上邊權最小值的和。所以問題就轉化為了求一個dfs序,使得最小割樹上邊權小的邊儘量少走。

先討論最小邊的貢獻的次數,只要經過最小邊則必然造成一次貢獻。因為要遍歷整棵樹,所以最小邊至少經過一次,即最少造成一次貢獻。
那麼有可能造成更多的貢獻嗎?這是不可能的。
令最小邊的兩端點為 \(u,v\),若斷掉最小邊,則一棵樹會變成兩顆樹,分別令 \(u,v\) 為根,那麼這兩棵樹分別為樹 \(u\) 和樹 \(v\)。一個非常顯然的事是,我們可以先遍歷完樹 \(u\),再遍歷樹 \(v\),這樣最小邊就只會造成 \(1\) 次貢獻,同理每條邊都會造成 \(1\) 次貢獻,答案的最大值就為等價流樹上邊權之和,方案可以 \(n^2\) 掃一遍。

點選檢視程式碼
#include<bits/stdc++.h>
#define ull unsigned long long
#define ll long long
#define pii pair<int,int>
#define pb push_back
#define mp make_pair
using namespace std;
namespace IO{
    template<typename T>
    inline void read(T &x){
        x=0;
        int f=1;
        char ch=getchar();
        while(ch>'9'||ch<'0'){
            if(ch=='-'){
                f=-1;
    }
            ch=getchar();
        }
        while(ch>='0'&&ch<='9'){
            x=x*10+(ch-'0');
            ch=getchar();
        }
        x=(f==1?x:-x);
    }
    template<typename T>
    inline void write(T x){
        if(x<0){
            putchar('-');
            x=-x;
        }
        if(x>=10){
            write(x/10);
        }
        putchar(x%10+'0');
    }
    template<typename T>
    inline void write_endl(T x){
        write(x);
        putchar('\n');
    }
    template<typename T>
    inline void write_space(T x){
        write(x);
        putchar(' ');
    }
}
using namespace IO;
const int N=300,M=1e4,inf=1e9;
int n,m,tot=1,tot1=1,head[N],hd[N],id[N],ans;
struct edge{
    int u,v,w,nxt;
}e[M],G[M];
void add(int u,int v,int w){
    e[++tot].v=v;
    e[tot].w=w;
    e[tot].nxt=head[u];
    head[u]=tot;
}
void add1(int u,int v,int w){
    G[++tot1].v=v;
    G[tot1].u=u;
    G[tot1].w=w;
    G[tot1].nxt=hd[u];
    ans+=w;
    hd[u]=tot1;
}
namespace Gusfield{
    int dep[N],cur[N],tmp1[N],tmp2[N];
    bool bfs(int S,int T){
        for(int i=1;i<=n;i++){
            dep[i]=0;
        }
        queue<int>q;
        q.push(S);
        dep[S]=1;
        while(!q.empty()){
            int u=q.front();
            q.pop();
            for(int i=head[u];i;i=e[i].nxt){
                int v=e[i].v,w=e[i].w;
                if(w&&!dep[v]){
                    dep[v]=dep[u]+1;
                    q.push(v);
                    if(v==T){
                        return 1;
                    }
                }
            }
        }
        return 0;
    }
    int dfs(int u,int flow,int T){
        if(u==T){
            return flow;
        }
        int s=0;
        for(int i=cur[u];i;i=e[i].nxt){
            cur[u]=i;
            int v=e[i].v,w=e[i].w;
            if(w&&dep[v]==dep[u]+1){
                int res=dfs(v,min(flow,w),T);
                e[i].w-=res;
                e[i^1].w+=res;
                flow-=res;
                s+=res;
            }
            if(!flow){
                break;
            }
        }
        if(!s){
            dep[u]=0;
        }
        return s;
    }
    void init(){
        for(int i=0;i<=tot;i+=2){
            e[i].w+=e[i^1].w;
            e[i^1].w=0;
        }
    }
    int dinic(int S,int T){
        init();
        int ans=0;
        while(bfs(S,T)){
            memcpy(cur,head,sizeof(head));
            ans+=dfs(S,inf,T);
        }
        return ans;
    }
    void build(int l,int r){
        if(l>=r){
            return;
        }
        int x=id[l],y=id[l+1];
        int flow=dinic(x,y);
        add1(x,y,flow);
        add1(y,x,flow);
        int cnt1=0,cnt2=0;
        for(int i=l;i<=r;i++){
            if(dep[id[i]]){
                tmp1[++cnt1]=id[i];
            }
            else{
                tmp2[++cnt2]=id[i];
            }
        }
        for(int i=1;i<=cnt1;i++){
            id[l+i-1]=tmp1[i];
        }
        for(int i=1;i<=cnt2;i++){
            id[cnt1+i+l-1]=tmp2[i];
        }
        build(l,l+cnt1-1);
        build(l+cnt1,r);
    }
}
bool vis[M];
vector<int>s;
void get(int u,int fa){
    s.pb(u);
    for(int i=hd[u];i;i=G[i].nxt){
        int v=G[i].v;
        if(vis[i]||v==fa){
            continue;
        }
        get(v,u);
    }
}
void dfs(int u){
    if(s.size()==1){
        write_space(s[0]);
        return;
    }
    int id=0,mn=inf;
    for(auto x:s){
        for(int i=hd[x];i;i=G[i].nxt){
            if(vis[i]){
                continue;
            }
            int w=G[i].w;
            if(w<mn){
                mn=w;
                id=i;
            }
        }
    }
    cerr<<id/2<<endl;
    vis[id]=vis[id^1]=1;
    vector<int>().swap(s);
    get(G[id].u,0);
    dfs(G[id].u);
    vector<int>().swap(s);
    get(G[id].v,0);
    dfs(G[id].v);
}
void solve(){
    read(n),read(m);
    for(int i=1,u,v,w;i<=m;i++){
        read(u),read(v),read(w);
        add(u,v,w);
        add(v,u,0);
        add(v,u,w);
        add(u,v,0);
    }
    for(int i=1;i<=n;i++){
        id[i]=i;
    }
    Gusfield::build(1,n);
    write_endl(ans>>1);
    for(int i=1;i<=n;i++){
        s.pb(i);
    }
    dfs(1);
}
signed main(){
    #ifndef ONLINE_JUDGE
        freopen("1.in","r",stdin);
        freopen("1.out","w",stdout);
    #endif
    int t=1;
    while(t--){
        solve();
    }
    return 0;
}

相關文章