最大流之Ford-Fulkerson方法详解及实现

最大流问题常常出现在物流配送中，可以规约为以下的图问题。最大流问题中，图中两个顶点之间不能同时存在一对相反方向的边。

边上的数字为该条边的容量，即在该条边上流过的量的上限值。最大流问题就是在满足容量限制条件下，使从起点s到终点t的流量达到最大。在介绍解决最大流问题的Ford-Fulkerson方法之前，先介绍一些基本概念。

1.  残存网络与增广路径

根据图和各条边上的流可以画出一幅图的残存网络如下所示。左图为流网络，右图为残存网络，其中流网络中边上的数字分别是流量和容量，如10/12，那么10为边上的流量，12为边的容量。残存网络中可能会存在一对相反方向的边，与流网络中相同的边代表的是流网络中该边的剩余容量，在流网络中不存在的边代表的则是其在流网络中反向边的已有流量，这部分流量可以通过“回流”减少。例如，右图残存网络中，边<s,v1>的剩余容量为4，其反向边<v1.s>的值为12，即左图流网络中的边<s,v1>的流量。在残存网络中，值为0的边不会画出，如边<v1,v2>。

残存网络描述了图中各边的剩余容量以及可以通过“回流”删除的流量大小。在Ford-Fulkerson方法中，正是通过在残存网络中寻找一条从s到t的增广路径，并对应这条路径上的各边对流网络中的各边的流进行修改。如果路径上的一条边存在于流网络中，那么对该边的流增加，否则对其反向边的流减少。增加或减少的值是确定的，就是该增广路径上值最小的边。

2. Ford-Fulkerson方法

Ford-Fulkerson方法的正确性依赖于这个定理：当残存网络中不存在一条从s到t的增广路径，那么该图已经达到最大流。这个定理的证明及一些与其等同的定理可以参考《算法导论》。

Ford-Fulkerson方法的伪代码如下。其中<u,v>代表顶点u到顶点v的一条边，<u,v>.f表示该边的流量，c是边容量矩阵，c(i,j)表示边<i,j>的容量，当边<i,j>不存在时，c(i,j)=0。e为残存网络矩阵，e(i,j)表示边<i,j>的值，当边<i,j>不存在时，e(i,j)=0。E表示边的集合。f表示流网络。

[cpp] view plain copy

1. Ford-Fulkerson  
2.   
3.     for <u,v> ∈ E  
4.   
5.         <u,v>.f = 0  
6.   
7.     while find a route from s to t in e  
8.   
9.         m = min(<u,v>.f, <u,v>  ∈ route)  
10.   
11.         for <u,v> ∈ route  
12.   
13.             if <u,v>  ∈ f  
14.   
15.                 <u,v>.f = <u,v>.f + m  
16.   
17.             else  
18.   
19.                 <v,u>.f = <v,u>.f - m  

Ford-Fulkerson方法首先对图中的所有边的流量初始化为零值，然后开始进入循环：如果在残存网络中可以找到一条从s到t的增广路径，那么要找到这条这条路径上值最小的边，然后根据该值来更新流网络。

Ford-Fulkerson有很多种实现，主要不同点在于如何寻找增广路径。最开始提出该方法的Ford和Fulkerson同学在其论文中都是使用广度优先搜索实现的，其时间复杂度为O(VE)，整个算法的时间复杂度为O(VE^2)。

下面我给出一个应用Bellman-Ford计算单源最短路径的算法实现寻找一条增广路径，对于用邻接矩阵表示的图来说，该实现的时间复杂度为O(V^3)，对于用邻接表表示的图来说，时间复杂度则为O(VE)。

[cpp] view plain copy

1. // 寻找增广路径  
2. int findRoute(int **e, int vertexNum, int *priorMatrix, int s,int t)  
3. {  
4.     s--; t--;  
5.     int *d = (int *)malloc(sizeof(int)*vertexNum);  
6.     // initialize  
7.     for (int i = 0; i < vertexNum; i++)  
8.     {  
9.         d[i] = 0;  
10.         priorMatrix[i] = -1;  
11.     }  
12.     d[s] = 1;  
13.     // 反复用边<i,j>做松弛操作，将<s,...,j>更新为<s,...,i,j>  
14.     for (int k = 0; k < vertexNum; k++)  
15.     {  
16.         for (int i = 0; i < vertexNum; i++)  
17.         {  
18.             for (int j = 0; j < vertexNum; j++)  
19.             {  
20.                 if (d[j] == 0)  
21.                 {  
22.                     d[j] |= (d[i] & (*((int*)e + i*vertexNum + j) > 0));  
23.                     if (d[j] == 1)  
24.                     {  
25.                         priorMatrix[j] = i;  
26.                     }  
27.                 }  
28.             }  
29.         }  
30.     }  
31.     if (d[t] == 0)  return 0;  
32.   
33.     int min = INT_MAX;  
34.     int pre = priorMatrix[t];  
35.     while (pre != -1)  
36.     {  
37.         if (min > *((int*)e + pre*vertexNum + t))  
38.         {  
39.             min = *((int*)e + pre*vertexNum + t);  
40.         }  
41.         t = pre;  
42.         pre = priorMatrix[t];  
43.     }  
44.     return min;  
45. }  

