From 1fb0708ee92cc1b37a91d2e50120ca1a3b21111a Mon Sep 17 00:00:00 2001 From: Yaoyao Date: Fri, 5 Jul 2019 13:33:32 +0800 Subject: [PATCH] =?utf8?q?=E6=8A=8A=E7=BD=91=E7=BB=9C=E5=92=8C=E7=BD=91?= =?utf8?q?=E7=BB=9C=E6=B5=81=E7=9A=84=E5=AE=9A=E4=B9=89=E5=81=9A=E5=8C=BA?= =?utf8?q?=E5=88=86?= MIME-Version: 1.0 Content-Type: text/plain; charset=utf8 Content-Transfer-Encoding: 8bit --- docs/graph/flow/max-flow.md | 438 ++------------------------------------------ 1 file changed, 19 insertions(+), 419 deletions(-) diff --git a/docs/graph/flow/max-flow.md b/docs/graph/flow/max-flow.md index f5838e9a..f57603cd 100644 --- a/docs/graph/flow/max-flow.md +++ b/docs/graph/flow/max-flow.md @@ -1,5 +1,9 @@ +网络流在 OI 中是显得尤为重要的。在《算法导论》中就用了 35 页来讲述网络流的知识，在这里，给大家介绍网络流中的一些基本知识。 + ## 网络 +首先，请分清楚**网络**与**网络流**的概念。 + 网络是指一个有向图 $G=(V,E)$。 每条边 $(u,v)\in E$ 都有一个权值 $c(u,v)$，称之为容量（Capacity），当 $(u,v)\notin E$ 时有 $c(u,v)=0$。 @@ -10,15 +14,13 @@ 设 $f(u,v)$ 定义在二元组 $(u\in V,v\in V)$ 上的实数函数且满足 -1. **容量限制**：$f(u,v)\leq c(u,v)$ -2. **斜对称性**：$f(u,v)=-f(v,u)$ -3. **流守恒性**：$\forall x\in V-\{s,t\},\sum_{(u,x)\in E}f(u,x)=\sum_{(x,v)\in E}f(x,v)$ - -那么 $f$ 称为网络 $G$ 的流函数。 +1. 容量限制：对于每条边，流经该边的流量不得超过该边的容量，即，$f(u,v)\leq c(u,v)$ +2. 斜对称性：每条边的流量与其相反边的流量之和为 0，即$f(u,v)=-f(v,u)$ +3. 流守恒性：从源点流出的流量等于汇点流入的流量，即$\forall x\in V-\{s,t\},\sum_{(u,x)\in E}f(u,x)=\sum_{(x,v)\in E}f(x,v)$ -对于 $(u,v)\in E$，$f(u,v)$ 称为边的**流量**，$c(u,v)-f(u,v)$ 称为边的**剩余容量**。 +那么 $f$ 称为网络 $G$ 的流函数。对于 $(u,v)\in E$，$f(u,v)$ 称为边的**流量**，$c(u,v)-f(u,v)$ 称为边的**剩余容量**。整个网络的流量为 $\sum_{(s,v)\in E}f(s,v)$，即**从源点发出的所有流量之和**。 -整个网络的流量为 $\sum_{(s,v)\in E}f(s,v)$，即**从源点发出的所有流量之和**。 +一般而言也可以把网络流理解为整个图的流量。而这个流量必满足上述三个性质。 *注*：流函数的完整定义为 @@ -30,424 +32,22 @@ f(u,v)=\left\{\begin{split} \end{split}\right. $$ -## 定义 - -我们有一张图，要求从源点流向汇点的最大流量（可以有很多条路到达汇点），就是我们的最大流问题。 - -求解最大流问题有三种常见算法：Edmonds-Karp 算法，Dinic 算法，以及 ISAP 算法。 - -## FF增广路算法 - -该方法通过寻找增广路来更新最大流，有 EK,dinic,SAP,ISAP 主流算法。 - -求解最大流之前，我们先认识以下增广路的概念。 - - **增广路** 指的是，从源点到汇点，只要有 $flow$ ($flow>0$) 流过去，这条路就是增广路。在一些最大流算法中，就是将这些路 **增广** （意思就是走掉这条路，带走的流量肯定就是这条路的最小流量），如图： - -![](./images/flow1.png) - -我们从 $4$ 到 $3$ ，肯定可以先从流量为 $20$ 的这条边先走。那么这条边就被走掉了，不能再选，总的流量为 $20$（现在）。然后我们可以这样选择： - -1. $4\rightarrow2\rightarrow3$ 这条 **增广路** 的总流量为 $20$ 。到 $2$ 的时候还是 $30$，到 $3$ 了就只有 $20$ 了。 - -2. $4\rightarrow2\rightarrow1\rightarrow3$ 这样子我们就很好的保留了 $30$ 的流量。 - -所以我们这张图的最大流就应该是 $20+30=50$ 。 - -求最大流是很简单的，接下来讲解求最大流的 3 种方法。 - -### Edmond-Karp 动能算法（EK 算法） - -这个算法很简单，就是 BFS **找增广路** ，然后对其进行 **增广** 。你可能会问，怎么找？怎么增广？ - -1. 找？我们就从源点一直 BFS 走来走去，碰到汇点就停，然后增广（每一条路都要增广）。我们在 BFS 的时候就注意一下流量合不合法就可以了。 - -2. 增广？其实就是按照我们找的增广路在重新走一遍。走的时候把这条路的能够成的最大流量减一减，然后给答案加上最小流量就可以了。 - -再讲一下 **反向边** 。增广的时候要注意建造反向边，原因是这条路不一定是最优的，这样子程序可以进行反悔。假如我们对这条路进行增广了，那么其中的每一条边的反向边的流量就是它的流量。 - -![](./images/flow2.png) - -讲一下一些小细节。如果你是用邻接矩阵的话，反向边直接就是从 $table[x,y]$ 变成 $table[y,x]$ 。如果是常用的链式前向星，那么在加入边的时候就要先加入反向边。那么在用的时候呢，我们直接 $i\operatorname{xor}1$ 就可以了 ( $i$ 为边的编号）。为什么呢？相信大家都是知道 $\operatorname{xor}$ 的，那么我们在加入正向边后加入反向边，就是靠近的，所以可以使用 $\operatorname{xor}$ 。我们还要注意一开始的编号要设置为 $tot=1$ ，因为边要从编号 $2$ 开始，这样子 $\operatorname{xor}$ 对编号 $2,3$ 的边才有效果。 - -EK 算法的时间复杂度为 $O(n^2m)$ （其中 $n$ 为点数， $m$ 为边数）。效率还有很大提升空间。 - -```cpp -#include -#include -#include -#include -#define INF 0x3f3f3f3f -using namespace std; -struct edge { - int v, w, next; -} e[200005]; -struct node { - int v, e; -} p[10005]; -int head[10005], vis[10005]; -int n, m, s, t, cnt = 1; -void addedge(int u, int v, int w) { - e[++cnt].v = v; - e[cnt].w = w; - e[cnt].next = head[u]; - head[u] = cnt; -} -bool bfs() { - queue q; - memset(p, 0, sizeof(p)); - memset(vis, 0, sizeof(vis)); - vis[s] = 1; - q.push(s); - while (!q.empty()) { - int cur = q.front(); - q.pop(); - for (int i = head[cur]; i; i = e[i].next) - if ((!vis[e[i].v]) && e[i].w) { - p[e[i].v].v = cur; - p[e[i].v].e = i; - if (e[i].v == t) return 1; - vis[e[i].v] = 1; - q.push(e[i].v); - } - } - return 0; -} -int main() { - scanf("%d%d%d%d", &n, &m, &s, &t); - for (int i = 1; i <= m; i++) { - int u, v, w; - scanf("%d%d%d", &u, &v, &w); - addedge(u, v, w); - addedge(v, u, 0); - } - int ans = 0; - while (bfs()) { - int minw = INF; - for (int i = t; i != s; i = p[i].v) minw = min(minw, e[p[i].e].w); - for (int i = t; i != s; i = p[i].v) { - e[p[i].e].w -= minw; - e[p[i].e ^ 1].w += minw; - } - ans += minw; - } - printf("%d\n", ans); - return 0; -} -``` - -### Dinic - -**Dinic 算法** 的过程是这样的：每次增广前，我们先用 BFS 来将图分层。设源点的层数为 0，那么一个点的层数便是它离源点的最近距离。 - -通过分层，我们可以干两件事情： - -1. 如果不存在到汇点的增广路（即汇点的层数不存在），我们即可停止增广。 -2. 确保我们找到的增广路是最短的。（原因见下文） - -接下来是 DFS 找增广路的过程。 - -我们每次找增广路的时候，都只找比当前点层数多 1 的点进行增广（这样就可以确保我们找到的增广路是最短的）。 - -Dinic 算法有两个优化： - -1. **多路增广** ：每次找到一条增广路的时候，如果残余流量没有用完怎么办呢？我们可以利用残余部分流量，再找出一条增广路。这样就可以在一次 DFS 中找出多条增广路，大大提高了算法的效率。 -2. **当前弧优化** ：如果一条边已经被增广过，那么它就没有可能被增广第二次。那么，我们下一次进行增广的时候，就可以不必再走那些已经被增广过的边。 - -设点数为 $n$ ，边数为 $m$ ，那么 Dinic 算法的时间复杂度是 $O(n^{2}m)$ ，在稀疏图上效率和 EK 算法相当，但在稠密图上效率要比 EK 算法高很多。 - -特别地，在求解二分图最大匹配问题时，可以证明 Dinic 算法的时间复杂度是 $O(m\sqrt{n})$ 。 - -```cpp -#include -#include -#include -#include -#define INF 0x3f3f3f3f -using namespace std; -struct edge { - int v, w, next; -} e[200005]; -int n, m, s, t, cnt = 1; -int head[100005], dep[100005], vis[100005], cur[100005]; -void addedge(int u, int v, int w) { - e[++cnt].v = v; - e[cnt].w = w; - e[cnt].next = head[u]; - head[u] = cnt; -} -bool bfs() { - queue q; - memset(dep, INF, sizeof(dep)); - memset(vis, 0, sizeof(vis)); - memcpy(cur, head, sizeof(head)); - dep[s] = 0; - vis[s] = 1; - q.push(s); - while (!q.empty()) { - int p = q.front(); - q.pop(); - vis[p] = 0; - for (int i = head[p]; i; i = e[i].next) - if (dep[e[i].v] > dep[p] + 1 && e[i].w) { - dep[e[i].v] = dep[p] + 1; - if (!vis[e[i].v]) { - vis[e[i].v] = 1; - q.push(e[i].v); - } - } - } - if (dep[t] == INF) return 0; - return 1; -} -int dfs(int p, int w) { - if (p == t) return w; - int used = 0; //已经使用的流量 - for (int i = cur[p]; i; i = e[i].next) //每条边都尝试找一次增广路 - { - cur[p] = i; //当前弧优化 - if (dep[e[i].v] == dep[p] + 1 && e[i].w) { - int flow = dfs(e[i].v, min(w - used, e[i].w)); - if (flow) { - used += flow; - e[i].w -= flow; - e[i ^ 1].w += flow; - if (used == w) break; //残余流量用尽了就停止增广 - } - } - } - return used; -} -int main() { - scanf("%d%d%d%d", &n, &m, &s, &t); - for (int i = 1; i <= m; i++) { - int u, v, w; - scanf("%d%d%d", &u, &v, &w); - addedge(u, v, w); - addedge(v, u, 0); - } - int ans = 0; - while (bfs()) ans += dfs(s, INF); - printf("%d\n", ans); - return 0; -} -``` - -### ISAP - -这个是 SAP 算法的加强版 (Improved)。 - -## PR 预留推进算法 - -该方法在求解过程中忽略流守恒性，并每次对一个结点更新信息，以求解最大流 - -有 HLPP 的主流算法 - -推送 - 重贴标签算法通过对单个结点的更新操作，直到没有结点需要更新来求解最大流 - -算法过程维护的流函数不一定保持流守恒性，对于一个结点，我们允许进入结点的流超过流出结点的流，超过的部分被称为结点 $u(u\in V-\{s,t\})$ 的**超额流**$e(u)$： - -$$ -e(u)=\sum_{(x,u)\in E}f(x,u)-\sum_{(u,y)\in E}f(u,y) -$$ - -若 $e(u)>0$，称结点 $u$**溢出**. - -推送 - 重贴标签算法维护每个结点的高度 $h(u)$，并且规定溢出的结点 $u$ 如果要推送超额流，只能向高度小于 $u$ 的结点推送；如果 $u$ 没有相邻的高度小于 $u$ 的结点，就修改 $u$ 的高度（重贴标签）。 - -#### 高度函数 - -准确地说，推送 - 重贴标签维护以下的一个映射 $h:V\to \mathbf{N}$： - -- $h(s)=|V|,h(t)=0$ -- $\forall (u,v)\in E_f,h(u)\leq h(v)+1$ - -则称 $h$ 是残存网络 $G_f=(V_f,E_f)$ 的高度函数。 - -引理 1：设 $G_f$ 上的高度函数为 $h$，对于任意两个结点 $u,v\in V$，如果 $h(u)>h(v)+1$，则 $(u,v)$ 不是 $G_f$ 中的边。 - -算法只会在 $h(u)=h(v)+1$ 的边执行推送。 - -#### 推送 -Push - -适用条件：结点 $u$ 溢出，且存在结点 $v((u,v)\in V_f,c(u,v)-f(u,v)>0,h(u)=h(v)+1)$，则 push 操作适用于 $(u,v)$。 - -于是，我们尽可能将超额流从 $u$ 推送到 $v$，推送过程中我们只关心超额流和 $c(u,v)-f(u,v)$ 的最小值，不关心 $v$ 是否溢出。 - -如果 $(u,v)$ 在推送完之后满流，将其从残存网络中删除。 - -#### 重贴标签 -Relabel - -适用条件：如果结点 $u$ 溢出，且 $\forall (u,v)\in E_f,h(u)\leq h(v)$，则 relabel 操作适用于 $u$。 - -则将 $h(u)$ 更新为 $min_{(u,v)\in E_f}h(v)+1$ 即可。 - -#### 初始化 - -$$ -\begin{split} -&\forall (u,v)\in E,&f(u,v)=\left\{\begin{split} -&c(u,v)&,u=s\\ -&0&,u\neq s\\ -\end{split}\right.. -\\ -&\forall u\in V,&h(u)=\left\{\begin{split} -&|V|&,u=s\\ -&0&,u\neq s\\ -\end{split}\right. -,e(u)=\sum_{(x,u)\in E}f(x,u)-\sum_{(u,y)\in E}f(u,y). -\end{split}. -$$ - -上述将 $(s,v)\in E$ 充满流，并将 $h(s)$ 抬高，使得 $(s,v)\notin E_f$，因为 $h(s)>h(v)$，而且 $(s,v)$ 毕竟满流，没必要留在残存网络中；上述还将 $e(s)$ 初始化为 $\sum_{(s,v)\in E}f(s,v)$ 的相反数。 - -#### 通用执行框架 - -无需按照特定顺序，执行以下过程： - -- 只要存在结点 $u$ 满足 push 或 relabel 的条件，就执行对应的操作。 - -如图，每个结点中间表示编号，左下表示高度值 $h(u)$，右下表示超额流 $e(u)$，结点颜色的深度也表示结点的高度；边权表示 $c(u,v)-f(u,v)$，绿色的边表示满足 $h(u)=h(v)+1$ 的边 $(u,v)$（即残存网络的边 $E_f$）： - -![p1](./images/2148.png) - -整个算法我们大致浏览一下过程，这里笔者使用的是一个暴力算法，即暴力扫描是否有溢出的结点，有就更新 - -![p2](./images/2149.gif) - -最后的结果 - -![p3](./images/2150.png) - -可以发现，最后的超额流一部分回到了 $s$，且除了源点汇点，其他结点都没有溢出；这时的流函数 $f$ 满足流守恒性，为最大流，即 $e(t)$。 +## 网络流的常见问题 -#### 核心代码 +网络流问题中常见的有以下三种：最大流，最小割，费用流。 -```cpp -const int N=1e4+4,M=1e5+5,INF=0x3f3f3f3f; -int n,m,s,t,maxflow,tot; -int ht[N],ex[N]; -void init(){// 初始化 - for(int i=h[s];i;i=e[i].nex){const int &v=e[i].t; - ex[v]=e[i].v,ex[s]-=ex[v],e[i^1].v=e[i].v,e[i].v=0; - } - ht[s]=n; -} -bool push(int ed){const int &u=e[ed^1].t,&v=e[ed].t; - int flow=min(ex[u],e[ed].v); - ex[u]-=flow,ex[v]+=flow,e[ed].v-=flow,e[ed^1].v+=flow; - return ex[u];// 如果 u 仍溢出，返回 1 -} -void relabel(int u){ - ht[u]=INF; - for(int i=h[u];i;i=e[i].nex)if(e[i].v)ht[u]=min(ht[u],ht[e[i].t]); - ++ht[u]; -} -``` +### 最大流 -### HLPP 算法 - -最高标号预流推进算法（High Level Preflow Push）是基于推送 - 重贴标签算法的优先队列实现，该算法优先推送高度高的溢出的结点，算法算法复杂度 $O(n^2\sqrt m)$.=。 - -具体地说，HLPP 维护以下过程： - -1. 初始化（基于推送 - 重贴标签算法）； -2. 选择溢出结点（除 $s,t$）中高度最高的结点 $u$，并对它所有可以推送的边进行推送； -3. 如果 $u$ 仍溢出，对它重贴标签，回到 2； -4. 如果没有溢出的结点，算法结束。 - -#### BFS 优化 - -HLPP 的上界为 $O(n^2\sqrt m)$，但在使用时卡得比较紧；我们可以在初始化高度的时候进行优化： - -具体来说，我们初始化 $h(u)$ 为 $u$ 到 $t$ 的最短距离；特别地，$h(s)=n$。 - -在 BFS 的同时我们顺便检查图的联通性，排除无解的情况。 - -#### GAP 优化 - -HLPP 推送的条件是 $h(u)=h(v)+1$，而如果在算法的某一时刻，$h(u)=t$ 的结点个数为 $0$，那么对于 $h(u)>t$ 的结点就永远无法推送超额流到 $t$，因此只能送回 $s$，那么我们就在这时直接让他们的高度变成 $n+1$，以尽快推送回 $s$，减少重贴标签的操作。 - -#### [LuoguP4722] 【模板】最大流 加强版 / 预流推进 +我们有一张图，要求从源点流向汇点的最大流量（可以有很多条路到达汇点），就是我们的最大流问题。 -```cpp -#include -#include -#include -using namespace std; -const int N=1e4+4,M=2e5+5,INF=0x3f3f3f3f; -int n,m,s,t; +### 最小费用最大流 -struct qxx{int nex,t,v;}; -qxx e[M*2]; -int h[N],cnt=1; -void add_path(int f,int t,int v){e[++cnt]=(qxx){h[f],t,v},h[f]=cnt;} -void add_flow(int f,int t,int v){add_path(f,t,v);add_path(t,f,0);} +最小费用最大流问题是这样的：每条边都有一个费用，代表单位流量流过这条边的开销。我们要在求出最大流的同时，要求花费的费用最小。 -int ht[N],ex[N],gap[N];// 高度；超额流；gap 优化 -bool bfs_init(){ - memset(ht,0x3f,sizeof(ht)); - queue q; - q.push(t),ht[t]=0; - while(q.size()){// 反向 BFS, 遇到没有访问过的结点就入队 - int u=q.front();q.pop(); - for(int i=h[u];i;i=e[i].nex){const int &v=e[i].t; - if(e[i^1].v&&ht[v]>ht[u]+1)ht[v]=ht[u]+1,q.push(v); - } - } - return ht[s]!=INF;// 如果图不连通，返回 0 -} -struct cmp{bool operator()(int a,int b)const{return ht[a],cmp> pq;// 将需要推送的结点以高度高的优先 -bool vis[N];// 是否在优先队列中 -int push(int u){// 尽可能通过能够推送的边推送超额流 - for(int i=h[u];i;i=e[i].nex){const int &v=e[i].t,&w=e[i].v; - if(!w||ht[u]!=ht[v]+1)continue; - int k=min(w,ex[u]);// 取到剩余容量和超额流的最小值 - ex[u]-=k,ex[v]+=k,e[i].v-=k,e[i^1].v+=k;//push - if(v!=s&&v!=t&&!vis[v])pq.push(v),vis[v]=1;// 推送之后，v 必然溢出，则入堆，等待被推送 - if(!ex[u])return 0;// 如果已经推送完就返回 - }return 1; -} -void relabel(int u){// 重贴标签（高度） - ht[u]=INF; - for(int i=h[u];i;i=e[i].nex)if(e[i].v)ht[u]=min(ht[u],ht[e[i].t]); - ++ht[u]; -} -int hlpp(){// 返回最大流 - if(!bfs_init())return 0;// 图不连通 - ht[s]=n; - memset(gap,0,sizeof(gap)); - for(int i=1;i<=n;i++)if(ht[i]!=INF)gap[ht[i]]++;// 初始化 gap - for(int i=h[s];i;i=e[i].nex){const int v=e[i].t,w=e[i].v;// 队列初始化 - if(!w)continue; - ex[s]-=w,ex[v]+=w,e[i].v-=w,e[i^1].v+=w;// 注意取消 w 的引用 - if(v!=s&&v!=t&&!vis[v])pq.push(v),vis[v]=1;// 入队 - } - while(pq.size()){ - int u=pq.top();pq.pop(),vis[u]=0; - while(push(u)){// 仍然溢出 - // 如果 u 结点原来所在的高度没有结点了，相当于出现断层 - if(!--gap[ht[u]])for(int i=1;i<=n;i++) - if(i!=s&&i!=t&&ht[i]>ht[u]&&ht[i]