Loading... # 1 引言 在现实生活中，我们可能会遇到这样一些问题：1. 有一堆任务待分配。对不同任务，不同员工所需费用不同。现在要求找到一种分配方式使得总费用最小。2. 给定一个网络，网络中每条边显示两地间潜在货运量，现在要求计算网络能支持的最大流量。 3.确定从工厂向零售商店运输商品性价比最高的方法。实际中，类似于这样的一些问题有很多，在数据量较小的情况下，我们可以使用遍历的方法解决。可如果，出现需要大量计算的问题，我们又该如何解决呢？ 其实，上面所讲述的这些问题都可以转化为网络流问题，因为他们都是从一或多个源点流向一或多个汇点。网络流问题中常见的有以下三种：最大流，最小割，费用流。 **最大流**：我们有一张图，要求从源点流向汇点的最大流量（可以有很多条路到达汇点），就是我们的最大流问题。 **最小费用最大流**：每条边都有一个费用，代表单位流量流过这条边的开销。我们要在求出最大流的同时，要求花费的费用最小。 **最小割**：割其实就是删边的意思，当然最小割就是割掉 $X$ 条边来让 $S$ 跟 $T$ 不互通。我们要求 $X$ 条边加起来的流量总和最小，这就是最小割问题。 本篇论文所研究的是网络流的最大流问题以及上下界网络流中的无源汇上下界可行流问题。 # 2 算法原理 ## 2.1 网络与网络流 首先我们介绍网络与网络流的基本概念。 网络是指一个有向图$G = (V,E)$。每条边 $(u,v)\in E$ 都有一个权值 $c(u,v)$，称之为容量（Capacity），当 $(u,v)\notin E$ 时有 $c(u,v)=0$ 。其中有两个特殊的点：源点（Source）$s\in V$ 和汇点（Sink）$t\in V,(s\neq t)$。 设 $f(u,v)$ 定义在二元组 $(u\in V,v\in V)$ 上的实数函数且满足 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)$ 那么 $f$ 称为网络 $G$ 的流函数。对于 称为边的 $(u,v)\in E$ , **流量**，$c(u,v)-f(u,v)$ 称为边的 **剩余容量**。整个网络的流量为 $\sum_{(s,v)\in E}f(s,v)$，即 **从源点发出的所有流量之和**。 一般而言也可以把网络流理解为整个图的流量。而这个流量必满足上述三个性质。 注：流函数的完整定义为 $$ f(u,v)=\left\lbrace \begin{aligned} &f(u,v),&(u,v)\in E\newline &-f(v,u),&(v,u)\in E\newline &0,&(u,v)\notin E,(v,u) \notin E \end{aligned} \right. $$ ## 2.2 残量网络 我们介绍一下一条边的剩余容量 $c_f(u,v)$（Residual Capacity），它表示的是这条边的容量与流量之差，即 $c_f(u,v)=c(u,v)-f(u,v)$。 对于流函数 $f$，残存网络 $G_f$（Residual Network）是网络 $G$ 中所有结点 和剩余容量大于 0的边构成的子图。形式化的定义，即 $G_f=(V_f=V,E_f=\{(u,v)\in E,c_f(u,v)>0\})$。 注意，剩余容量大于 0 的边可能不在原图 $G$ 中（根据容量、剩余容量的定义以及流函数的斜对称性得到）。可以理解为，残量网络中包括了那些还剩了流量空间的边构成的图，也包括虚边（即反向边）。 ## 2.3 增广路 在原图 $G$ 中若一条从源点到汇点的路径上所有边的 **剩余容量都大于 0**，这条路被称为增广路（Augmenting Path）。或者说，在残存网络 $G_f$ 中，一条从源点到汇点的路径被称为增广路。如图：  ## 2.4 Dinic 算法 求解网络流的的算法现在有很多种，本篇论文主要介绍的是Dinic 算法。 Dinic 算法的过程是这样的：每次增广前，我们先用 BFS 来将图分层。设源点的层数为 $0$，那么一个点的层数便是它离源点的最近距离。 通过分层，我们可以干两件事情： 1. 如果不存在到汇点的增广路（即汇点的层数不存在），我们即可停止增广。 2. 确保我们找到的增广路是最短的。（原因见下文） 接下来是 DFS 找增广路的过程。 我们每次找增广路的时候，都只找比当前点层数多 $1$ 的点进行增广（这样就可以确保我们找到的增广路是最短的）。 Dinic 算法有两个优化： 1. **多路增广**：每次找到一条增广路的时候，如果残余流量没有用完怎么办呢？我们可以利用残余部分流量，再找出一条增广路。这样就可以在一次 DFS 中找出多条增广路，大大提高了算法的效率。 2. **当前弧优化**：如果一条边已经被增广过，那么它就没有可能被增广第二次。那么，我们下一次进行增广的时候，就可以不必再走那些已经被增广过的边。 ## 2.5 基于 Dinic 算法实现的无源汇上下界可行流 给定无源汇流量网络 $G$。询问是否存在一种标定每条边流量的方式，使得每条边流量满足上下界同时每一个点流量平衡。 不妨假设每条边已经流了 $b(u,v)$ 的流量，设其为初始流。同时我们在新图中加入 $u$ 连向 $v$ 的流量为 $c(u,v) - b(u,v)$ 的边。考虑在新图上进行调整。 由于最大流需要满足初始流量平衡条件（最大流可以看成是下界为 $0$ 的上下界最大流），但是构造出来的初始流很有可能不满足初始流量平衡。假设一个点初始流入流量减初始流出流量为 $M$。 若 $M=0$，此时流量平衡，不需要附加边。 若 $M>0$，此时入流量过大，需要新建附加源点 $S'$，$S'$ 向其连流量为 $M$ 的附加边。 若 $M<0$，此时出流量过大，需要新建附加汇点 $T'$，其向 $T'$ 连流量为 $-M$ 的附加边。 如果附加边满流，说明这一个点的流量平衡条件可以满足，否则这个点的流量平衡条件不满足。（因为原图加上附加流之后才会满足原图中的流量平衡。） 所以我们只需要在建图完毕之后，使用 Dinic 算法跑 $S'$ 到 $T'$ 的最大流，若 $S'$ 连出去的边全部满流，则存在可行流，否则不存在。 ## 2.6 图的基本存储结构 我们所使用的算法基本储存结构为邻接表与链式前向星，图的结构用C++描述如下： ```C++ struct Edge { int u, v, down, up; } E[maxm]; struct edge { int to, nxt, flow; } e[(maxn + maxm) << 1]; ``` 在代码中，$Edge$ 采取邻接表的储存结构，储存内容是我们所构建的图，其中$u$，$v$为顶点，$down$，$up$为下确界与上确界；$edge$ 采取链式前向星的储存结构，储存内容是我们所构建的附加边。 ## 2.6 图的构建 为解决无源汇上下界可行流问题，Dinic 算法需要为图添加两种边：初始边与附加边。 初始边为我们所构建的图的边，附加边是为了满足初始流量平衡而构造出的边。 初始边的添加通过直接为邻接表赋值，代码如下： ```C++ E[num].u = u, E[num].v = v, E[num].down = low, E[num].up = up; ``` 附加边的添加通过函数，代码如下： ```c++ inline void add(int u, int v, int w) { e[++cnt_edge].nxt = head[u]; head[u] = cnt_edge; e[cnt_edge].to = v; e[cnt_edge].flow = w; } inline void addflow(int u, int v, int w) { add(u, v, w); add(v, u, 0); } ``` ## 2.7 寻找增广路 在Dinic算法中，每次增广我们先使用BFS来将图分层，如果不存在到汇点的增广路（即汇点的层数不存在），我们即可停止增广。代码如下： ```c++ inline bool bfs() { memset(dep, 0, sizeof(dep)); for (int i = S; i <= T; i++) cur[i] = head[i]; queue<int> q; q.push(S); dep[S] = 1; while (!q.empty()) { int x = q.front(); q.pop(); for (int i = head[x]; i; i = e[i].nxt) { int y = e[i].to; if (dep[y] || e[i].flow <= 0) continue; dep[y] = dep[x] + 1; q.push(y); } } return dep[T] > 0; } ``` 接下来是 DFS 找增广路的过程，我们每次找增广路的时候，都只找比当前点层数多 $1$ 的点进行增广。代码如下： ```c++ int dfs(int x, int lim) { if (x == T || lim <= 0) return lim; int res = lim; for (int i = cur[x]; i; i = e[i].nxt) { cur[x] = i; int y = e[i].to; if (dep[y] != dep[x] + 1 || e[i].flow <= 0) continue; int tmp = dfs(y, min(res, e[i].flow)); if (tmp <= 0) dep[y] = 0; res -= tmp; e[i].flow -= tmp, e[i ^ 1].flow += tmp; if (res <= 0) break; } return lim - res; } ``` ## 2.8 Dinic 算法的正确性和时间复杂性分析 设点数为 $n$，边数为 $m$，那么 Dinic 算法的时间复杂度（在应用上面两个优化的前提下）是 $O(n^{2}m)$。 首先考虑单轮增广的过程。在应用了 **当前弧优化** 的前提下，对于每个点，我们维护下一条可以增广的边，而当前弧最多变化 $m$ 次，从而单轮增广的最坏时间复杂度为 $O(nm)$。 接下来我们证明，最多只需 $n-1$ 轮增广即可得到最大流。 我们先回顾下 Dinic 的增广过程。对于每个点，Dinic 只会找比该点层数多 $1$ 的点进行增广。 首先容易发现，对于图上的每个点，一轮增广后其层数一定不会减小。而对于汇点 $t$，情况会特殊一些，其层数在一轮增广后一定增大。 对于后者，我们考虑用反证法证明。如果 $t$ 的层数在一轮增广后不变，则意味着在上一次增广中，仍然存在着一条从 $s$ 到 $t$ 的增广路，且该增广路上相邻两点间的层数差为 $1$。这条增广路应该在上一次增广过程中就被增广了，这就出现了矛盾。 从而我们证明了汇点的层数在一轮增广后一定增大，即增广过程最多进行 $n-1$ 次。 综上 Dinic 的最坏时间复杂度为 $O(n^{2}m)$。事实上在一般的网络上，Dinic 算法往往达不到这个上界。 特别地，在求解二分图最大匹配问题时，Dinic 算法的时间复杂度是 $O(m\sqrt{n})$。接下来我们将给出证明。 首先我们来简单归纳下求解二分图最大匹配问题时，建立的网络的特点。我们发现这个网络中，所有边的流量均为 $1$，且除了源点和汇点外的所有点，都满足入边最多只有一条，或出边最多只有一条。我们称这样的网络为 **单位网络**。 对于单位网络，一轮增广的时间复杂度为 $O(m)$，因为每条边只会被考虑最多一次。 接下来我们试着求出增广轮数的上界。假设我们已经先完成了前 $\sqrt{n}$ 轮增广，因为汇点的层数在每次增广后均严格增加，因此所有长度不超过 $\sqrt{n}$ 的增广路都已经在之前的增广过程中被增广。设前 $\sqrt{n}$ 轮增广后，网络的流量为 $f$，而整个网络的最大流为 $f'$，设两者间的差值 $d=f'-f$。 因为网络上所有边的流量均为 $1$，所以我们还需要找到 $d$ 条增广路才能找到网络最大流。又因为单位网络的特点，这些增广路不会在源点和汇点以外的点相交。因此这些增广路至少经过了 $d\sqrt{n}$ 个点（每条增广路的长度至少为 $\sqrt{n}$），且不能超过 $n$ 个点。因此残量网络上最多还存在 $\sqrt{n}$ 条增广路。也即最多还需增广 $\sqrt{n}$ 轮。 综上，对于包含二分图最大匹配在内的单位网络，Dinic 算法可以在 $O(m\sqrt{n})$ 的时间内求出其最大流。 # 3 实验 ## 3.1 实验题目 网络流-矩阵计算 C时间限制：3000 毫秒 | C内存限制：3000 Kb **题目内容：** 有一个n行m列的整数矩阵A， 知道每行的和以及每列的和，还知道一些矩阵元素的约束如A[i][j]<x, 或者A[i][j]>y等， 判断该是否存在满足上述条件的可行矩阵。 **输入描述：** 第一行是测试用例的数目c 每个测试用例的第一行是n,m 表示行和列 接下来一行是n个行和 接下来一行是m个列和 然后一行是约束个数k 接下来k行是约束，每个约束如 a b c d, 其中a,b是某个元素的行列坐标，c是一个字符（>,=,<）, d是一个整数， 2 3 > 4 表示的意思是A[2][3]>4。矩阵左上角坐标规定为（1,1），所以一个约束的a为0，则表示b列所有的元素， 而如果b为0，则表示a行所有的元素。 **输出描述：** 如果存在，则输出这个矩阵；否则输出“不存在” **输入样例：** 2 2 3 8 10 5 6 7 4 0 2 > 2 2 1 = 3 2 3 > 2 2 3 < 5 2 2 4 5 6 7 1 1 1 > 10 **输出样例：** 2 3 3 3 3 4 不可能 ## 3.2 问题分析 这一道题是网络流中的**无源汇上下界可行流问题**，我们可以通过Dinic算法来解决。 首先问题给出一个矩形的长和宽，同时要求矩形的每一行行、列之和分别为一个定值；之后给出 $k$ 个约束，代表某一个矩阵的值必须满足这个约束条件（大于、小于或等于一个值）。 那么我们建图就设置两列，一列表示行，一列表示列，连边代表限制，行与列之间的连边代表点的限制。他们是一个矩阵，本身就是有源有汇的，我们可以使用 $x_i$ 表示行，$y_i$ 表示列；我们再新建两个点 $s$ 和 $t$ ，建立 $s$ 连向行 $x_i$ 、 $t$ 连向列 $y_i$ 的边；而他们被赋予一个定值，所以他们的上下界便设置为这个定值。之后 $k$ 个约束，我们将 $x_i$ 向 $y_i$ 根据题目意思建立边，并且添加上界与下界，然后再建立附加源汇 $ss$ 和 $tt$。 这样，这个问题就被转化为**无源汇上下界可行流问题**，我们求出行 $x_i$ 向列 $y_i$ 的可行流，便解决了这个问题。 ## 3.3 核心算法描述 构图思路如3.2所示，我们构建出这样一个图：  接下来我们只需要为他添加附加边，使其满足初始流量平衡，即可解决问题。 ## 3.4 测试数据 一共给出6个测试数据 (1) 2 (2) 128 (3) 200 (4) 256 (5) 512 (6) 9216 ## 3.5 性能测试 | | 测试数据1 | 测试数据2 | 测试数据3 | 测试数据4 | 测试数据5 | 测试数据6 | | -------------------- | ----------- | ----------- | ----------- | ----------- | ----------- | ----------- | | 数据规模 | N = 2 | N = 128 | N = 200 | N = 256 | N = 512 | N = 9216 | | 算法运行时间（ms） | 603 | 2954 | 3482 | 4654 | 13427 | 122815 | 可以看出，算法在数据量较小时，运行时间增长较缓慢。而数据量达到很高时，算法运行时间增长较快，但仍然可以看出该算法解决无源汇上下界可行流问题的迅速。 # 4 小结 Dinic 算法的最坏时间复杂度是 $O(n^{2}m)$，而解决二分图最大匹配问题时，Dinic 算法的时间复杂度是 $O(m\sqrt{n})$。 Dinic 算法可以通过添加附加边的方式解决无源汇上下界可行流问题，并且在这一类问题的时间复杂性上，拥有优越性。 # 5 参考文献 [1] 徐义春,万书振,解德祥. 算法设计与分析[M]. 北京：清华大学出版社，2016. [2] Cherkassky B V, Goldberg A V. On implementing push-relabel method for the maximum flow problem[C]. International Conference on Integer Programming and Combinatorial Optimization. Springer, Berlin, Heidelberg, 1995: 157-171. [3] Ahuja R K, Kodialam M, Mishra A K, et al. Computational investigations of maximum flow algorithms[J]. European Journal of Operational Research, 1997, 97(3): 509-542. [4] Derigs U, Meier W. Implementing Goldberg's max-flow-algorithm—A computational investigation[J]. Zeitschrift für Operations Research, 1989, 33(6): 383-403. [5] Even, Shimon; Tarjan, R. Endre. Network Flow and Testing Graph Connectivity[J]. SIAM Journal on Computing,1975, 4 (4): 507–518. # 附录1：完整代码 ```c++ #include <string.h> #include <algorithm> #include <cstdio> #include <iostream> #include <queue> #include <vector> using namespace std; const int maxn = 2010; const int maxm = 10210; const int inf = 1e9; int epoch; int n, num, cnt_edge = 1, sum, S, T; int x, y, k; int head[maxn], cur[maxn], in[maxn], out[maxn], dep[maxn]; int low_[maxn][maxn]; int up_[maxn][maxn]; struct Edge { int u, v, down, up; } E[maxm]; struct edge { int to, nxt, flow; } e[(maxn + maxm) << 1]; inline void add(int u, int v, int w) { e[++cnt_edge].nxt = head[u]; head[u] = cnt_edge; e[cnt_edge].to = v; e[cnt_edge].flow = w; } inline void addflow(int u, int v, int w) { add(u, v, w); add(v, u, 0); } inline bool bfs() { memset(dep, 0, sizeof(dep)); for (int i = S; i <= T; i++) cur[i] = head[i]; queue<int> q; q.push(S); dep[S] = 1; while (!q.empty()) { int x = q.front(); q.pop(); for (int i = head[x]; i; i = e[i].nxt) { int y = e[i].to; if (dep[y] || e[i].flow <= 0) continue; dep[y] = dep[x] + 1; q.push(y); } } return dep[T] > 0; } int dfs(int x, int lim) { if (x == T || lim <= 0) return lim; int res = lim; for (int i = cur[x]; i; i = e[i].nxt) { cur[x] = i; int y = e[i].to; if (dep[y] != dep[x] + 1 || e[i].flow <= 0) continue; int tmp = dfs(y, min(res, e[i].flow)); if (tmp <= 0) dep[y] = 0; res -= tmp; e[i].flow -= tmp, e[i ^ 1].flow += tmp; if (res <= 0) break; } return lim - res; } inline int Dinic() { int res = 0; while (bfs()) res += dfs(S, inf); return res; } int build(int x, int y) { for (int i = 1; i <= x; i++) for (int j = 1; j <= y; j++) if (low_[i][j] > up_[i][j]) return 0; else { num++; int u = i, v = j + x; int low = low_[i][j], up = up_[i][j]; E[num].u = u, E[num].v = v, E[num].down = low, E[num].up = up; } return 1; } void init() { n = 0, num = 0, cnt_edge = 1, sum = 0; memset(head, 0, sizeof(head)); memset(cur, 0, sizeof(cur)); memset(in, 0, sizeof(in)); memset(out, 0, sizeof(out)); memset(dep, 0, sizeof(dep)); for (int i = 0; i < maxn; i++) { for (int j = 0; j < maxn; j++) { low_[i][j] = 0; up_[i][j] = inf; } } } int main() { scanf("%d", &epoch); while (epoch--) { init(); scanf("%d%d", &x, &y); n = x + y + 2; S = 0, T = n + 1; int s = n - 1, t = n; int sum1 = 0, sum2 = 0; for (int i = 1; i <= x; i++) { num++; int u = s, v = i; int up, low; int flow; scanf("%d", &flow); up = low = flow; sum1 += flow; E[num].u = u, E[num].v = v, E[num].down = low, E[num].up = up; } for (int i = 1; i <= y; i++) { num++; int u = x + i, v = t; int up, low; int flow; scanf("%d", &flow); up = low = flow; sum2 += flow; E[num].u = u, E[num].v = v, E[num].down = low, E[num].up = up; } scanf("%d", &k); for (int i = 1; i <= k; i++) { int u, v; int up, low; char c; int flow; scanf("%d%d", &u, &v); cin >> c; scanf("%d", &flow); int s1, t1, s2, t2; s1 = t1 = u, s2 = t2 = v; if (u == 0) s1 = 1, t1 = x; if (v == 0) s2 = 1, t2 = y; for (int i = s1; i <= t1; i++) { for (int j = s2; j <= t2; j++) { if (c == '>') { low_[i][j] = max(flow + 1, low_[i][j]); } else if (c == '<') { up_[i][j] = min(flow - 1, up_[i][j]); } else { low_[i][j] = max(flow, low_[i][j]); up_[i][j] = min(flow, up_[i][j]); } } } } if (!build(x, y)) { printf("不可能\n"); continue; } for (int i = 1; i <= num; i++) in[E[i].v] += E[i].down, out[E[i].u] += E[i].down; for (int i = 1; i <= num; i++) addflow(E[i].u, E[i].v, E[i].up - E[i].down); for (int i = 1; i <= n - 2; i++) { if (in[i] >= out[i]) addflow(S, i, in[i] - out[i]), sum += in[i] - out[i]; else addflow(i, T, out[i] - in[i]); } if (Dinic() != sum|| sum1 != sum2) { printf("不可能\n"); continue; }; for (int i = 2 * (x + y + 1); i <= 2 * num + 1; i += 2) { if ((i - 2 * (x + y)) % (2 * y) == 0) printf("%d\n", E[i / 2].down + e[i ^ 1].flow); else printf("%d ", E[i / 2].down + e[i ^ 1].flow); } } return 0; } ``` 最后修改：2021 年 12 月 04 日 04 : 07 PM © 允许规范转载