[模板] 最大流和费用流分别的两种做法

注意:以下所有说明均以帮助理解模板为目的,不保证正确性。

最大流

dinic

考虑每次找一条S到T的不满流的路径并进行增广,但需要解决转圈圈的问题

所以首先用bfs给你的网络分层,之后再按照层(每个点只能走到下一层的点)来dfs,尽可能地把流量都占满

当前弧优化:在一次dfs中,每个点已经被访问完的出边再被访问没有意义,所以可以记录下来已经做到了哪个边

LOJ101
 1 #include<bits/stdc++.h>
 2 #define pa pair<int,int>
 3 #define CLR(a,x) memset(a,x,sizeof(a))
 4 #define MP make_pair
 5 #define fi first
 6 #define se second
 7 using namespace std;
 8 typedef long long ll;
 9 typedef unsigned long long ull;
10 typedef unsigned int ui;
11 typedef long double ld;
12 const int maxn=105,maxm=5005;
13 
14 inline char gc(){
15     return getchar();
16     static const int maxs=1<<16;static char buf[maxs],*p1=buf,*p2=buf;
17     return p1==p2&&(p2=(p1=buf)+fread(buf,1,maxs,stdin),p1==p2)?EOF:*p1++;
18 }
19 inline ll rd(){
20     ll x=0;char c=gc();bool neg=0;
21     while(c<'0'||c>'9'){if(c=='-') neg=1;c=gc();}
22     while(c>='0'&&c<='9') x=(x<<1)+(x<<3)+c-'0',c=gc();
23     return neg?(~x+1):x;
24 }
25 
26 struct Edge{
27     int b,l,ne;
28 }eg[maxm*2];
29 int egh[maxn],ect=1,cur[maxn];
30 int N,M,S,T;
31 int dep[maxn];
32 
33 inline void adeg(int a,int b,int l){
34     eg[++ect].b=b,eg[ect].l=l,eg[ect].ne=egh[a],egh[a]=ect;
35 }
36 
37 inline bool bfs(){
38     static queue<int> q;
39     for(int i=1;i<=N;i++) dep[i]=1e9;
40     dep[S]=0;q.push(S);
41     while(!q.empty()){
42         int p=q.front();q.pop();
43         for(int i=egh[p];i;i=eg[i].ne){
44             int b=eg[i].b;if(!eg[i].l||dep[b]<=dep[p]+1) continue;
45             dep[b]=dep[p]+1;q.push(b);
46         }
47     }
48     return dep[T]<1e9;
49 }
50 
51 ll dinic(int x,ll y){
52     if(x==T) return y;
53     ll tmp=y;
54     for(int &i=cur[x];i;i=eg[i].ne){
55         int b=eg[i].b;if(!eg[i].l||dep[b]!=dep[x]+1) continue;
56         ll r=dinic(b,min(tmp,(ll)eg[i].l));
57         tmp-=r,eg[i].l-=r,eg[i^1].l+=r;
58         if(!tmp) break;
59     }return y-tmp;
60 }
61 
62 int main(){
63     //freopen("","r",stdin);
64     N=rd(),M=rd(),S=rd(),T=rd();
65        for(int i=1;i<=M;i++){
66            int a=rd(),b=rd(),l=rd();
67            adeg(a,b,l),adeg(b,a,0);
68        }
69        ll ans=0;
70        while(bfs()){
71            memcpy(cur,egh,sizeof(cur));
72            ans+=dinic(S,1e15);
73        }
74        printf("%lld\n",ans);
75     return 0;
76 }
View Code

 

预流推进(HLPP)

考虑对某个点,尽可能地把它目前有的流量都推出去,当然这会导致流量不平衡,也就是说,每个点都会有一个额外流(我记做exc)

假设S的exc是inf,那么当除S和T外所有点的exc都为0时,T的exc的最大值就是最大流量

我们考虑给每个点附上一个高度h,并钦定只能从h推向h-1,以及钦定S的高度是N,T的高度是0不会改变

于是就可以每次取出最高的有exc的高度,然后推流,之后如果它的exc还有剩,就把它的h改成它出边中最小的h+1

所以用一个优先队列维护一下就好了

优化1:先用一次T到S的bfs,给每个点附上一个初始的h

优化2:如果某个高度x没有点了,那么所有高度大于x,小于N+1的点都不会再有机会流向T,所以直接把它们的高度设成N+1,所以拿一个cnt记一下每个高度的点的个数即可。记得cnt要开两倍大小

