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

小咪咪 2021-11-10 12:36 363阅读 0赞

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

## 最大流 ##

### dinic ###

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

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

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

##### LOJ101 #####

![ContractedBlock.gif][] ![ExpandedBlockStart.gif][]

### 预流推进(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 #####

![ContractedBlock.gif][] ![ExpandedBlockStart.gif][]

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 }

## 费用流 ##

### spfa ###

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

##### LOJ102 #####

![ContractedBlock.gif][] ![ExpandedBlockStart.gif][]

### 原始对偶 ###

首先，我们刚刚的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 #####

![ContractedBlock.gif][] ![ExpandedBlockStart.gif][]

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 }

转载于:https://www.cnblogs.com/Ressed/p/11136619.html

[ContractedBlock.gif]: https://images.cnblogs.com/OutliningIndicators/ContractedBlock.gif
[ExpandedBlockStart.gif]: /images/20211109/d914aac673e14e6db1e3f147285f29e5.png