最短路径算法

发表于 | 分类于 | 阅读次数:

综述

本文主要介绍最短路径的相关概念、单源最短路径和所有结点对的最短路径。本文的完整代码可以在我的github找到。

概念

对于一个带有权重的有向图$G=(V, E)$和权重函数$w: E \rightarrow R$，该权重函数将每条边映射到实数值的权重上。图中一条路径$p=(v_0, v_1, v_2, …, v_k)$的权重$w(p)$是构成该路径的所有边的权重之和：

$$w(p)=\sum_{i=1}^{k}w(v_{i-1}, v_i)$$
定义从结点u到结点v的最端路径权重$\delta(u, v)$如下：

$$\delta(u, v) = \begin{cases}
\min\lbrace w(p): u \leadsto v\rbrace & 如果存在一条从u到v的路径 \\
\infty & 其他
\end{cases}$$

从结点u到结点v的最短路径定义为任何一条权重为$w(p)=\delta(u, v)$的从u到v的路径p。

最短路径问题由几种 变体，接下来我们主要讨论单源最短路径问题和所有结点对的最短路径问题。

单源最短路径

单源最短路径问题：给定$G=(V, E)$，找到从源结点$s \in V$到每个结点$v \in V$的最短路径。有时单源最短路径问题可能会包含具有负权重的边。但如果图$G=(V, E)$不包含从源结点s可以到达的权重为负值的环路，则对于所有的结点$v \in V$都有最短路径（可能为无穷）。如果环路的权重为正值，则此环路必定不存在于最短路径上；如果环路的权重为负值，则最短路径没有定义，因为我们总是可以尽可能多的包含此环路，使其最短路径为$-\infty$；如果环路的权重为0， 则必然有不包含环路的最短路径。

最短路径通常同前驱子图来表示，我们记录每个结点v的前驱结点$v.\pi$。我们定义前驱子图$G_\pi=(V_\pi, E_\pi)$，其中$V_\pi=\lbrace v \in V: v.\pi \neq NIL\rbrace \bigcup \lbrace s\rbrace $，有向边集合$E_\pi=\lbrace (v.\pi, v) \in E: v \in V.\pi -\lbrace s\rbrace \rbrace $。

我们定义图如下：

1
2
3
4
5
6
INF = float('inf')
class Graph:
    def __init__(self, n):
        self.n = n
        self.e = []
        self.w = [[0 if i == j else INF for j in range(n)] for i in range(n)]

Bellman-Ford算法

Bellman-Ford算法逐渐通过松弛操作降低从源结点s到结点v的最短路径的估计值$d[v]$，直到达到最短路径$\delta(s, v)$为止。松弛操作的python的代码如下：

1
2
3
4
def relax(g, p, d, u, v):
    if d[v] > d[u] + w[u][v]:
        d[v] = d[u] + w[u][v]
        p[v] = u

Bellman-Ford算法的python代码示例如下：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
def bellman_ford(g, s):
    d, p = [], []
    for i in range(g.n):
        d.append(INF)
        p.append(None)
    d[s] = 0
    for i in range(g.n - 1):
        for u, v in g.e:
            # 松弛操作
            if d[v] > d[u] + g.w[u][v]:
                d[v] = d[u] + g.w[u][v]
                p[v] = u
    for u, v in g.e:
        if d[v] > d[u] + g.w[u][v]:
            return False, p, d
    return True, p, d

松弛操作的原理也很简单，因为边(u, v)存在，所以最短路径的估计值d[v]必定小于或者等于d[u] + w[u][v]。d[v]相当于是最短路径$\delta(s, v)$的一个上界。Bellman-Ford算法每次对所有的边进行松弛操作， 共循环|V|-1次，使得d[v]最终成为最短路径。至于其准确性简单说明如下。首先我们引入定理：

最短路径的子路径也为最短路径。