该实现应用了计算图的最短路径方法中的思想，对图中的边反复在松弛操作，从而计算得到一个源点到其它所有点的路径。这里我们不需要计算最短路径，只要找到一条可行路径即可。上述findRoute方法的实现原理可以参考我前面的一篇文章 单源最短路径之Bellman-Ford算法 。在寻找路径的同时，我们还要记录一个前驱子图priorMatrix，其本质上是一个一位数组，其记录了从顶点s到其它顶点的一条可行路径上的终点的前一个顶点。于是我们就可以从前驱子图中找到从s到t的一条完整路径。其正确性由图的最短路径的计算方法思想保证。具体可以参考我另一篇博客 结点对最短路径之Floyd算法详解及实现 。

下面给出根据图和流网络计算残存网络的代码。

[cpp] view plain copy

1. // 计算残存网络  
2. void calculateENet(int **c, int vertexNum, int **f, int **e)  
3. {  
4.     for (int i = 0; i < vertexNum; i++)  
5.     {  
6.         for (int j = 0; j < vertexNum; j++)  
7.         {  
8.             int a = *((int*)c + i*vertexNum + j);  
9.             if (a != 0)  
10.             {  
11.                 *((int*)e + i*vertexNum + j) = a - *((int*)f + i*vertexNum + j);  
12.                 *((int*)e + j*vertexNum + i) = *((int*)f + i*vertexNum + j);  
13.             }  
14.             else  
15.             {  
16.                 *((int*)e + i*vertexNum + j) = 0;  
17.             }  
18.         }  
19.     }  
20. }  

下面给出Ford-Fulkerson方法的实现代码。

[cpp] view plain copy

1. /** 
2. * Ford-Fulkerson方法的一种实现 
3. * @param c 二维矩阵，记录每条边的容量 
4. * @param vertexNum 顶点个数，包括起点和终点 
5. * @param s 起点编号，编号从1开始 
6. * @param t 终点编号 
7. * @param f 输出流网络矩阵，二维矩阵，记录每条边的流量 
8. */  
9. void Ford_Fulkerson(int **c, int vertexNum, int s, int t, int **f)  
10. {  
11.     int *e = (int *)malloc(sizeof(int)*vertexNum*vertexNum);    // 残存网络  
12.     int *priorMatrix = (int *)malloc(sizeof(int)*vertexNum);    // 增广路径的前驱子图  
13.   
14.     // initialize  
15.     for (int i = 0; i < vertexNum;i++)  
16.     {  
17.         for (int j = 0; j < vertexNum; j++)  
18.         {  
19.             *(f + i*vertexNum + j) = 0;  
20.         }  
21.     }  
22.   
23.     while (1)  
24.     {  
25.         calculateENet(c, vertexNum, (int **)f, (int **)e);  // 计算残存网络  
26.         int min;  
27.         if ((min = findRoute((int **)e, vertexNum, priorMatrix, s, t)) == 0)    // 寻找增广路径及其最小流值  
28.         {  
29.             break;  
30.         }  
31.         int pre = priorMatrix[t - 1];  
32.         int next = t - 1;  
33.         while (pre != -1)       // 按增广路径更新流网络  
34.         {  
35.             if (*((int*)c + pre * vertexNum + next) != 0)  
36.             {  
37.                 *((int*)f + pre * vertexNum + next) += min;  
38.             }  
39.             else  
40.             {  
41.                 *((int*)f + next * vertexNum + pre) -= min;  
42.             }  
43.             next = pre;  
44.             pre = priorMatrix[pre];  
45.         }  
46.     }  
47. }  

3. 测试及效果

下面给出用于测试的代码。

[cpp] view plain copy

1. void testFord()  
2. {  
3.     int c[6][6] = { 0,      16,     13,     0,      0,      0,  
4.                     0,      0,      0,      12,     0,      0,  
5.                     0,      4,      0,      0,      14,     0,  
6.                     0,      0,      9,      0,      0,      20,  
7.                     0,      0,      0,      7,      0,      4,  
8.                     0,      0,      0,      0,      0,      0   };  
9.     int f[6][6];  
10.     Ford_Fulkerson((int **)c, 6, 1, 6, (int **)f);  
11.     for (int i = 0; i < 6; i++)  
12.     {  
13.         for (int j = 0; j < 6; j++)  
14.         {  
15.             int flow = f[i][j];  
16.             if (flow != 0)  
17.             {  
18.                 printf("%d -> %d : %d\n", i + 1, j + 1, flow);  
19.             }  
20.         }  
21.     }  
22. }  

上述代码构造的图如下所示。

运行结果如下，其中1为顶点s，5为顶点t，2~5依次为顶点v1、v2、v3和v4。

流网络和残存网络如下所示，其中左图为流网络，右图为残存网络。

我们可以看到残存网络中的确已经不存在一条从s到t的路径了。此时Ford-Fulkerson方法的循环应该终止，最大流量为各边的流量相加之和，即76。

完整的程序可以看到我的github项目 数据结构与算法

这个项目里面有本博客介绍过的和没有介绍的以及将要介绍的《算法导论》中部分主要的数据结构和算法的C实现，有兴趣的可以fork或者star一下哦~ 由于本人还在研究《算法导论》，所以这个项目还会持续更新哦~ 大家一起好好学习~