LOJ101
  1 #include<bits/stdc++.h>
  2 #define pa pair<int,int>
  3 #define CLR(a,x) memset(a,x,sizeof(a))
  4 #define MP make_pair
  5 #define fi first
  6 #define se second
  7 using namespace std;
  8 typedef long long ll;
  9 typedef unsigned long long ull;
 10 typedef unsigned int ui;
 11 typedef long double ld;
 12 const int maxn=105,maxm=5005;
 13 
 14 inline char gc(){
 15     return getchar();
 16     static const int maxs=1<<16;static char buf[maxs],*p1=buf,*p2=buf;
 17     return p1==p2&&(p2=(p1=buf)+fread(buf,1,maxs,stdin),p1==p2)?EOF:*p1++;
 18 }
 19 inline ll rd(){
 20     ll x=0;char c=gc();bool neg=0;
 21     while(c<'0'||c>'9'){if(c=='-') neg=1;c=gc();}
 22     while(c>='0'&&c<='9') x=(x<<1)+(x<<3)+c-'0',c=gc();
 23     return neg?(~x+1):x;
 24 }
 25 
 26 struct Edge{
 27     int b,l,ne;
 28 }eg[maxm*2];
 29 int egh[maxn],ect=1;
 30 int N,M,S,T;
 31 int h[maxn],cnt[maxn*2];
 32 ll exc[maxn];
 33 bool flag[maxn];
 34 struct cmp{
 35     inline bool operator ()(int a,int b) const{return h[a]<h[b];}
 36 };
 37 priority_queue<int,vector<int>,cmp> q;
 38 
 39 inline void adeg(int a,int b,int l){
 40     eg[++ect].b=b,eg[ect].l=l,eg[ect].ne=egh[a],egh[a]=ect;
 41 }
 42 
 43 inline bool bfs(){
 44     for(int i=1;i<=N;i++) h[i]=1e9;
 45     queue<int> q;
 46     q.push(T);h[T]=0;
 47     while(!q.empty()){
 48         int p=q.front();q.pop();
 49         for(int i=egh[p];i;i=eg[i].ne){
 50             int b=eg[i].b;if(!eg[i^1].l||h[b]<=h[p]+1) continue;
 51             h[b]=h[p]+1,q.push(b);
 52         }
 53     }return h[S]<1e9;
 54 }
 55 
 56 inline ll HLPP(){
 57     if(!bfs()) return 0;
 58     h[S]=N;flag[S]=flag[T]=1;
 59     for(int i=1;i<=N;i++) if(h[i]<1e9) cnt[h[i]]++;
 60     for(int i=egh[S];i;i=eg[i].ne){
 61         int b=eg[i].b,l=eg[i].l;if(!l) continue;
 62         exc[b]+=l,eg[i].l-=l,eg[i^1].l+=l;
 63         if(!flag[b]) q.push(b),flag[b]=1;
 64     }
 65     while(!q.empty()){
 66         int p=q.top();q.pop();
 67         flag[p]=0;
 68         for(int i=egh[p];i;i=eg[i].ne){
 69             int b=eg[i].b;if(!eg[i].l||h[b]!=h[p]-1) continue;
 70             int l=min(exc[p],(ll)eg[i].l);
 71             exc[p]-=l,exc[b]+=l,eg[i].l-=l,eg[i^1].l+=l;
 72             if(!flag[b]) q.push(b),flag[b]=1;
 73             if(!exc[p]) break;
 74         }
 75         if(exc[p]){
 76             if(!--cnt[h[p]]){
 77                 for(int i=1;i<=N;i++){
 78                     if(!flag[i]&&h[i]>h[p]&&h[i]<N+1) h[i]=N+1;
 79                 }
 80             }
 81             int mi=1e9;
 82             for(int i=egh[p];i;i=eg[i].ne){
 83                 int b=eg[i].b;if(!eg[i].l) continue;
 84                 mi=min(mi,h[b]+1);
 85             }
 86             cnt[h[p]=mi]++;
 87             q.push(p),flag[p]=1;
 88         }
 89     }
 90     return exc[T];
 91 }
 92 
 93 int main(){
 94     //freopen("","r",stdin);
 95     N=rd(),M=rd(),S=rd(),T=rd();
 96        for(int i=1;i<=M;i++){
 97            int a=rd(),b=rd(),l=rd();
 98            adeg(a,b,l),adeg(b,a,0);
 99        }
100        ll ans=HLPP();
101        printf("%lld\n",ans);
102     return 0;
103 }
View Code

 

费用流

spfa

考虑使用spfa来找到一条从S到T的费用和最小的路径,然后增广这条路径。记一下每个点最短距离取的那个入边即可