即如果从$v_0$到$v_k$有一条最短路径$p=( v_0, v_1, …, v_{k-1}, v_{k})$，则$p_{ij}=(v_i, v_{i+1}, …, v_{j-1}, v_j)$也是从$v_i$到$v_j$的最短路径。那么在第一次循环时，对于$p_{sv}=(s, v)$的结点v，即最短路径只有一条边的结点，d[v]即为其最小距离。对于最短路径只有两条边的结点在第二次循环时即可求出，依次类推，因为s到v的简单路径上最多包含|V|-1条边，则循环|V|-1次即可求解。若最短路径不包含权值为负的环路，则最短路径算法返回True。容易知道，此算法的时间复杂度为O(VE)

Dijkstra算法

Dijkstra算法同样可以解决带权重的有向图的单源最短路径问题，不过它要求所有边的权重都为非负值。Dijkstra算法维护一个结点集S，从源结点s到达该集合的每一个结点的最短路径均已找到，每次从V - S中取出最短路径估计最小的结点u将其加入S。其python代码示例如下：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
def find(d, flag):
    index, min_dist = -1, INF
    for i in range(len(d)):
        if flag[i] is False and d[i] < min_dist:
            index = i
            min_dist = d[i]
    if index == -1:
        return None
    return index

def dijkstra(g, s):
    d, p, flag = [], [], []
    for _ in range(g.n):
        d.append(INF)
        p.append(None)
        flag.append(False)

    d[s] = 0
    while True:
        u = find(d, flag)
        if u is None:
            break
        flag[u] = True
        for v, weight in enumerate(g.w[u]):
            if weight < INF:
                if d[v] > d[u] + weight:
                    d[v] = d[u] + weight
                    p[v] = u

    return p, d

从代码中，易知其时间复杂度为$O(V^2)$。但如果像Prim算法那样使用设计良好的最小优先队列，并使用邻接链表法，则可以使其使见复杂度降低至O((V+E)logV)，即O(ElogV)。

有向无环图的最短路径算法

对于有向无环图我们可以利用拓扑排序来优化最短路径的时间，其伪代码示例如下：

1
2
3
4
5
6
7
8
DAG-SHORTEST_PATH(G, w, s)
    利用深度优先搜索对G进行拓扑排序
    初始化父亲数组p，距离数组d
    for u in 拓扑排序后的结点
        for v in G.adj[u]:
            if d[v] > d[u] + w[u][v]
                d[v] = d[u] + w[u][v]
                p[v] = u

按照拓扑排序进行遍历的原因如下：如果存在一条边(u, v)，则结点u的拓扑排序在结点v的前面。当前循环的结点到达v时，d[v]已经是最短路径了，因为在此之前所有能到达v的结点u已经被遍历了。依此类推，当遍历结束，最短路径就都求出来了。当使用邻接链表的表示方法时，其时间复杂度度为O(V+E)。但当使用邻接矩阵的表示方法时，其时间复杂度为$O(V^2)$。使用邻接矩阵法的有向无环图的最短路径算法的python代码如下：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
def dfs(g, u, f, flag=None):
    flag = flag if flag else g.n * [False]
    flag[u] = True
    for v, weight in enumerate(g.w[u]):
        if flag[v] is False and weight < INF:
            dfs(g, v, f, flag)

    f.insert(0, u)
def dag_short_path(g, s):
    d, p = [], []
    for i in range(g.n):
        d.append(INF)
        p.append(None)
    d[s] = 0
    f = []
    dfs(g, s, f)
    for u in f:
        for v, weight in enumerate(g.w[u]):
            if weight < INF:
                if d[v] > d[u] + weight:
                    d[v] = d[u] + weight
                    p[v] = u
    return p, d

从代码中可以看到我们并没有对整个图进行深度优先遍历，我们只对从源结点s出发的结点进行了深度优先遍历。因此我们所得到的拓扑排序f为以s为根的深度优先树的拓扑排序。因为如果源结点s可以到达结点v，则v必定在拓扑排序中。上述实现因为使用邻接矩阵，因此时间复杂度为$O(V^2)$。

