这篇文章其实主要想说的是如何解决最短路径的问题。
其实最短路径问题,我们在生活中都在不知不觉的使用。比如我们在上网的时候,互联网传输采用了各种各样的数据包路由方法。这些路由算法都在幕后工作。
还有一些图算法的寻路操作,比如游戏中让游戏角色自动寻路。
其实寻找最短路径一般都是其他非图类算法中的一个重要的子程序。
下面开始说一下最短路径问题:
最短路径问题可以分为以下几种形式。例如,在有向图和无向图中找出最短路径,二者最重要的区别在于起点和目的地。
这个问题有很多种形式,比如从一个节点到其他所有节点的最短路径的源节点唯一的问题,从一个节点到另一个节点的最短路径的一对一问题,还有从其他所有节点到一个节点的最短路径的目标节点唯一问题,还有从所有节点到其他节点的最短路径的所有节点两两组合问题。
进去正文之前,我们先说一下松弛化和逐步改进的概念,具体的例子不再举了,大家知道是一种求解方法就好,这些方法会通过逐步接近的方式来获得相关问题的最佳解法。
然后大家看一下这个代码:
# 松弛技术
inf = float('inf')
def relax(W, u, v, D, P):
d = D.get(u, inf) + W[u][v]
if d < D.get(v, inf):
D[v], P[v] = d, u
return True
把图表示为字典的字典,用字典D存放距离值估计。P作为前导节点字典。前导指针构成了最短路径树。然后用公共代码的relax分解出松弛技术。D中不存在的项目都可以看做无穷大。也在主算法对它们初始化定义为无穷大。
这段代码,我们通过u来观察是否可以缩短路径,从而改进目前已知的到达v点的最短路径。不是最短路径的话,不再理会。是最短路径的话,就记住经过了哪些节点。
然后我们看第一个比较正式的算法,Bellman-Ford算法。这个算法适用于任意有向或者无向图的单源最短路径算法。如果图包含负环,算法将会输出信息,然后放弃查找。
# Bellman-Ford算法
def bellman_ford(G, s):
D, P = dict(s=0), dict()
for rnd in G:
changed = False
for u in G:
for v in G[u]:
if relax(G, u, v, D, P):
changed = True
if not changed:
break
else:
raise ValueError('negative cycle')
return D, P
这个算法包含了检查变化的charged,所以不需要多次迭代,程序可以提前终止。哈可以通过判断多余的迭代是否带来变化,来检测负环是否存在。
然后我们看一个权重图:
权重图改写为字典的字典,可以用如下代码表示:
a, b, c, d, e, f, g, h = range(8)
G = {
a: {b: 2, c: 1, d: 3, e: 9, f: 4},
b: {c: 4, e: 3},
c: {d: 8},
d: {e: 7},
e: {f: 5},
f: {c: 2, g: 2, h: 2},
g: {f: 1, h: 6},
h: {f: 9, g: 8}
}
算法具体的实现,可以用调试器,增加两条打印命令来显示松弛操作的边和分配给D的值。
如果测试几次,就会发现,g、h和f的距离值估计不断建议,到了第8轮还有减一的情况发生,这就提示了有负环的存在。
解决方法可以用del G[f][g]去掉(f, g)边。这样就不会出现负环的问题。
然后我们给单源DAG图最短路径问题加一个约束条件,环可以有,但边的权值不能为负值。
然后我们引入这个Dijkstra算法,这个算法的结构和Prim算法很类似,都使用优先级队列进行遍历。在Dijkstra算法中,优先级是距离值估计。然后我们看具体的算法实现:
# Dijkstra算法
from heapq import heappush, heappop
def dijkstra(G, s):
D, P, Q, S = {s: 0}, {}, [(0, s)], set()
while Q:
_, u = heappop(Q)
if u in S:
continue
S.add(u)
for v in G[u]:
relax(G, u, v, D, P)
heappush(Q, (D[v], v))
return D, P
这个算法的运行时间是Q((m+n)lgn),其中m为边数,n为节点数。
和Dijkstra算法关系比较密切的还有广度优先搜索(BFS)算法。我们可以考虑一下边权值为正整数的情况,将权值为w的一条边替换为w-1条无权边,连成一条由虚拟节点组成的路径,比如下图:
BFS总会找到正确答案,但是效率会变得很低。这时候我们会发现BFS解决问题的方式和Dijkstra算法十分类似:在各条边上华为的时间与边权值成正比,所以会按照距离大小从起始节点依序到达各节点。
然后下面说多对多问题:
这里引入Johnson算法,算法动机其实很简单,解决稀疏图矩阵所有节点对之间的最短路径问题,对各个节点相同的情况下使用Dijkstra算法。但是Dijkstra算法不允许负权边的存在。对于单源最短路径的问题,除了改用Bellman-Ford算法,并没有其他办法。
然后我们可以这样设想,增加一个新的节点s,它到所有现在的节点的边权值为0.然后对从s出发的情况运行Bellman-Ford算法。这样可以计算出s到图中每个节点的距离。我们称为h(v)。然后我们可以用h调节各边的权值。可以如下定义:w'(u, v) = w(u, v) + h(u) - h(v)。
这样定义可以保证新的权值非负值。然后也不会引入新的干扰项。
然后再用Dijkstra算法发现最短路径,再逆向变换所有的路径长度。具体的Johnson算法实现如下:
# 求解稀疏矩阵图
from copy import deepcopy
def johnson(G):
G = deepcopy(G)
s = object()
G[s] = {v: 0 for v in G}
h, _ = bellman_ford(G, s)
del G[s]
for u in G:
for v in G[u]:
G[u][v] += h[u] - h[v]
D, P = dict(), dict()
for u in G:
D[u], P[u] = dijkstra(G, u)
for v in G:
D[u][v] += h[v] - h[u]
return D, P
Johnson算法效率也还不错,不过它的的运行时间是Dijkstra算法运行时间的n倍。
然后开始说下一个问题,求解所有节点对之间最短距离的方法:
这个算法名字叫做Floyd-Warshall算法。基于动态编程的原理,其实Dijkstra算法也是基于动态编程。
Floyd-Warshall基于缓存式递归分解,实现的过程一般都具有迭代性。
我们需要寻找一组递归相关的子问题。
我们随意对节点排序,并限制允许用于构成最短路径的中间节点的数量,即前k个。
这里直接用三个参数对子问题进行参数化:起始节点
终止节点
允许经过的最大节点编号
然后设节点u到节点v的最短路径的长度为d(u, v, k):
d(u, v, k) = min(d(u, v, k-1), d(u, k, k-1) + d(k, v, k-1))。
这个和背包问题一样,要考虑是否包括节点k。不包括,就是d(u, v, k-1)。包括就必须使用到达k的最短路径d(u, k, k-1)和从k出来的最短路径d(k, v, k-1)。
然后下面看代码:
# Floyd-Warshall算法的缓存式递归实现
def rec_floyd_warshall(G):
@memo
def d(u, v, k):
if k == 0:
return G[u][k]
return min(d(u, v, k-1), d(u, k, k-1) + d(k, v, k-1))
return {(u, v): d(u, v, len(G)) for u in G for v in G}
# 记忆体化的装饰器的函数
from functools import wraps
def memo(func):
cache = dict()
@wraps(func)
def wrap(*args):
if args not in cache:
cache[args] = func(*args)
return cache[args]
return wrap
然后我们可以尝试一下迭代的版本:
三个参数,就需要三个for循环。下面直接上代码:
# Floyd-Warshall算法,仅考虑距离
def floyd_warshall(G):
D = deepcopy(G)
for k in G:
for u in G:
for v in G:
D[u][v] = min(D[u][v], D[u][k] + D[k][v])
return D
D为当前的距离图,先前的距离图为C。全程只用一个距离图换算之后,公式由:
D[u][v] = min(D[u][v], C[u][k] + C[k][v])。
变换为:
D[u][v] = min(D[u][v], D[u][k] + D[k][v])。
如果再加入一个P矩阵的话,P[u][v]将被替换为P[k][v]。代码就会变为:
# Floyd-Warshall算法
def floyd_warshall(G):
D, P = deepcopy(G), dict()
for u in G:
for v in G:
if u == v or G[u][v] == inf:
P[u, v] = None
else:
P[u, v] = u
for k in G:
for u in G:
for v in G:
shortcut = D[u][k] + D[k][v]
if shortcut < D[u][v]:
D[u][v] = shortcut
P[u, v] = P[k, v]
return D, P
这个地方应为shortcut < D[u][v],而不是shortcut <= D[u][v],因为在某些情况,最后一步为D[v][v],这时候将导致P[u, v] = None。
然后说下一个问题,说一下“中途相遇的问题”:
Dijkstra算法子问题的解,如果节点为s和t的话,将由s到t在图上不断扩散。但是如果从起点和终点同时出发,展开遍历,这样就会减少很多工作量。具体抽象出来就像下面的图:
把原来的稍作修改就可以变为,Dijkstra算法的双向图版本。可以变成一个子解生成器,让我们尽可能提取更多的子解。
这样必须抛开距离表,仅仅依靠优先队列中保存的距离值。于是放上示例代码:
# Dijkstra算法作为解决方案生成器的实现
from heapq import heappush, heappop
def idijkstra(G, s):
Q, S = [(0, s)], set()
while Q:
d, u = heappop(Q)
if u in S:
continue
S.add(u)
yield u, d
for v in G[u]:
heappush(Q, (d+G[u][v], v))
没有引入前导节点信息,但是可以通过向堆元组中添加前导节点来扩展这个解决方案。获取距离表,只需调用dict(idijkstra(G, s))。
但是从两个节点同时出发,可能会遇到下面的问题:
首次相遇的节点(高亮的节点)并不一定位于最短路径(高亮的边)上。
其实这个问题只需要对终止条件做限定就好了,只要看他们已经走了多远,就是目前获得的最新的距离,这个是不能减小的。二者综合至少和我们目前发现的最佳路径相等,这样就不可能再找到更好的方案。
如果G是无向图的话,并且任意节点u满足G[u][u] = 0,则可以使用下面的代码:
# Dijkstra算法的双向图版本
from itertools import cycle
def bidir_dijkstra(G, s, t):
Ds, Dt = dict(), dict()
forw, back = idijkstra(G, s), idijkstra(G, t)
dirs = (Ds, Dt, forw), (Dt, Ds, back)
try:
for D, other, step in cycle(dirs):
v, d = next(step)
D[v] = d
if v in other:
break
except StopAsyncIteration:
return inf
m = inf
for u in Ds:
for v in G[u]:
if not v in Dt:
continue
m = min(m, Ds[u] + G[u][v] + Dt[v])
return m
当然,你也可以很轻松的拓展这段代码。
下面我们会引入A*算法:
如果真的知道哪一条边可以离目标更近,我们就可以用贪婪算法解决问题。直接沿着最短的路径移动,不用理会其他的旁支路线。
A*算法有点像人工智能中启发式搜索的概念。而不是像DFS和BFS那样子盲目搜索。也不像Dijkstra算法一样对未来走向一无所知。
因为A*加入了一个潜在势函数,也可以叫做启发式函数h(v)。
就像上面介绍的Johnson算法一样,我们可以定义修正后的边权:
w'(u, v) = w(u, v) - h(u) + h(v)
然后你会发现,这样子调整之后,我们可以奖励正确的边,惩罚不正确的边。然后给各个边的权值加上了启发式番薯。这个算法将沿导致剩余距离下降的方向发展。
A*算法相当于针对修正图的Dijkstra算法,h可行的话,算法就是正确的。
h(s)是一个常数,所以我们只用将h(v)加入现有的优先级的队列中。这个和是我们对从s到t的最佳估计。如果w'(u, v)可行,那么h(v)也会是d(v, t)的下界。
然后这只直接给出代码,里面调用了上面代码的idijkstra函数:
# A*算法
from heapq import heappush, heappop
inf = float('inf')
def a_star(G, s, t, h):
P, Q = dict(), [(h(s), None, s)]
while Q:
d, p, u = heappop(Q)
if u in P:
continue
P[u] = p
if u == t:
return d - h(t), P
for v in G[u]:
w = G[u][v] - h(u) + h(v)
heappush(Q, (d + w, u, v))
return inf, None
A*算法的优势就在于启发式函数,至于具体的启发式函数是什么取决于要解决的问题。
另外,A*算法也可以搜索解空间。可以解决魔方问题或者单词梯问题。这里说一下后一个问题。
单词梯是从第一个起始单词开始构建的,比如 lead,然后以另外一个单词作为结尾,比如gold。单词梯每一步搭建都要用到实际单词,从一个单词推进到下一个单词,只能更换一个字母。比如这个题目的一个解法就可以通过 load 和 goad 这两个单词,到达 lead 和 gold。如果每个单词看做一个节点,我们可以加上边。可能没有必要这样建立,不过这样假设,便可以做出如下代码:
# 单词梯路径的隐式图
from string import ascii_letters as chars
def variants(wd, words):
wasl = list(wd)
for i, c in enumerate(wasl):
for oc in chars:
if c == oc:
continue
wasl[i] = oc
ow = ''.join(wasl)
if ow in words:
yield ow
words[i] = c
class WordSpace:
def __init__(self, words):
self.words = words
self.M = dict()
def __getitem__(self, wd):
if wd not in self.M:
self.M[wd] = dict.fromkeys(self.variants(wd, self.words), 1)
return self.M[wd]
def heuristic(self, u, v):
return sum(a!=b for a, b in zip(u, v))
def ladder(self, s, t, h=None):
if h is None:
def h(v):
return self.heuristic(v, t)
_, p = a_star(self, s, t, h)
if p is None:
return [s, None, t]
u, p = t, []
while u is not None:
p.append(u)
u = P[u]
p.reverse()
return p
WordSpace可以作为加权图使用,配合a_star使用。这里的启发式函数只是统计单词不同的字符位数。G是WordSpace的对象,G['lead'] 就是字典,其他单词是键值。各个边的权值是1。
文章写到了这里,差不多就说完了。
A*算法一般在人工智能的书里面才会涉及到。其他的一些算法在一般的算法书中也能找得到。
设计新的算法的时候,比如Dijkstra算法的双向版本和A*的启发式结合在一起,很可能因为一些陷阱而导致算法无效。
以上。谢谢大家关注。
今天立冬,祝大家冬天愉快,所有的不开心伴随着秋天的过去都过去了。新的季节来了,希望大家都有好心情。
天气寒冷,大家注意保暖。