EK不够快？再学个Dinic吧

[上次的网络流EK讲解](https://www.luogu.org/blog/ONE-PIECE/wang-lao-liu-di-zong-jie)应该都能听懂吧。 上次讲了EK，有好多人想让我讲一下Dinic，那我就讲解一下Dinic吧~~（顺便看看还能不能再上洛谷日报）~~ # 特别感谢：@ComeIntoPower提出了许多修改建议 ------------ ## 什么是Dinic？ 百度百科上的定义是这样的： Dinic算法（又称Dinitz算法）是一个在网络流中计算最大流的强多项式复杂度的算法，设想由以色列（前苏联）的计算机科学家Yefim (Chaim) A. Dinitz在1970年提出。 ~~定义没什么用。~~ 按我的理解，Dinic是这样一种算法： 目的是找最大流，并且很高效。~~（废话）~~ ------------ ## 时间复杂度： 相比于EK的**O(nm^2**)(n为点数，m为边数)，DInic可以达到**O(n^2m)**，在稀疏图上，两者相差不大，但在稠密图上（二分匹配之类的）Dinic的优势就很明显了。 对于一些特殊情况， 在具有单位容量的网络中，可以证明每次寻找增广路的复杂度时间为O（m），并且分层次数不超过sqrt（m）和n^(2/3),因此复杂度为O(min(sqrt(m),n^(2/3))m) 在二分图匹配中，复杂度则为O（sqrt（n）m） ~~这些特殊情况的分析是网上找到，我也不会证明，总之Dinic很快，一般不会被卡，除了某些[恶意造数据的题](https://www.luogu.org/problemnew/show/P4722)，如果你的Dinic被卡了，可能是图建的不恰当~~ ------------ # Dinic讲解 有了之前EK的讲解，现在看Dinic的讲解应该很好理解，但以防万一，我还是做了图片 不知道你之前学习EK时有没有想过这样的问题：寻找增广路为什么只能一条一条找？可不可以一次找多条？一次找多条增广路可以实现吗？如果能，效率怎么样？ 没错，这就是Dinic的主要思想： **多路增广** 具体怎么实现呢？ Dinic算法分为两个步骤： 1. bfs分层(在EK中bfs是用于寻找增广路的) 2. dfs增广~~(dfs？EK中貌似没有这玩意啊，确定能高效？)~~ 3. ~~咦！刚才不是说两个步骤吗？~~重复执行1.2.直到图中无增广路为止 什么意思呢？ 以下图为例：  ~~反向边已经加好了，不用讲你也知道反向边是干嘛的吧。~~ 与EK一样，我们仍要通过bfs来判断图中是否还存在增广路，但是DInic算法里的bfs略有不同，这次，我们不用记录路径，而是给每一个点分层，对于任意点i，从s到i每多走过一个点，就让层数多一。 代码与EK略有不同： ``` bool bfs(){ memset(dep,0x3f,sizeof(dep)); memset(inque,0,sizeof(inque)); dep[s]=0; queue<int>q; q.push(s); while(!q.empty()){ int u=q.front(); q.pop(); inque[u]=0; for(int i=head[u];i;i=node[i].next){ int d=node[i].v; if(dep[d]>dep[u]+1&&node[i].val){ dep[d]=dep[u]+1;//注意与EK的区别 if(inque[d]==0){ q.push(d); inque[d]=1; } } } } if(dep[t]!=0x3f3f3f3f)return 1; return 0; }//给增广路上的点分层 ``` 分完层效果是这样的：  蓝色的数字是每个点层数 分层有什么用？ 有了每个点的层数编号，对任意点u到点d的路径如果有dep[d]==dep[u]+1，我们就可以判断该路径在一条**最短增广路**上。 为什么要找**最短增广路**？ 举个极端例子：  有了分层，我们就不会选s->1->2->4->5->3->t了 刚才说了的，分完层下一步是dfs增广。 在Dinic中，我们找增广路是用深搜： ~~由于这部分代码应该一看就懂，~~就先放代码，主要讲思想 ``` int dfs(int u,int flow){//u是当前点，low是当前增广路上的最小边权 int rlow=0; if(u==t)return flow;//已经到达t可以用当前路径 for(int i=head[u];i;i=node[i].next){ int d=node[i].v; if(node[i].val&&dep[d]==dep[u]+1){//寻找最短可增广路径 if(rlow=dfs(d,min(flow,node[i].val))){ node[i].val-=rlow; node[i^1].val+=rlow; //正向边加流，反向边减流，思想同EK return rlow; } } } return 0;//无法到达t，无增广路 }//寻找增广路 ``` ~~感觉这玩意应该会很慢，怎么可能比EK快？别急，讲完你就知道了~~ 再放一次刚才的图：  分完了层就要开始找增广路了。 比如说，第一次我们找s->1->4->t：  路径已经标出来了，由于作者画图时手比较抖~~（主要原因）~~，外加电脑自带的画图软件功能差，略丑。 再仔细看一看图，发现了什么？ 还有增广路（显然），标号还可以继续利用！！！那我们可以再执行一次dfs 继续利用第一次bfs出来的标号，再找第二条增广路： s->1->5->t:  再找找 竟然还能继续利用标号找第三条增广，再执行dfs s->1->3->t:  还有第四条！再执行dfs s->2->3->t:  发现了吗？一次bfs我们找了4条增广路！ 多么高效，这就是Dinic算法。 思想讲完了就附上我丑陋的代码： ``` //Dinic网络最大流模板 #include<iostream> #include<cstdio> #include<cstring> #include<algorithm> #include<queue> using namespace std; const int inf=1<<30; int dep[10101],head[10101],inque[10101]; int top=1; int maxflow=0; int n,m,s,t; struct Node{ int v; int val; int next; }node[200100]; inline void addedge(int u,int v,int val){ node[++top].v=v; node[top].val=val; node[top].next=head[u]; head[u]=top; } inline int Read(){ int x=0;char c=getchar(); while(c>'9'||c<'0')c=getchar(); while(c>='0'&&c<='9')x=x*10+c-'0',c=getchar(); return x; } bool bfs(){ memset(dep,0x3f,sizeof(dep)); memset(inque,0,sizeof(inque)); dep[s]=0; queue<int>q; q.push(s); while(!q.empty()){ int u=q.front(); q.pop(); inque[u]=0; for(int i=head[u];i;i=node[i].next){ int d=node[i].v; if(dep[d]>dep[u]+1&&node[i].val){ dep[d]=dep[u]+1;//注意与EK的区别 if(inque[d]==0){ q.push(d); inque[d]=1; } } } } if(dep[t]!=0x3f3f3f3f)return 1; return 0; }//给增广路上的点分层 inline int min(int x,int y){ return x<y?x:y; } int dfs(int u,int flow){ int rlow=0; if(u==t)return flow; for(int i=head[u];i;i=node[i].next){ int d=node[i].v; if(node[i].val&&dep[d]==dep[u]+1){ if(rlow=dfs(d,min(flow,node[i].val))){ node[i].val-=rlow; node[i^1].val+=rlow; return rlow; } } } return 0; }//寻找增广路 int Dinic(){ int lowflow; while(bfs()){ while(lowflow=dfs(s,inf))maxflow+=lowflow; } return maxflow; }//Dinic寻找最大流 int main(){ n=Read(),m=Read(),s=Read(),t=Read(); register int i; int u,v,val; for(i=1;i<=m;i++)u=Read(),v=Read(),val=Read(),addedge(u,v,val),addedge(v,u,0); printf("%d",Dinic()); return 0; } ``` 补充说明：Dinic在跑二分图匹配时比匈牙利快很多。 ------------ ## 优化： 上面讲的是普通的Dinic（一般在网上博客找到的版本），其实还可以优化。 其实只要加一点小小的改变，我们就可以在*dfs*时直接实现*多路增广*！ 具体怎么实现呢？ 我们定义一个变量vis，表示是否能到达t，若能到达t则可以再次dfs，否则重新bfs。（看起来没什么用） 再定义一个变量used，用来表示这个点的流量用了多少。 然后在dfs时我们可以在找到一条增广路时不直接返回，而是改变used的值，如果used还没达到该点流量上限fow，可以继续找别的增广路。 至于代码： ``` int dfs(int u,int flow){ int rlow=0; if(u==t){ vis==1;//可以到达t maxflow+=flow;//其实可以直接在这里累加最大流 return flow; } int used=0;//该点已经使用的流量 for(int i=head[u];i;i=node[i].next){ int d=node[i].v; if(node[i].val&&dep[d]==dep[u]+1){ if(rlow=dfs(d,min(flow-used,node[i].val))){ used+=rlow;//该点使用的流量增加 node[i].val-=rlow; node[i^1].val+=rlow; if(used==flow)break;//该点流量满了，没必要再找了 } } } return used;//返回该点已使用流量 }//寻找增广路 ``` 与之前不同的部分打了注释，请放心食用 当然主函数也要略作修改： ``` int Dinic(){ while(bfs()){ vis=1; while(vis==1){ vis=0; dfs(s,inf); } } return maxflow; }//Dinic寻找最大流 ``` 其实也可以直接这样：（上面的只是为了便于理解） ``` int Dinic(){ while(bfs()){ dfs(s,inf); } return maxflow; }//Dinic寻找最大流 ``` ------------ ## 优化2.0：（传说中的当前弧优化） 我们定义一个数组cur记录当前边（弧）（**功能类比邻接表中的head数组，只是会随着dfs的进行而修改**）， 每次我们找过某条边（弧）时，修改cur数组，改成该边（弧）的编号， 那么下次到达该点时，会直接从cur对应的边开始（**也就是说从head到cur中间的那一些边（弧）我们就不走了**）。 有点抽象啊，感觉并不能加快，然而实际上确实快了很多。 在代码中的实现，bfs中： 只需将原来的 ``` memset(dep,0x3f,sizeof(dep)); memset(inque,0,sizeof(inque)); ``` 修改为 ``` for(register int i=0;i<=n;i++)cur[i]=head[i],dep[i]=0x3f3f3f3f,inque[i]=0; //两个memset一个memcpy还不如直接循环 ``` ## 在dfs中 将这一部分： ``` for(int i=head[u];i;i=node[i].next){ int d=node[i].v; if(node[i].val&&dep[d]==dep[u]+1){ if(rlow=dfs(d,min(flow-used,node[i].val))){ used+=rlow; node[i].val-=rlow; node[i^1].val+=rlow; if(used==flow)break; } } } ``` 修改为： ``` for(int i=cur[u];i;i=node[i].next){ cur[u]=i;//修改当前弧 int d=node[i].v; if(node[i].val&&dep[d]==dep[u]+1){ if(rlow=dfs(d,min(flow-uesd,node[i].val))){ used+=rlow; node[i].val-=rlow; node[i^1].val+=rlow; if(used==flow)break; } } } ``` 完整的弧优化代码我就不放了，按上面的改优化一就可以了。 一点点小小的改动就能~~极大程度~~提高你的代码速度 ## 这里还要讲一下当前弧优化的原理： 首先，我们在按顺序dfs时，先被遍历到的边肯定是已经增广过了（或者已经确定无法继续增广了），那么这条边就可以视为“**废边**” 那么下次我们再到达该节点时，就可以直接无视掉所有废边，只走还有用的边，也就是说，每次dfs结束后，下次dfs可以更省时间。 # 感谢观看 ------------ 等等，讲完了吗？ 其实还没有， 这种思想也可以运用到*费用流*上，而且肯定比EK快啊。（大多数情况） 第一步当然和EK一样，要跑一遍spfa（要不然怎么保证是最小费用呢） 这里用了一个思想：只要一个从u来的点d满足 dist[d]==dist[u]+node[i].w 就可以表示该点在最短路上（~~应该很好理解，~~和上面的dep[d]==dep[u]+1类似），那我们在增广时就可以走这条边。所以只要把dfs时的条件略作修改即可。 然后与dinic一样，跑dfs进行增广，同时累加最大流，最小费用，同时正向边减流，反向边加流。~~是不是很简单啊，~~直接上代码： ``` int dfs(int u,int flow){ if(u==t){vis[t]=1;maxflow+=flow;return flow;}//可以到达t，加流 int used=0;//该点使用了的流量 vis[u]=1; for(int i=head[u];i;i=node[i].next){ int d=node[i].v; if((vis[d]==0||d==t)&&node[i].val!=0&&dist[d]==dist[u]+node[i].w){ int minflow=dfs(d,min(flow-used,node[i].val)); if(minflow!=0)cost+=node[i].w*minflow,node[i].val-=minflow,node[i^1].val+=minflow,used+=minflow; //可以到达t，加费用，扣流量 if(used==flow)break;//小小的优化 } } return used; } ``` 关于used的解释在上文已经说过了，至于为什么这里的**vis**要用数组，那是因为这是**费用流**，可能存在某条边的费用为0，**如果不给每一个点打访问标记，可能会出现两个点来回跑的情况，**卡死你的程序。 **当然如果t被访问过了，还是可以访问的。** 主函数部分： ``` int mincostmaxflow(){//最小费用最大流 while(spfa()){ vis[t]=1; while(vis[t]){//能到达t就继续 memset(vis,0,sizeof(vis));//由于是数组，要memset dfs(s,inf); } } return maxflow; } ``` 完整代码： ``` #include<iostream> #include<cstring> #include<algorithm> #include<queue> #include<cstdio> using namespace std; const int inf=1<<30; int top=1,head[5010]; int n,m,s,t; int cost,maxflow; int vis[5010];//是否到达过该点 int dist[5010];//到t的单位费用 struct Node{ int val; int v; int next; int w; }node[101010]; inline int Read(){ int x=0; char c=getchar(); while(c>'9'||c<'0')c=getchar(); while(c>='0'&&c<='9')x=x*10+c-'0',c=getchar(); return x; } inline void addedge(int u,int v,int val,int w){ node[++top].v=v; node[top].w=w; node[top].val=val; node[top].next=head[u]; head[u]=top; } bool spfa(){ memset(vis,0,sizeof(vis)); memset(dist,0x3f,sizeof(dist)); dist[s]=0; vis[s]=1; queue<int>q; q.push(s); while(!q.empty()){ int u=q.front(); q.pop(); vis[u]=0; for(int i=head[u];i;i=node[i].next){ int d=node[i].v; if(dist[d]>dist[u]+node[i].w&&node[i].val){ dist[d]=dist[u]+node[i].w; if(vis[d]==0){ q.push(d); vis[d]=1; } } } } return dist[t]!=0x3f3f3f3f; }//spfa同EK inline int min(int x,int y){ return x<y?x:y; } int dfs(int u,int flow){ if(u==t){vis[t]=1;maxflow+=low;return flow;}//可以到达t，加流 int used=0;//该条路径可用流量 vis[u]=1; for(int i=head[u];i;i=node[i].next){ int d=node[i].v; if((vis[d]==0||d==t)&&node[i].val!=0&&dist[d]==dist[u]+node[i].w){ int minflow=dfs(d,min(flow-used,node[i].val)); if(minflow!=0)cost+=node[i].w*minflow,node[i].val-=minflow,node[i^1].val+=minflow,used+=minflow; //可以到达t，加费用，扣流量 if(used==flow)break; } } return used; } int mincostmaxflow(){ while(spfa()){ vis[t]=1; while(vis[t]){ memset(vis,0,sizeof(vis)); dfs(s,inf); } } return maxflow; } int main(){ n=Read(),m=Read(),s=Read(),t=Read(); int u,v,w,val; register int i; for(i=1;i<=m;i++){ u=Read(),v=Read(),val=Read(),w=Read(); addedge(u,v,val,w); addedge(v,u,0,-w); } mincostmaxflow(); printf("%d %d",maxflow,cost); return 0; } ``` ### 总结：在图稀疏时此算法与EK差不多，但当图较稠密时就明显比EK优，而且代码也没比EK长多少（只多了一个dfs），根据情况决定使用哪个吧。 究竟能快到什么程度呢？ 我做[洛谷P4003](https://www.luogu.org/problemnew/show/P4003) 时两种都用了一次~~（主要是第一次EK没吸氧T了）~~ 可以对比一下时间： [EK（还吸了氧气的）](https://www.luogu.org/record/show?rid=9615262) [这种算法（没吸氧气）](https://www.luogu.org/record/show?rid=9874702) ## 感谢观看 ## 欢迎各位dalao提出修改建议 [ISAP和HLPP已经写好了（作者快累死了）](https://www.luogu.org/blog/ONE-PIECE/jiu-ji-di-zui-tai-liu-suan-fa-isap-yu-hlpp)