所有结点对的最短路径

平方重复法

所有结点对的最短路径算法可以使用动态规划的方法，其基本思路如下：假设从结点i到j存在一条最短路径p，且p最多包含m条边，假定没有权重为负值的环路。设路径p上结点j的父结点为k，即$p=(i, …, k, j)$，则$p_1=(i, …, k)$为从结点i到k的最短路径。则递推式为

$$\delta(i, j)=min(\delta(i, k) + w[k][j])　　\forall(k, j) \in E$$
则按照这个思路实现的伪代码如下：

1
2
3
4
5
6
7
8
9
10
EXTEND-SHORTEST-PATHS(L, W)
    n = L.rows
    let L' be a new n*n matrix
    for i =1 to n
        for j = 1 to n
            L'[i][j]=INF
            for k = 1 to n
                L'[i][j] = min(L'[i][j], L[i][k]+W[k][j])

    return L'

L初始化为权重矩阵W。因此需要进入这样的函数|V| - 1次，因此其时间复杂度为$O(V^4)$。这样的过程类似于矩阵的乘法，对于矩阵的乘法我们有以下的伪代码：

1
2
3
4
5
6
7
8
9
10
SQURE-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] = C[i][j] + A[i][k] * A[k][j]
            
    return C

比较之下，发现两者非常相似。将min替换成+，将+替换成*，即完全一致。对于矩阵乘法我们有:

$$L^{(1)} = W \\
L^{(2)} = W * W \\
L^{(4)} = W^2 * W^2 \\
L^{(8)} = W^4 * W^4 \\
…\\
L^{(2m)} = W^m * W^m$$
因此我们在求解最短路径时，可以利用上述性质，因为只要运行EXTEND-SHORTEST-PATHS函数|V|-1次必然得到最短路径，即$W^{n-1}, n=|V|$必然为最短路径，且之后再次运行此函数，最短路径矩阵将不会改变。因此我们只需要利用上述的原理即可得到如下的伪代码：

1
2
3
4
5
6
7
8
9
10
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 = 2m
    
    return L[m]

因为我们只需要得到$W^{n-1}$，而$W^{m} = W^{n-1} \forall m>n-1$。此时时间复杂度降低至$O(V^3logV)$。其python代码示例如下：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
def extend_paths(l, w, p):
    n = len(l)
    y = [[INF for _ in range(n)] for _ in range(n)]
    pp = [[None for _ in range(n)] for _ in range(n)]
    for i in range(n):
        for j in range(n):
            pp[i][j] = p[i][j]
            for k in range(n):
                if i == j:
                    y[i][j] = 0
                    pp[i][j] = i
                else:
                    if y[i][j] > l[i][k] + w[k][j] and k != j:
                        y[i][j] = l[i][k] + w[k][j]
                        pp[i][j] = p[k][j]

    return y, pp
    
def all_short_paths(g):
    n = g.n
    l = g.w
    p = [[i if g.w[i][j] < INF else None for j in range(n)] for i in range(n)]
    m = 1
    while m < n - 1:
        l, p = extend_paths(l, l, p)
        m = 2 * m

    return l, p

需要指出的是以上代码虽然能正确的得到最短路径d，但是记录前驱结点的数组p不一定准确。$p[i][j]=k$表示从源结点i到j的路径中j的前驱结点为k。数组p只有在图中存在权重为0的环时才可能不准确，大部分情况下都是准确的。

Floyd-Warshall算法

Floyd-Warshall算法使用了一种不同的动态规划方法，这种动态规划方法更加抽象，如果让我子集设计一个所有结点对的最短路径算法，我可能倾向于会想到上面的重复平方法，而非这种方法。但事实证明，设计更有抽象的最优子问题，往往会得到更高效的动态规划算法。

