文章目录
- 0 前言
- 1 最短路径和矩阵乘法
- 2 Floyd-Warshall算法
- 3 用于稀疏图的Johnson算法
0 前言
如何找到一个图中所有结点之间的最短路径?
给定带权重的有向图G=(V,E),其权重函数为w:E→R,该函数将边映射到实数值的权重上。对于所有的结点对u,v∈V,希望找到一条从结点u到结点v的最短路径,其中,一条路径的权重为组成该路径的所有边的权重之和。
通过运行|V|次单源最短路径算法可以解决所有结点对之间的最短路径问题,只需要每一次使用一个不同的结点作为源结点即可。其中:
1、如果所有的边的权重为非负值,那么可以使用Dijkstra算法。如果使用线性数组来实现最小优先队列,那么该算法的运行时间将是O(V3+V·E)=O(V3)。
2、如果图中有权重为负值的边,可运行Bellman-Ford算法,同样地,只需要每次使用一个不同的结点作为源结点即可,算法的运行时间将是O(V2·E)。
本博文中的多数算法使用邻接矩阵来表示图,约定结点的编号为1,2,…,|V|,因此算法的输入将是一个n×n的矩阵W,该矩阵代表的是一个有n个结点的有向图G=(V,E)的边的权重,即W=(wij):
允许存在负权重的边,但图中不包括权重为负值的环路。
本博文讨论的所有结点对最短路径算法的表格输出也是一个n×n的矩阵D=(dij),其中dij代表的是从结点i到结点j的一条最短路径的权重。如果用δ(i,j)来表示从结点i到结点j之间的最短路径权重,则在算法终结时有dij=δ(i,j)。
为了解决所有结点对最短路径问题,不仅需要计算出最短路径权重,还需要计算出前驱结点矩阵Ⅱ=(πij),其中πij在i=j或从i到j不存在路径时为NIL,在其他情况下给出的是从结点i到结点j的某条最短路径上结点j的前驱结点。由矩阵Ⅱ的第i行所诱导的子图应当是一棵根结点为i的最短路径树。对于每个结点i∈V,定义图G对于结点i的前驱子图为Gπ,i=(Vπ,i,Eπ,i),其中:
如果Gπ,i是一棵最短路径树,则下面的过程将打印出从结点i到结点j的一条最短路径:
PRINT-ALL-PAIRS-SHORTEST-PATH(Ⅱ,i,j)
if i==j
print i
else if π[i][j]==NIL
print "no path from" i "to" j "exists"
else PRINT-ALL-PAIRS-SHORTEST-PATH(Ⅱ,i,π[i][j])
print j
下面将首先讨论如何用基于矩阵乘法的动态规划算法来解决所有结点对最短路径问题,之后会给出另一种动态规划算法,即Floyd-Warshall算法,同时还会讨论如何在有向图中找出传递闭包的问题,最后会讨论Johnson算法。
下面给出邻接矩阵表示的一些约定:
1、假定输入图G=(V,E)有n个结点,因此n=|V|。
2、使用大写字母来表示矩阵,如W、L或D,而用带下标的小写字母来表示矩阵中的元素,如wij、lij或dij。
3、一些矩阵将有带括号的上标,如L(m)=(lij(m))或D(m)=(dij(m))用来表示迭代。
4、对于一个给定的n×n矩阵A,假定矩阵的维度n存储在属性A.rows中。
1 最短路径和矩阵乘法
本节讨论有向图G=(V,E)上所有结点对最短路径问题的一种动态规划算法。在动态规划的每个大循环里,将调用一个与矩阵乘法非常相似的操作,因此该算法看上去就像是重复的矩阵乘法。设计动态规划算法的步骤如下:
1、分析最优解的结构。
2、递归定义最优解的值。
3、自底向上计算最优解的值。
先对最优解的结构开始进行分析。一条最短路径的所有子路径都是最短路径。假定用邻接矩阵来表示输入图,即W=(wij)。考虑从结点i到结点j的一条最短路径p,假定p至多包含m条边,还假定没有权重为负值的环路,且m为有限值。如果i=j,则p的权重为0且不包含任何边。如果结点i和结点j不同,则将路径p分解为:
其中路径p’至多包含m-1条边。根据引理:
可知p’是从结点i到结点k的一条最短路径,因此δ(i,j)=δ(i,k)+wkj。现在设lij(m)为从结点i到结点j的至多包含m条边的任意路径中的最小权重。当m=0时,从结点i到结点j之间存在一条没有边的最短路径当且仅当i=j。因此:
上图中m=0。对于m≥1,需要计算的lij(m)是【lij(m-1),即从i到j最多由m-1条边组成的最短路径的权重】和【从i到j最多由m条边组成的任意路径的最小权重】的最小值,通过对j的所有可能前驱k进行检查来获得该值。因此递归定义:
因为对于所有的j有wjj=0,所以上述式子中后面的等式成立(需要仔细想想)。但真正的最短路径权重δ(i,j)到底是多少呢?如果图G不包含权重为负值的环路,则对于每一对结点i和j,如果δ(i,j)<∞,则从i到j之间存在一条最短路径。由于该路径是简单路径,其包含的边最多为n-1条。从结点i到结点j的由多于n-1条边构成的路径不可能有比从i到j的最短路径权重更小的权重。因此真正的最短路径权重可以由下面的公式给出:
根据输入矩阵W=(wij),现在可以计算出矩阵序列L(1),L(2),…,L(n-1),其中对于m=1,2,…,n-1,有L(m)=(lij(m))。最后的矩阵L(n-1)包含的是最短路径的实际权重。注意,对于所有的结点i和j,L(1)=(wij),因此L(1)=W。
该算法的核心如下面的伪代码程序所示。该伪代码程序可以在给定W和L(m-1)的情况下,计算出L(m),即该伪代码将最近计算出的最短路径扩展了一条边:
EXTEND-SHORTEST-PATHS(L,W)
n=L.rows
let L'=(l'[i][j]) be a new n×n matrix
for i=1 to n
for j=1 to n
l'[i][j]=∞
for k=1 to n
l'[i][j]=min(l'[i][j],l[i][k]+w[k][j])
return L'
该过程计算在算法结束时返回的矩阵L’=(l’ij),计算该矩阵的方式是对于所有的i和j来计算:
使用L作为L(m-1),L’作为L(m)。这里在写法上没有注明上标的目的是让输入和输出矩阵独立于变量m。由于有3层嵌套的for循环,该算法的运行时间为Θ(n3)。
现在可以看到该算法与矩阵乘法的关系了。
假定希望计算矩阵乘积C=A×B,其中A和B均为n×n的矩阵。那么对于i,j=1,2,…,n,计算:
如果在:
中进行如下替换:
则将获得上上上一张图中cij的公式(你可以试试)。因此如果对算法EXTEND-SHORTEST-PATHS进行上述替换,并用0(加法操作+
的不变量)来替换∞(求最小值操作min的不变量),则获得的就是与下面的平方矩阵乘法同为时间Θ(n3)的矩阵乘法。
SQUARE-MATRIX-MULTIPLY(A,B)
n=A.rows
let C be a new n×n matrix
for i=1 to n
for j=1 to n
c[i][j]=0
for k=1 to n
c[i][j]+=a[i][k]*b[k][j]
return C
回到所有结点对最短路径问题,通过对最短路径一条边一条边地扩展来计算最短路径权重。设A·B表示由算法EXTEND-SHORTEST-PATHS(A,B)所返回的矩阵“乘积”,可以计算出下面由n-1个矩阵所构成的矩阵序列:
如上所述,矩阵L(n-1)=Wn-1包含的是最短路径权重。下面的伪代码程序在Θ(n4)时间内计算出该矩阵序列:
SLOW-ALL-PAIRS-SHORTEST-PATHS(W)
n=W.rows
L(1)=W
for m=2 to n-1
let L(m) be a new n×n matrix
L(m)=EXTEND-SHORTEST-PATHS(L(m-1),W)
return L(n-1)
下图(记为图1)所示为一个图和在图上运行SLOW-ALL-PAIRS-SHORTEST-PATHS所计算出的矩阵序列L(m):
算法的目标并不是要计算所有的L(m)矩阵,感兴趣的仅仅是矩阵L(n-1)。在没有权重为负值的环路的情况下,公式:
意味着对于所有的m≥n-1,有L(m)=L(n-1),可仅用⌈lg(n-1)⌉个矩阵乘积来计算矩阵L(n-1),计算的方法如下:
由于2⌈lg(n-1)⌉>n-1,因此最后的乘积L(2⌈lg(n-1)⌉)等于L(n-1)。下面的过程使用重复平方技术来计算上述矩阵序列:
FASTER-ALL-PAIRS-SHORTEST-PATHS(W)
n=W.rows
L(1)=W
m=1
while m<n-1
let L(2m) be a new n×n matrix
L(2m)=EXTEND-SHORTEST-PATHS(L(m),L(m))
m=2*m
return L(m)
在算法第5-8行的while循环的每一次迭代,计算L(2m)=(L(m))2,整个计算从m=1开始。在每次迭代的末尾,对m的取值进行加倍。最后的迭代通过实际计算L(2m)所计算出的是L(n-1),其中,n-1≤2m<2n-2。根据公式:
L(2m)=L(n-1)。下一次再执行算法第5行的测试时,m已经加倍,因此现在m>n-1,该测试将不通过,整个算法返回的是其最后所计算的矩阵。
因为⌈lg(n-1)⌉个矩阵中的每个矩阵的计算时间为Θ(n3),因此FASTER-ALL-PAIRS-SHORTEST-PATHS的运行时间为Θ(n3·lgn)。由于该代码非常紧凑,没有包含任何精巧的数据结构,隐藏在Θ记号中的常数应该较小。
2 Floyd-Warshall算法
本节将使用一种不同的动态规划公式来解决所有结点对最短路径问题,所产生的算法称为Floyd-Warshall算法,其运行时间为Θ(V3)。与前面的假设一样,负权重的边可以存在,但不能存在权重为负值的环路。如前所述,仍然按照动态规划过程的通常步骤来阐述算法。在对算法进行研究后,将提供一种类似的方法来找出有向图的传递闭包。
在Floyd-Warshall算法中,对一条最短路径的结构特征做出的描述与前面有所不同。Floyd-Warshall算法考虑的是一条最短路径上的中间结点,这里,简单路径p=<v1,v2,…,vl>上的中间结点指的是路径p上除v1和vl之外的任意结点,也就是处于集合{v2,v3,…,vl-1}中的结点。
Floyd-Warshall算法依赖于下面的观察。
假定图G的所有结点为V={1,2,…,n},考虑其中的一个子集{1,2,…,k},这里k是某个小于n的整数。对于任意结点对i,j∈V,考虑从结点i到结点j的所有中间结点均取自集合{1,2,…,k}的路径,并且设p为其中权重最小的路径(路径p是简单路径)。Floyd-Warshall算法利用了路径p和从i到j之间中间结点均取自集合{1,2,…,k-1}的最短路径之间的关系。该关系依赖于结点k是否是路径p上的一个中间结点:
1、如果结点k不是路径p上的中间结点,则路径p上的所有中间结点都属于集合{1,2,…,k-1}。因此从结点i到结点j的中间结点取自集合{1,2,…,k-1}的一条最短路径也是从结点i到结点j的中间结点取自集合{1,2,…,k}的一条最短路径。
2、如果结点k是路径上的中间结点,则将路径分解为:
如下图(记为图3)所示:
根据引理:
可知p1是从结点i到结点k的中间结点全部取自集合{1,2,…,k}的一条最短路径。又因为结点k不是路径p1上的中间结点,因此路径p1上的所有中间结点都属于集合{1,2,…,k-1}。因此p1是从结点i到结点k的中间结点全部取自集合{1,2,…,k-1}的一条最短路径。类似地,p2是从结点k到结点j的中间结点全部取自集合{1,2,…,k-1}的一条最短路径。基于上面的观察,可以定义一个不同于【1 最短路径和矩阵乘法】一节所描述的最短路径估计的递归公式。设dij(k)为从结点i到结点j的所有中间结点全部取自集合{1,2,…,k}的一条最短路径的权重。当k=0时,从结点i到结点j的一条不包括编号大于0的中间结点的路径将没有任何中间结点,这样的路径最多只有一条边,因此dij(0)=wij。根据上面的讨论,递归定义dij(k)如下:
因为对于任何路径来说,所有的中间结点都属于集合{1,2,…,n},因此矩阵D(n)=(dij(n))给出的就是最后答案:对于所有的i,j∈V,dij(n)=δ(i,j)。根据递归上图的公式,可以使用下下一张图中的自底向上的算法FLOYD-WARSHALL以递增次序来计算dij(k)的值。该算法的输入为一个n×n的矩阵W,该W由公式:
所定义。下面的算法(由于我不会在文章中打出hat或者上标符号,下面包括后面的许多内容是图片格式的,请见谅)返回的是最短路径权重矩阵D(n):
下图(记为图4)描述的是在图1上运行Floyd-Warshall算法所计算出的矩阵D(k):
Floyd-Warshall算法的运行时间由算法第3-7行的3层嵌套的for循环所决定。因为算法第7行的每次执行时间为O(1),因此该算法的运行时间为Θ(n3)。如【1 最短路径和矩阵乘法】一节中所描述的最终算法一样,该代码也非常紧凑,没有使用精巧的数据结构,隐藏在表述后面的常数比较小。因此即使对于输入规模为中等的图,Floyd-Warshall算法的效率也相当好。
在Floyd-Warshall算法中,可以有多种不同的方法来构建最短路径。一种办法是先计算最短路径权重矩阵D,然后从D矩阵来构造前驱矩阵Ⅱ。一旦给定了前驱矩阵Ⅱ,本博文最开始的PRINT-ALL-PAIRS-SHORTEST-PATH过程将打印出给定最短路径上的所有结点。
可以在计算矩阵D(k)的同时计算前驱矩阵Ⅱ。具体来说,将计算一个矩阵序列Ⅱ(0),Ⅱ(1),…,Ⅱ(n),这里Ⅱ=Ⅱ(n)并且定义πij(k)为从结点i到结点j的一条所有中间结点都取自集合{1,2,…,k}的最短路径上j的前驱结点。
可以给出πij(k)的一个递归公式:
1、当k=0时,从i到j的一条最短路径没有中间结点,因此:
2、对于k≥1,如果考虑路径i→k→j,这里k≠j,则所选择的结点j的前驱与选择的从结点k到结点j的一条中间结点全部取自集合{1,2,…,k-1}的最短路径上的前驱是一样的。否则所选择的结点j的前驱与选择的从结点i到结点j的一条中间结点全部取自集合{1,2,…,k-1}的最短路径上的前驱是一样的。即对于k≥1,有:
图4(往上面翻)描述的是进行此种合并后的算法在图1上运行时所生成的Ⅱ(k)矩阵序列。
给定有向图G=(V,E),结点集合为V={1,2,…,n},希望判断对于所有的结点对i和j,图G是否包含一条从结点i到结点j的路径。定义图G的传递闭包为图G*=(V,E*),其中E*={(i,j):如果图G中包含一条从结点i到结点j的路径}。
一种时间复杂度为Θ(n3)的计算图G的传递闭包的办法是给E中的每条边赋予权重1,然后运行Floyd-Warshall算法。如果存在一条从结点i到结点j的路径,则有dij<n;否则dij=∞。
还有另一种类似的办法,其时间复杂度也是Θ(n3),但在实际场景中能够节省时间和空间。该办法以逻辑或操作(∨)和逻辑与操作(∧)来替换Floyd-Warshall算法中的算术操作min和+。对于i,j,k=1,2,…,n,定义:如果图G中存在一条从结点i到结点j的所有中间结点都取自集合{1,2,…,k}的路径,则tij(k)为1;否则tij(k)为0。构建传递闭包G*=(V,E*)的方法为:将边(i,j)置于集合E*当且仅当tij(k)为1。递归定义tij(k)如下:
1、当k=0时,有:
2、对于k≥1,有:
下面以k递增的次序来计算矩阵T(k)=(tij(k)):
下图(记为图5):
描述的是在一个样本图上运行TRANSITIVE-CLOSURE(G)过程所计算出的矩阵序列T(k)。如Floyd-Warshall算法一样,TRANSITIVE-CLOSURE(G)过程的运行时间也是Θ(n3)。
3 用于稀疏图的Johnson算法
Johnson算法可以在O(V2·lgV+V·E)的时间内找到所有结点对之间的最短路径,算法要么返回一个包含所有结点对的最短路径权重的矩阵,要么报告输入图包含一个权重为负值的环路,算法需要调用此博文中的Dijkstra算法和Bellman-Ford算法作为自己的子程序。Johnson算法使用的技术称为重新赋予权重,该技术的工作原理如下。
如果图G=(V,E)中所有的边权重w皆为非负值,可以通过对每个结点运行一次Dijkstra算法来找到所有结点对之间的最短路径。如果图G包含权重为负值的边,但没有权重为负值的环路,那么只要计算出一组新的非负权重值,然后使用同样的方法即可。
注意:
下面的引理1描述的是可以很容易地对边的权重进行重新赋值来满足上面的第一个性质:引理1:
证明:首先来证明:
有:
因此:
最后:
根据:
有:
引理1得证。下一个目标是确保第二个性质,即对于所有的边(u,v)∈E,有:
为非负值。给定带权重的有向图G=(V,E),其权重函数为w:E→R。令一幅新图G’=(V’,E’),这里V’=V∪{s},s是一个新结点,s∉V,E’=E∪{(s,v):v∈V}。对权重函数w进行扩展,使得对于所有的结点v∈V,w(s,v)=0。注意,因为结点s没有进入的边,因此除了以s为源结点的最短路径外,图G’中没有其他包含结点s的最短路径。而且,图G’不包含权重为负值的环路当且仅当图G不包含权重为负值的环路。下图(记为图6)的(a)描述的是【对应于图1中的图G的】图G’:
下面是对上图的解释:
现在假定图G和图G’都不包含权重为负值的环路。定义对于所有的结点v∈V’,h(v)=δ(s,v)。根据三角不等式:
对于所有的边(u,v)∈E’,有h(v)≤h(u)+w(u,v)。因此如果根据公式:
来定义新的权重,则有:
至此满足了第二条性质。图6(b)描述的是对图6(a)的图进行权重重新赋值后的图G’。Johnson算法在执行过程中需要使用Bellman-Ford算法和Dijkstra算法作为子程序来计算所有结点对之间的最短路径。假定所有的边都保存在邻接链表里,算法要么返回一个|V|×|V|的矩阵D=(dij),其中dij=δ(i,j),要么报告输入图包含一个权重为负值的环路。对于所有结点对最短路径算法来说,通常假定结点的编号为1-|V|。
上述算法第1行生成图G’,第2行在图G’上运行Bellman-Ford算法,使用权重函数w和源结点s。如果图G’包含一条权重为负值的环路,算法的第3行将报告这个问题。算法第4-12行假定图G’不包含权重为负值的环路。第4-5行将h(v)的值设置为由Bellman-Ford算法所计算出来的最短路径权重δ(s,v)。算法的第6-7行计算新的权重:
对于每一对结点u,v∈V,算法的第9-12行的for循环通过调用Dijkstra算法来计算最短路径权重:
算法第12行将根据公式:
所计算出来的最短路径权重δ(u,v)保存在矩阵元素duv里。最后,算法的第13行返回构造完毕的矩阵D。图6描述的便是Johnson算法的执行过程。
END