LOJ102
 1 #include<bits/stdc++.h>
 2 #define pa pair<int,int>
 3 #define CLR(a,x) memset(a,x,sizeof(a))
 4 #define MP make_pair
 5 #define fi first
 6 #define se second
 7 using namespace std;
 8 typedef long long ll;
 9 typedef unsigned long long ull;
10 typedef unsigned int ui;
11 typedef long double ld;
12 const int maxn=405,maxm=15005,inf=(1u<<31)-1;
13 
14 inline char gc(){
15     return getchar();
16     static const int maxs=1<<16;static char buf[maxs],*p1=buf,*p2=buf;
17     return p1==p2&&(p2=(p1=buf)+fread(buf,1,maxs,stdin),p1==p2)?EOF:*p1++;
18 }
19 inline ll rd(){
20     ll x=0;char c=gc();bool neg=0;
21     while(c<'0'||c>'9'){if(c=='-') neg=1;c=gc();}
22     while(c>='0'&&c<='9') x=(x<<1)+(x<<3)+c-'0',c=gc();
23     return neg?(~x+1):x;
24 }
25 
26 struct Edge{
27     int b,l,c,ne;
28 }eg[maxm*2];
29 int N,M,egh[maxn],ect=1,S,T;
30 int dis[maxn],ine[maxn];
31 bool flag[maxn];
32 queue<int> q;
33 
34 inline void adeg(int a,int b,int c,int l){
35     eg[++ect]=Edge{b,l,c,egh[a]},egh[a]=ect;
36 }
37 
38 inline bool spfa(){
39     for(int i=1;i<=N;i++) dis[i]=inf,flag[i]=0;
40     q.push(S);dis[S]=0;
41     while(!q.empty()){
42         int p=q.front();q.pop();
43         flag[p]=0;
44         for(int i=egh[p];i;i=eg[i].ne){
45             int b=eg[i].b;if(!eg[i].l||dis[b]<=dis[p]+eg[i].c) continue;
46             dis[b]=dis[p]+eg[i].c,ine[b]=i;
47             if(!flag[b]) q.push(b),flag[b]=1;
48         }
49     }return dis[T]<inf;
50 }
51 
52 int main(){
53     //freopen("","r",stdin);
54     N=rd(),M=rd(),S=1,T=N;
55     for(int i=1;i<=M;i++){
56         int a=rd(),b=rd(),l=rd(),c=rd();
57         adeg(a,b,c,l),adeg(b,a,-c,0);
58     }
59     int flw=0,cst=0;
60     while(spfa()){
61         int mi=inf;
62         for(int x=T;x!=S;x=eg[ine[x]^1].b) mi=min(mi,eg[ine[x]].l);
63         for(int x=T;x!=S;x=eg[ine[x]^1].b) eg[ine[x]].l-=mi,eg[ine[x]^1].l+=mi;
64         flw+=mi,cst+=mi*dis[T];
65     }
66     printf("%d %d\n",flw,cst);
67     return 0;
68 }
View Code

 

原始对偶

首先,我们刚刚的spfa实际上已经像dinic的bfs一样给图分了个层,所以我们可以仿照dinic的做法,按照这些层来选增广路(就是只能走最短路)。但要注意因为可能有环,已经走过的点就不能再走了,要打个标记,也正是因为这个标记,一次增广可能并不完全,所以我需要做多次增广直到找不到增广路。而且如果spfa是从S到T求的的话,要从T到S来增广,因为你只知道从S到x的距离,你再从S开始去增广的话,它有可能根本到不了T。

然而这还不够,我们想试着用更高效的方法来分层。

考虑dijkstra,然而图中有负边。但我们可以先进行一次spfa,之后每次都根据上次算出来的dis改变边权,使得它们(还有剩余流量的)全是非负的。

考虑连接i和j的边(i,j),一定有dis[i]+cost[(i,j)]>=dis[j],所以我们可以把(i,j)的边权变成dis[i]+cost[(i,j)]-dis[j]。之后我们会发现,S到T的路径的和就变成了原来的和加上dis[S]减去dis[T]。于是就可以dijkstra了

但要注意的是,每次dijkstra完都要更新边权,因为会出现新边(就是说 原来没有剩余流量的边现在有剩余流量了,而因为它原来没有剩余流量 所以它的边权并不一定是非负的)