Floyd-Warshall算法考虑的是最短路径上的中间结点。一条简单路径$p=\lbrace v_1, v_2, …, v_l\rbrace $上的中间结点即指p上除了$v_1$和$v_l$的任意结点。Floyd-Warshall算法定义最有子问题如下：

对于任意结点i和j，从结点i到结点j的最短路径p的中间结点全部取自$\lbrace 1, 2, …, k\rbrace $, 其中结点的全集为$V=\lbrace 1, 2, 3, …, n\rbrace $

如果从i到j的最短路径的中间结点均取自$\lbrace 1, 2, …, k\rbrace $，则有两种情况：

  1. 如果k不在最短路径上，则从i到j的中间结点取自结点集$\lbrace 1, 2, …, k-1\rbrace $的最短路径也为从i到j的中间结点取自结点集$\lbrace 1, 2, …, k\rbrace $的最短路径。
  2. 如果k在最短路径上，则将路径分解为从i到k再到j，即$i \rightsquigarrow k \rightsquigarrow j$，则从i到k的一条最短路径取自$\lbrace 1, 2, …, k-1\rbrace $，同样从k到j的一条最短路径也取自$\lbrace 1, 2, …, k-1\rbrace $。

子问题分解完成。当k=0时，表示从结点i到j一条最短路径不包含任何的中间结点。记$d_{ij}^{(k)}$表示从结点i到j的所有中间结点取自$\lbrace 1, 2, …, k\rbrace $的最短路径的权重。则$d_{ij}^{(0)}=w_{ij}$。其$d_{ij}^{(k)}$的递归式如下：

$$
d_{ij}^{(k)} = \begin{cases}
w_{i, j} & k = 0 \\
min(d_{ij}^{(k-1)}, d_{ik}^{(k-1)}+d_{kj}^{(k-1)}) & k\geq 1 \\
\end{cases}
$$

当k=n时，矩阵$D^{(n)}=(d_{ij}^{(n)})$接触所有结点对的最短路径。我们可以再计算矩阵$D^{(k)}=(d_{ij}^{(k)})$的同时计算前驱矩阵$P^{(k)}$。与前面类似，$p_{ij}^{(k)}$表示从结点i到j的一条所有中间结点都取自集合$\lbrace 1, 2, …, k\rbrace $的最短路径上j的前驱结点。则当k=0时，从i到j的最短路径没有中间结点，则有

$$
p_{ij}^{(0)} = \begin{cases}
NIL & i=j或w_{ij}=\infty\\
i & i\neq j或w_{ij}\lt\infty \\
\end{cases}
$$

其递归式如下：

$$
p_{ij}^{(k)} = \begin{cases}
p_{ij}^{(k-1)} & d_{ij}^{(k-1)} \leq d_{ik}^{(k-1)} + d_{kj}^{(k-1)} \\
p_{kj}^{(k-1)} & d_{ij}^{(k-1)} \gt d_{ik}^{(k-1)} + d_{kj}^{(k-1)} \\
\end{cases}
$$

根据以上的递归式，其python实现如下：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
def floyd_warshall(g):
    n = g.n
    d, p = [], []
    for i in range(n):
        d.append([])
        p.append([])
        for j in range(n):
            d[i].append(g.w[i][j])
            if i != j and g.w[i][j] < INF:
                p[i].append(i)
            else:
                p[i].append(None)

    for k in range(n):
        for i in range(n):
            for j in range(n):
                if d[i][j] > d[i][k] + d[k][j]:
                    d[i][j] = d[i][k] + d[k][j]
                    p[i][j] = p[k][j]

    return d, p

相比于重复平方法，floyd_warshall算法更加简洁和高效。其时间复杂度为$O(V^3)$，空间复杂度为$O(V^2)$。而且其能够克服上述所说的权值为0的环路的缺点。

Reference

本文主要参考《算法导论》。