网络流(最大流和费用流)小结
网络流
(持续复健中)
网络(流网络):一个有向图,每条边都有一定的容量,图上还有一个源点以及汇点
流从源点流向汇点,满足如下性质
- 每条边流过的流量大小不得超过边的容量
- 从源点流出的总流量等于汇点流入的总流量
最大流
最大流问题就是求解从源点最多能流多少流量到汇点的问题
FF(Ford-Fulkerson)方法
残量网络:在一个网络中所有剩余容量大于 \(0\) 的边构成的以及所有节点构成的子图
增广路:一条从源点到汇点的路径如果满足这条路径上所有边的剩余容量都大于 \(0\),那这条边就被称为增广路,顾名思义,流量从这条边上流过可以增加流到汇点的总流量,也就是可以实现增广,从定义不难看出,残量网络上的所有从源点到汇点的路径都是增广路
FF 方法就是通过不断地寻找增广路来得到最大流
下面来看常用的求解最大流的具体算法
EK 算法
首先对于每条有向边我们建立一条反向边,容量一开始为 \(0\),始终保证流量与这条边的流量的和为 \(0\),也就是说正向边流过一定的流量,反向边就增加这么多的容量表示可以“反悔”
从源点开始,我们不断地 BFS,带着流量进行增广,每次只要找到了汇点,就把从源点到汇点的这条增广路上的所有边都增加这么多流量就可以了
这样每次 BFS 只能找到一条增广路,效率有点低下,复杂度为 \(O(nm^2)\),\(n\) 为节点数,\(m\) 为边数(证明可以看看 这个)
实现我就不写了(太懒了)
Dinic 算法
在每次增广前,我们都用 BFS 将残量网络分层一次,一个点的层数就是他到源点最近的距离,这样如果一次 BFS 后发现无法从残量网络上到达汇点,我们就直接停止增广,在找增广路时,我们进行 DFS,每次一定只是往层数比当前点层数大 \(1\) 的点进行增广,这样就可以保证我们每次找到的增广路是最短的(至于为什么要这样做和复杂度分析有关)
Dinic 算法可以实现多路增广,如果在找一条增广路的过程中我们发现当前节点还有一些剩余的流量没有用完,那就可以接着从它的出边继续找增广路,一次 DFS 就可以实现同时找到多条增广路
还可以在此基础上再添加当前弧优化,如果一条边已经被增广过,那么它就不可能被再次增广,那么,我们下一次进行增广的时候如果找到这个边,就可以不必再走那些已经被增广过的边,直接从接下来的边开始继续增广
这样的复杂度就是 \(O(n^2m)\),\(n\) 为节点数,\(m\) 为边数,稠密图上的效率就被大大提高了,特别地,如果用 Dinic 来解决二分图最大匹配问题,复杂度为 \(O(m\sqrt n)\)
复杂度证明可以看看 这个
下面是代码(Luogu P3376 【模板】网络最大流),里面还有一些小剪枝可以优化一点常数,比如一个节点如果当前所有流量都用完了就不能再从它继续多路增广了,具体可以看代码实现
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N=211,M=1e4+5;
const ll INF=2e18;
int n,m,s,t;
int head[N],cur[N],nxt[M],to[M],cnt=1;//邻接表存图,cnt设置为1便于找反向边
ll cap[M],fl[M];//边的容量、边已有的流量
int dep[N];
inline void create(int ff,int tt,ll vv){nxt[++cnt]=head[ff],head[ff]=cnt,to[cnt]=tt,cap[cnt]=vv;}
inline bool bfs()//BFS分层
{
queue<int>q;
for(int i=1;i<=n;++i)cur[i]=head[i],dep[i]=0;
q.push(s);
dep[s]=1;
int tmp;
while(!q.empty())
{
tmp=q.front();
q.pop();
for(int i=head[tmp];i;i=nxt[i])
{
if(dep[to[i]]||cap[i]==fl[i])continue;//已经找到过或者没有流量了就不能走
dep[to[i]]=dep[tmp]+1;
q.push(to[i]);
}
}
return dep[t];
}
inline ll dfs(int pos,ll flow)//DFS找增广路,返回流到终点的流量
{
if(pos==t||!flow)return flow;//流到终点或者没有流量了就返回
ll suc=0,delta;
for(int &i=cur[pos];i;i=nxt[i])//i为引用,实现当前弧优化
{
if(dep[pos]+1==dep[to[i]]&&cap[i]>fl[i])//只有层数+1且还有流量才能流
{
delta=dfs(to[i],min(cap[i]-fl[i],flow));//delta表示流出了多少
suc+=delta;//成功流过流量增加
fl[i]+=delta;//当前边已有流量增加
fl[i^1]-=delta;//反向边已有流量减少
flow-=delta;//当前还剩的可流流量减少(后面的剩余流量可以继续多路增广)
if(!flow)break;//不能多路增广就可以退出了
}
}
return suc;
}
inline ll dinic()
{
ll ans=0;
while(bfs())ans+=dfs(s,INF);//每次分层,从起点以无限流量往下流
return ans;
}
signed main()
{
ios::sync_with_stdio(false);
cin.tie(0),cout.tie(0);
cin>>n>>m>>s>>t;
ll vv;
for(int i=1,ff,tt;i<=m;++i)cin>>ff>>tt>>vv,create(ff,tt,vv),create(tt,ff,0);
cout<<dinic();
return 0;
}
ISAP 算法
在 Dinic 算法中每次 DFS 都需要 BFS 分层一次,非常麻烦,ISAP 就在此方面进行了优化
ISAP 同样要执行分层,不过只执行一次,它是在 DFS 时更新层数
ISAP 在反图上先 BFS 分层一次(即从汇点开始 BFS),在增广时就只选择比当前点层数少 \(1\) 的节点来增广,和 Dinic 的区别就是在结束点 \(i\) 的增广后,我们遍历残量网络上 \(i\) 的所有出边,令 \(i\) 的层数变为层数最小的出点的层数加 \(1\),若没有出边那就直接让 \(i\) 的层数变为 \(n+1\)(大于最大值即可),不难发现若源点的层数大于 \(n\) 就说明已经无法通过残量网络到达汇点了
ISAP 不支持多路增广,每次找到一条增广路后才会回溯回源点并修改路径上的流量,但是它和 Dinic 一样有当前弧优化,此外 ISAP 还可以用一个数组记录每一个层数所拥有的点数,一旦中间有一种点数出现次数变为 \(0\),那就说明图上已经有了断层,一定无法找到增广路,可以直接终止算法
复杂度和 Dinic 是一样的,只是常数会小很多
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N=211,M=1e4+5;
const ll INF=2e18;
int n,m,s,t;
int head[N],fr[M],cur[N],nxt[M],to[M],cnt=1;
ll cap[M],fl[M];
int dep[N],num[N],pre[N];
inline void create(int ff,int tt,ll vv){nxt[++cnt]=head[ff],fr[cnt]=ff,head[ff]=cnt,to[cnt]=tt,cap[cnt]=vv;}
inline void bfs()//一遍BFS分层
{
queue<int>q;
for(int i=1;i<=n;++i)cur[i]=head[i],dep[i]=0;
q.push(t);
dep[t]=1;
int tmp;
while(!q.empty())
{
tmp=q.front();
q.pop();
for(int i=head[tmp];i;i=nxt[i])
{
if(dep[to[i]])continue;
dep[to[i]]=dep[tmp]+1;
q.push(to[i]);
}
}
}
inline ll getflow()//找到增广路后回溯回源点并修改增广路上的路径流量
{
int pos=t,edge;
ll ans=INF;
while(pos!=s)
{
edge=pre[pos];
ans=min(ans,cap[edge]-fl[edge]);//计算出这条增广路能够流过多少流量
pos=fr[edge];
}
pos=t;
while(pos!=s)
{
fl[pre[pos]]+=ans;//修改流量
fl[pre[pos]^1]-=ans;
pos=fr[pre[pos]];
}
return ans;
}
inline ll isap()
{
ll ans=0;
bfs();
for(int i=1;i<=n;++i)++num[dep[i]];//记录每种层数对应的点数
int pos=s;
while(dep[s]<=n)//源点不可到达其他节点就停止
{
if(pos==t) ans+=getflow(),pos=s;//找到增广路就流动并返回起点继续下一次DFS
bool isok=0;
for(int &i=cur[pos];i;i=nxt[i])
{
if(dep[to[i]]+1==dep[pos]&&cap[i]>fl[i])//可以流动
{
isok=1;//可以增广
pre[to[i]]=i;//记录从哪条边来的,以支持回溯
pos=to[i];
break;
}
}
if(!isok)//如果这个点完全不能再增广就就更新层数
{
int newdep=n;
for(int i=head[pos];i;i=nxt[i])
if(cap[i]>fl[i])
newdep=min(newdep,dep[to[i]]);
if(!(--num[dep[pos]]))break;
++num[dep[pos]=newdep+1];//深度变为到点深度+1,如果残量网络上没有出边,那就是n+1
cur[pos]=head[pos];
if(pos!=s)pos=fr[pre[pos]];//回到上一条边继续寻找增广路
}
}
return ans;
}
signed main()
{
ios::sync_with_stdio(false);
cin.tie(0),cout.tie(0);
cin>>n>>m>>s>>t;
ll vv;
for(int i=1,ff,tt;i<=m;++i)cin>>ff>>tt>>vv,create(ff,tt,vv),create(tt,ff,0);
cout<<isap();
return 0;
}
其他算法
(最高标号)预流推进等算法能做到更加优秀的时间复杂度,但是我们知道,卡 Dinic 的出题人都不是好出题人,而最大流或者费用流问题主要的考察点就是如何建图,对于算法的要求并不高,所以这里就不说了(更何况 CSP,NOIP 都不考)
费用流
顾名思义,费用流就是有费用了,现在对于网络上的每一条边我们都增加一个费用,这条边上的花费就是流过的流量乘以这个费用,要求在保证流是最大流的情况下使花费的费用最少,这就是最小费用最大流
找费用最小增广路法
怎么做呢?考虑贪心,正常的最大流就是直接找增广路,每次找到增广路就沿着增广路流动,现在要求费用最小,那么我们就在找增广路的时候先保证它是流过每单位流量费用最小的增广路就可以了(之一)
使用 SPFA (可能会出现负权环,Dijkstra 不支持)就可以了,这里给出使用 SPFA 寻找增广路的 EK 算法供参考
#include<bits/stdc++.h>
#define INF 2100000000
#define MAXN 5003
#define MAXM 100005
using namespace std;
int n,m,s,t,max_flow,min_cost;
int head[MAXN],nxt[MAXM],to[MAXM],fl[MAXM],cost[MAXM],cnt=1;
int minfl[MAXN],dis[MAXN],pre_point[MAXN],pre_edge[MAXN];
bool vis[MAXN];
inline void create(int f,int t,int cap,int c)
{
nxt[++cnt]=head[f],head[f]=cnt,to[cnt]=t,fl[cnt]=cap,cost[cnt]=c;
nxt[++cnt]=head[t],head[t]=cnt,to[cnt]=f,fl[cnt]=0,cost[cnt]=-c;
}
bool spfa()//SPFA找单位费用最小增广路
{
for(int i=0;i<=n;++i) dis[i]=INF,vis[i]=0;
queue<int>q;
int cur;
q.push(s);
vis[s]=1;
dis[s]=0;
minfl[s]=INF;
while(!q.empty())
{
cur=q.front();
q.pop();
vis[cur]=0;
for(int i=head[cur];i;i=nxt[i])
if(fl[i]>0&&cost[i]+dis[cur]<dis[to[i]])//单位费用更小并且能流动才进入
{
dis[to[i]]=cost[i]+dis[cur];
minfl[to[i]]=min(minfl[cur],fl[i]);//修改最小流量
pre_point[to[i]]=cur,pre_edge[to[i]]=i;//记录来源以支持回溯
if(!vis[to[i]]) vis[to[i]]=1,q.push(to[i]);
}
}
return dis[t]!=INF;
}
void find()
{
int cur=t,i;
while(cur!=s)//回溯并修改
{
i=pre_edge[cur];
fl[i]-=minfl[t];
fl[i^1]+=minfl[t];
cur=pre_point[cur];
}
max_flow+=minfl[t];//更新答案
min_cost+=dis[t]*minfl[t];
}
signed main()
{
ios::sync_with_stdio(false);
cin.tie(0),cout.tie(0);
cin>>n>>m>>s>>t;
int f,t,cap,c;
while(m--) cin>>f>>t>>cap>>c,create(f,t,cap,c);
while(spfa())find();
cout<<max_flow<<' '<<min_cost;
return 0;
}
原始对偶算法
如果我非要用 Dijkstra 呢?要使用 Dijkstra 并优化复杂度需要保证没有负权环,要想办法把每条边的边权变为非负数
原始对偶(Primal-Dual)算法就是这样一个方法,它给每个点增加了一个“势能”,利用势能和边权来把边权全部变为非负数
下面看看具体实现,首先用最短路算法求出源点到每个点的最短距离 \(h_i\) 作为一个点的初始势能,把的从点 \(u\) 到点 \(v\) 的单位费用为 \(w_i\) 边 \(i\) 的边权设置为 \(w_i+h_u-h_v\),每次增广后重新将每个点的势能更新并计算新的边权,势能 \(h_i\) 更新为 \(h_i+d_i\) 即可(\(d_i\) 表示在源点到当前节点最短路径上的更新后的边权的和,也就是 \(i\) 到源点的新的最短距离)
为什么这样是对的呢?还是看看 大佬的 Blog 吧(逃)
#include<bits/stdc++.h>
using namespace std;
const int N=5e3+5,M=1e5+5,INF=2e9;
int n,m,s,t;
int maxflow=0,mincost=0;
int head[N],fr[M],nxt[M],to[M],cnt=1;
int cap[M],w[M],fl[M],v[M],h[N]/*势能*/,dis[N];
int pre[N];
bool vis[N];
inline void create(int ff,int tt,int cc,int ww){nxt[++cnt]=head[ff],fr[cnt]=ff,head[ff]=cnt,to[cnt]=tt,cap[cnt]=cc,w[cnt]=ww;}
void spfa()//SPFA算初始势能
{
static queue<int>q;
for(int i=1;i<=n;++i)h[i]=INF;
h[s]=0,vis[s]=1;
q.push(s);
int tmp;
while(q.size())
{
tmp=q.front();
q.pop();
vis[tmp]=0;
for(int i=head[tmp];i;i=nxt[i])
{
if(cap[i]&&h[to[i]]>h[tmp]+w[i])
{
h[to[i]]=h[tmp]+w[i];
if(!vis[to[i]]) vis[to[i]]=1,q.push(to[i]);
}
}
}
}
bool dij()//找增广路
{
typedef pair<int,int> pii;
static priority_queue<pii,vector<pii>,greater<pii>>q;
for(int i=1;i<=n;++i)dis[i]=INF,vis[i]=0;
dis[s]=0;
q.push(make_pair(0,s));
int tmp;
while(!q.empty())
{
tmp=q.top().second;
q.pop();
if(vis[tmp])continue;
vis[tmp]=1;
for(int i=head[tmp],neww;i;i=nxt[i])
{
neww=w[i]+h[tmp]-h[to[i]];//更新边权
if(cap[i]>fl[i]&&dis[to[i]]>dis[tmp]+neww)
{
dis[to[i]]=dis[tmp]+neww;
pre[to[i]]=i;
if(!vis[to[i]])q.push(make_pair(dis[to[i]],to[i]));
}
}
}
return dis[t]!=INF;
}
signed main()
{
ios::sync_with_stdio(false);
cin.tie(0),cout.tie(0);
cin>>n>>m>>s>>t;
for(int i=1,t1,t2,t3,t4;i<=m;++i)cin>>t1>>t2>>t3>>t4,create(t1,t2,t3,t4),create(t2,t1,0,-t4);
spfa();
while(dij())
{
int flow=INF;
for(int i=t;i!=s;i=fr[pre[i]])flow=min(flow,cap[pre[i]]-fl[pre[i]]);//回溯
for(int i=t;i!=s;i=fr[pre[i]])
{
fl[pre[i]]+=flow;
fl[pre[i]^1]-=flow;
}
maxflow+=flow;
mincost+=flow*(h[t]+dis[t]);
for(int i=1;i<=n;++i)h[i]+=dis[i];
}
cout<<maxflow<<' '<<mincost;
cout<<maxflow<<' '<<mincost;
return 0;
}
代码中更新答案最小花费时用的是 mincost+=flow*(h[t]+dis[t]);
很奇怪,为什么呢?
这里使用 \(h_t\) 更新的原因是:假设此时汇点到起点在原本的网络上的距离(单位花费)是 \(d_t\),那么
所以单位花费 \(d_t\) 就是 \(\text{dis}_t+h_t\)