LOJ102
  1 #include<bits/stdc++.h>
  2 #define pa pair<int,int>
  3 #define CLR(a,x) memset(a,x,sizeof(a))
  4 #define MP make_pair
  5 #define fi first
  6 #define se second
  7 using namespace std;
  8 typedef long long ll;
  9 typedef unsigned long long ull;
 10 typedef unsigned int ui;
 11 typedef long double ld;
 12 const int maxn=405,maxm=15005,inf=(1u<<31)-1;
 13 
 14 inline char gc(){
 15     return getchar();
 16     static const int maxs=1<<16;static char buf[maxs],*p1=buf,*p2=buf;
 17     return p1==p2&&(p2=(p1=buf)+fread(buf,1,maxs,stdin),p1==p2)?EOF:*p1++;
 18 }
 19 inline ll rd(){
 20     ll x=0;char c=gc();bool neg=0;
 21     while(c<'0'||c>'9'){if(c=='-') neg=1;c=gc();}
 22     while(c>='0'&&c<='9') x=(x<<1)+(x<<3)+c-'0',c=gc();
 23     return neg?(~x+1):x;
 24 }
 25 
 26 struct Edge{
 27     int b,l,c,ne;
 28 }eg[maxm*2];
 29 int N,M,egh[maxn],ect=1,S,T;
 30 int dis[maxn],off;
 31 bool flag[maxn];
 32 
 33 inline void adeg(int a,int b,int c,int l){
 34     eg[++ect]=Edge{b,l,c,egh[a]},egh[a]=ect;
 35 }
 36 
 37 inline bool spfa(){
 38     static queue<int> q;
 39     for(int i=1;i<=N;i++) dis[i]=inf,flag[i]=0;
 40     q.push(S);dis[S]=0;
 41     while(!q.empty()){
 42         int p=q.front();q.pop();
 43         flag[p]=0;
 44         for(int i=egh[p];i;i=eg[i].ne){
 45             int b=eg[i].b;if(!eg[i].l||dis[b]<=dis[p]+eg[i].c) continue;
 46             dis[b]=dis[p]+eg[i].c;
 47             if(!flag[b]) q.push(b),flag[b]=1;
 48         }
 49     }
 50     for(int i=2;i<=ect;i++){
 51         eg[i].c+=dis[eg[i^1].b]-dis[eg[i].b];
 52     }
 53     off+=dis[T]-dis[S];
 54     return dis[T]<inf;
 55 }
 56 
 57 inline bool dijk(){
 58     static priority_queue<pa,vector<pa>,greater<pa> > q;
 59     for(int i=1;i<=N;i++) dis[i]=inf,flag[i]=0;
 60     q.push(MP(0,S)),dis[S]=0;
 61     while(!q.empty()){
 62         int p=q.top().se;q.pop();
 63         if(flag[p]) continue;
 64         flag[p]=1;
 65         for(int i=egh[p];i;i=eg[i].ne){
 66             int b=eg[i].b;if(!eg[i].l||dis[b]<=dis[p]+eg[i].c) continue;
 67             dis[b]=dis[p]+eg[i].c;
 68             q.push(MP(dis[b],b));
 69         }
 70     }
 71     for(int i=2;i<=ect;i++){
 72         eg[i].c+=dis[eg[i^1].b]-dis[eg[i].b];
 73     }
 74     off+=dis[T]-dis[S];
 75     return dis[T]<inf;
 76 }
 77 
 78 int dinic(int x,int y){
 79     if(x==S) return y;
 80     int tmp=y;flag[x]=1;
 81     for(int i=egh[x];i;i=eg[i].ne){
 82         int b=eg[i].b;if(flag[b]||!eg[i^1].l||eg[i^1].c) continue;
 83         int r=dinic(b,min(tmp,eg[i^1].l));
 84         tmp-=r,eg[i^1].l-=r,eg[i].l+=r;
 85         if(!tmp) break;
 86     }return y-tmp;
 87 }
 88 
 89 int main(){
 90     //freopen("","r",stdin);
 91     N=rd(),M=rd(),S=1,T=N;
 92     for(int i=1;i<=M;i++){
 93         int a=rd(),b=rd(),l=rd(),c=rd();
 94         adeg(a,b,c,l),adeg(b,a,-c,0);
 95     }
 96     int flw=0,cst=0;
 97     spfa();
 98     while(dijk()){
 99         int a=-1;
100         while(a){
101             CLR(flag,0);
102             a=dinic(T,inf);
103             flw+=a,cst+=a*off;
104         }
105     }
106     printf("%d %d\n",flw,cst);
107     return 0;
108 }
View Code

 

posted @ 2019-07-05 10:08  Ressed  阅读(556)  评论(3编辑  收藏  举报