第 25 章 所有结点对的最短路径

本章对应原书第 25 章。第 24 章从一个源点出发；本章一次求出每一对顶点之间的最短路径，输出一张 \(n\times n\) 的距离表——就像地图册里的城市距离表，或者交易台上"任意两种货币之间最优兑换率"的矩阵。本章有三条路线：把最短路径递推看成一种特殊的矩阵乘法；Floyd-Warshall 动态规划；以及用"势函数"重新赋权、把负权变成非负权的 Johnson 算法。最后一条和第 24 章的无套利条件、第 29 章的对偶变量是同一个思想。

学习目标

读完本章，你应当能够：

1. 把"再多走一条边"的递推写成 \((\min,+)\) 半环上的矩阵乘法，并用重复平方在 \(\Theta(n^3\lg n)\) 内求全源最短路径。
2. 写出 Floyd-Warshall 的递推式 \(d_{ij}^{(k)}=\min(d_{ij}^{(k-1)},d_{ik}^{(k-1)}+d_{kj}^{(k-1)})\)，理解"按允许的中间顶点集合递推"，并能同时维护前驱矩阵、求传递闭包。
3. 用引理 25.1 解释 Johnson 算法的重新赋权为什么不改变最短路径，以及为什么必须用势函数而不能"统一加一个常数"。
4. 用对角线 \(d_{ii}<0\) 识别负权环（在汇率图中即识别处于套利环上的资产）。
5. 用 numpy 把 Floyd-Warshall 的内两层循环向量化，并在多资产最优兑换表、因子流水线依赖分析中应用。

读前导读

这一章在解决什么问题。 第 24 章回答"从一个起点出发，到每个地方最便宜要花多少"；本章一次回答"任意两个地方之间最便宜要花多少"，输出一张 \(n\times n\) 的表。你最熟悉的类比是交易台上的交叉汇率表：行是卖出的货币，列是买入的货币，格子里是"最优可得兑换率"。直接报价不一定最优，绕道美元或 BTC 可能更划算，本章算法会把所有绕道都考虑进去。这张表的对角线也有意义：从 USD 出发绕一圈回到 USD，正常应该"不赚不亏"（成本 0）；若对角线出现负数，说明这种货币处在一个套利环上。

本章有三条路线，各有一个值得带走的思想。第一，"再多走一步"的递推形式上和矩阵乘法一模一样，只是把"乘"换成"加"、把"求和"换成"取最小"；所以矩阵乘法的技巧（结合律、快速幂、向量化）都能借用。第二，Floyd-Warshall 用"逐步开放中转站"的方式递推，三重循环就写完，是工程中最常用的全源算法。第三，Johnson 算法给每个顶点配一个"价格" \(h\)，用 \(\hat w(u,v)=w(u,v)+h(u)-h(v)\) 把负权都变成非负，再借用快速的 Dijkstra。这个"价格"与第 24 章套利检测中的势函数、CFA 里"无套利 ⟺ 存在一组自洽的定价"是同一个东西：约化后的权重就是"按公允价计算，这笔交易多付了多少"。

需要先想起来的数学。

• 矩阵乘法与下标。 \(C=AB\) 的元素是 \(c_{ij}=\sum_k a_{ik}b_{kj}\)：\(A\) 的第 \(i\) 行与 \(B\) 的第 \(j\) 列逐项相乘再相加。例：\(\begin{pmatrix}1&2\\3&4\end{pmatrix}\begin{pmatrix}5\\6\end{pmatrix}=\begin{pmatrix}1\cdot5+2\cdot6\\3\cdot5+4\cdot6\end{pmatrix}=\begin{pmatrix}17\\39\end{pmatrix}\)。矩阵乘法满足结合律 \((AB)C=A(BC)\)，但一般不满足交换律。见 第 00 册第 06 章 线性代数速成。
• 上标不是乘方。 本章 \(l_{ij}^{(m)}\)、\(d_{ij}^{(k)}\) 里带括号的上标是"第几轮"的编号，不是 \(m\) 次方。只有 \(W^2,W^4\) 这种不带括号的写法才是"矩阵的幂"，而且是 \((\min,+)\) 意义下的幂。
• 望远镜求和。 若各项是"前后之差"，相加后中间全部抵消：\((h_0-h_1)+(h_1-h_2)+(h_2-h_3)=h_0-h_3\)。这是引理 25.1 的全部秘密。见 第 00 册第 07 章 概率中的分析工具。
• 对数与 \(\lg\)。 \(\lg n\) 是以 2 为底的对数，"把 \(n\) 对半分几次到 1"；重复平方只需 \(\lceil\lg(n-1)\rceil\) 次（\(\lceil x\rceil\) 表示向上取整，\(\lceil 2.3\rceil=3\)）。汇率图中 \(w=-\ln R\)、最优兑换率 \(=e^{-d}\) 的来回换算，见 第 00 册第 04 章 级数与收敛。
• 归纳法与"动态规划"。 动态规划就是"把大问题拆成同类小问题，按从小到大的顺序把答案填进表里"；正确性靠归纳法：小问题的答案对了，大问题也就对了。见 第 00 册第 08 章 读懂数学证明与符号 和本册第 15a 章。

怎么读这一章。 必读：25.1（表示法）、25.3.1 与 25.3.3（Floyd-Warshall 递推和对角线检测负环）、25.4.2（引理 25.1）、25.4.5（金融解读）、25.6.2（多资产兑换表）。25.2 的 \((\min,+)\) 矩阵乘法建议读懂 25.2.2 的替换表即可，重复平方第一次可以只看结论。25.3.2 构造路径、25.3.4 传递闭包、25.5 思考题与注记第一次可以略读。建议先把第 24 章的 24.9 节读完再来，本章 25.4 和 25.6.2 是它的直接延续。

25.1 问题与表示

给定带权有向图 \(G=(V,E)\)，\(w:E\to\mathbb R\)，求所有 \(u,v\) 间的最短路径权重 \(\delta(u,v)\)。

朴素做法：从每个顶点各跑一次单源算法。

• 边权非负时跑 \(|V|\) 次 Dijkstra：数组实现 \(O(V^3)\)；二叉堆 \(O(VE\lg V)\)；斐波那契堆 \(O(V^2\lg V+VE)\)。
• 有负权边时只能跑 \(|V|\) 次 Bellman-Ford：\(O(V^2E)\)，稠密图上是 \(O(V^4)\)。

本章要做得比这更好。

输入：与第 24 章用邻接表不同，本章多数算法用邻接矩阵。顶点编号 \(1..n\)，输入 \(n\times n\) 矩阵 \(W=(w_{ij})\)：

允许负权边，暂时假定没有负权环。

白话解释：邻接矩阵就是一张"报价表"：第 \(i\) 行第 \(j\) 列写着"从 \(i\) 直接到 \(j\) 的代价"。对角线写 0（原地不动不花钱），没有直接通道的格子写 \(\infty\)（没有报价）。在汇率图里，\(w_{ij}=-\ln R[i,j]\)，\(\infty\) 对应 \(R=0\)，即没有这个盘口。本章算法的输出 \(D\) 是同样形状的表，但格子里换成了"允许任意绕道后的最低代价"。

输出：距离矩阵 \(D=(d_{ij})\)，\(d_{ij}=\delta(i,j)\)；以及前驱矩阵（predecessor matrix）\(\Pi=(\pi_{ij})\)：\(i=j\) 或无路径时为 NIL，否则是 \(i\) 到 \(j\) 某条最短路径上 \(j\) 的前驱。\(\Pi\) 的第 \(i\) 行诱导出一棵以 \(i\) 为根的最短路径树。打印路径的递归很简单：

def print_all_pairs_shortest_path(Pi, i, j):
    if i == j:
        print(i)
    elif Pi[i][j] is None:
        print(f"no path from {i} to {j} exists")
    else:
        print_all_pairs_shortest_path(Pi, i, Pi[i][j])
        print(j)

25.2 最短路径与矩阵乘法

25.2.1 按"边数"做动态规划

按第 15a 章动态规划的步骤来。

最优解的结构：由引理 24.1，最短路径的子路径也是最短路径。设 \(i\) 到 \(j\) 的最短路径 \(p\) 至多含 \(m\) 条边。若 \(i=j\)，权重为 0；否则把 \(p\) 拆成 \(i\overset{p'}{\leadsto}k\to j\)，\(p'\) 至多 \(m-1\) 条边且是 \(i\) 到 \(k\) 的最短路径，于是 \(\delta(i,j)=\delta(i,k)+w_{kj}\)。

递归式：令 \(l_{ij}^{(m)}\) 为 \(i\) 到 \(j\) 至多含 \(m\) 条边的路径的最小权重。

（原式还有一项 \(l_{ij}^{(m-1)}\)，但因为 \(w_{jj}=0\)，它已包含在 \(k=j\) 的项中。）无负权环时最短路径是简单路径，至多 \(n-1\) 条边，因此

推导拆解：(25.2) 的意思是"最多走 \(m\) 步到 \(j\) 的最好方案 = 在所有可能的倒数第二站 \(k\) 中，挑'最多 \(m-1\) 步到 \(k\) 的最好方案 + 最后一步 \(k\to j\)'最便宜的那个"。 小例子：三种货币 1、2、3，成本 \(w_{12}=1\)，\(w_{23}=1\)，\(w_{13}=5\)，其余无边。\(l_{13}^{(1)}=w_{13}=5\)（只能直接换）。算 \(l_{13}^{(2)}\) 时，倒数第二站 \(k\) 有三种选择：\(k=1\) 得 \(0+5=5\)；\(k=2\) 得 \(l_{12}^{(1)}+w_{23}=1+1=2\)；\(k=3\) 得 \(l_{13}^{(1)}+w_{33}=5+0=5\)。取最小，\(l_{13}^{(2)}=2\)。\(k=3\) 那一项就是括号里说的"原地不动"，因为 \(w_{33}=0\)，它自动包含了"不多走一步"的选项。 (25.3) 说递推到 \(n-1\) 步以后就不再变化：没有负权环时最短路径不会重复经过顶点，\(n\) 个顶点的简单路径最多 \(n-1\) 条边，允许更多步也用不上。

自底向上：\(L^{(1)}=W\)，依次算 \(L^{(2)},\dots,L^{(n-1)}\)。核心过程把已有的最短路径再延伸一条边：

def extend_shortest_paths(L, W):            # Θ(n^3)
    n = len(L)
    Lp = [[INF] * n for _ in range(n)]
    for i in range(n):
        for j in range(n):
            for k in range(n):
                Lp[i][j] = min(Lp[i][j], L[i][k] + W[k][j])
    return Lp

25.2.2 这就是矩阵乘法

对比普通矩阵乘法 \(c_{ij}=\sum_k a_{ik}\cdot b_{kj}\)。在 (25.2) 中做替换

EXTEND-SHORTEST-PATHS 就变成了 SQUARE-MATRIX-MULTIPLY。换句话说，它是 \((\min,+)\) 半环（也叫热带代数，tropical algebra）上的矩阵乘法。记这种乘法为 "\(\cdot\)"，则

\(L^{(0)}\)（对角线为 0、其余为 \(\infty\)）扮演单位矩阵的角色（习题 25.1-3）。逐次相乘 \(n-2\) 次的 SLOW-ALL-PAIRS-SHORTEST-PATHS 用时 \(\Theta(n^4)\)。

推导拆解：把普通乘法和 \((\min,+)\) 乘法并排写，取 \(A=\begin{pmatrix}0&1\\4&0\end{pmatrix}\)，算 \(A\cdot A\) 的左下角元素（从 2 到 1）。 普通乘法：\(a_{21}a_{11}+a_{22}a_{21}=4\cdot0+0\cdot4=0\)，这个数没有意义。 \((\min,+)\) 乘法：\(\min(a_{21}+a_{11},\ a_{22}+a_{21})=\min(4+0,\ 0+4)=4\)，含义是"两步以内从 2 到 1 的最低代价"。 两者的骨架完全一样（遍历中间下标 \(k\)，把两个元素结合，再汇总），只是"结合"用加法、"汇总"用取最小。单位元也要跟着换：普通乘法里"什么都不加"是 0，取最小时"什么都没有"是 \(\infty\)，因为 \(\min(x,\infty)=x\)。

例（原书图 25.1）　5 个顶点，边：\(1\to2=3\)，\(1\to3=8\)，\(1\to5=-4\)，\(2\to4=1\)，\(2\to5=7\)，\(3\to2=4\)，\(4\to1=2\)，\(4\to3=-5\)，\(5\to4=6\)。最终

且 \(L^{(5)}=L^{(4)}\cdot W=L^{(4)}\)，此后不再变化。

25.2.3 重复平方

我们只需要 \(L^{(n-1)}\)，并且对所有 \(m\ge n-1\) 都有 \(L^{(m)}=L^{(n-1)}\)。\((\min,+)\) 乘法满足结合律（习题 25.1-4），所以可以像快速幂那样计算 \(W,W^2,W^4,W^8,\dots\)，只要 \(\lceil\lg(n-1)\rceil\) 次乘法：

def faster_all_pairs_shortest_paths(W):     # Θ(n^3 lg n)
    n = len(W)
    L, m = W, 1
    while m < n - 1:
        L = extend_shortest_paths(L, L)     # L^(2m) = (L^(m))^2
        m *= 2
    return L

最后算出的 \(L^{(2m)}\) 满足 \(n-1\le2m<2n-2\)，由 (25.3) 它等于 \(L^{(n-1)}\)。只用两个矩阵交替存储，空间是 \(\Theta(n^2)\)（习题 25.1-8）。

金融直觉：重复平方和计算 \((1+r)^{16}\) 的方法相同：不必乘 16 次，先算 \((1+r)^2\)，平方得 4 次方，再平方得 8 次方、16 次方，4 次乘法就够。这里能"多算"是安全的：\(n=6\) 时需要 \(L^{(5)}\)，重复平方会算到 \(L^{(8)}\)，但 (25.3) 保证 \(L^{(8)}=L^{(5)}\)，多走的步数不会让答案变化。这一切依赖结合律：\(W^4=(W^2)(W^2)\) 必须等于 \(W\cdot W\cdot W\cdot W\)。

负权环检测（习题 25.1-9）：再多乘一次，若结果变了，或者对角线出现负值，就说明有负权环。

这个视角的价值不只是一个 \(\Theta(n^3\lg n)\) 的算法。它说明：凡是形如"\(\min_k\{a_{ik}+b_{kj}\}\)"的递推都是矩阵乘法，因而可以用结合律做重复平方、可以向量化、可以并行。多期最优执行、库存控制之类的动态规划常常可以写成 \((\min,+)\) 或 \((\max,+)\) 的矩阵乘积，长期限时就能用重复平方加速。

25.3 Floyd-Warshall 算法

25.3.1 按"中间顶点"做动态规划

Floyd-Warshall 换了一种刻画最优子结构的方式。简单路径 \(p=\langle v_1,\dots,v_l\rangle\) 的中间顶点（intermediate vertex）是除两端以外的顶点。考虑所有"中间顶点都取自 \(\{1,\dots,k\}\)"的 \(i\to j\) 路径，设 \(p\) 是其中权重最小的一条：

• 若 \(k\) 不是 \(p\) 的中间顶点，则 \(p\) 的中间顶点都在 \(\{1,\dots,k-1\}\) 中，\(p\) 也是"中间顶点限于 \(\{1..k-1\}\)"时的最短路径。
• 若 \(k\) 是中间顶点，把 \(p\) 拆成 \(i\overset{p_1}{\leadsto}k\overset{p_2}{\leadsto}j\)。\(k\) 在简单路径上只出现一次，所以 \(p_1\)、\(p_2\) 的中间顶点都在 \(\{1..k-1\}\) 中，并且它们分别是该限制下 \(i\to k\)、\(k\to j\) 的最短路径。

于是得到递推式，\(d_{ij}^{(k)}\) 表示中间顶点都在 \(\{1..k\}\) 中的 \(i\to j\) 最短路径权重：

\(k=0\) 时路径没有中间顶点，至多一条边。\(D^{(n)}\) 就是答案。

白话解释：Floyd-Warshall 的思路是"逐步开放中转站"。第 0 轮：不许中转，只能直接兑换，\(D^{(0)}=W\)。第 1 轮：允许经过顶点 1 中转，对每一对 \((i,j)\) 只问一个问题："先到 1、再从 1 到 \(j\)，是否比现在的方案便宜？"第 2 轮：再允许经过顶点 2 中转（顶点 1 仍然可用，因为 \(d_{i2}^{(1)}\)、\(d_{2j}^{(1)}\) 里已经包含了经过 1 的方案）……第 \(n\) 轮后，所有顶点都可以当中转站，表中就是真正的最短距离。 例：USD、EUR、CNH 三种货币，先开放 USD 作为中转站。若 EUR→CNH 直接成本 20 bp，EUR→USD 成本 1 bp、USD→CNH 成本 1 bp，则第一轮就把 EUR→CNH 更新为 2 bp。这正是第 24 章 24.9.5 节"EUR→CNH 绕道美元更划算"的结论。 和 25.2 的区别：25.2 每一轮都要在 \(n\) 个中转点里找最好的一个；这里每一轮只考虑一个新中转点，所以每轮工作量少一个 \(n\) 倍。

def floyd_warshall(W):
    n = len(W)
    D = [row[:] for row in W]               # D^(0) = W
    for k in range(n):
        for i in range(n):
            for j in range(n):
                if D[i][k] + D[k][j] < D[i][j]:
                    D[i][j] = D[i][k] + D[k][j]
    return D

时间 \(\Theta(n^3)\)，三重循环，没有复杂的数据结构，常数很小，中等规模的图上非常实用。原书为每个 \(k\) 新建一个矩阵，空间 \(\Theta(n^3)\)；习题 25.2-4 证明去掉上标、原地更新仍然正确，空间 \(\Theta(n^2)\)。原因是：第 \(k\) 轮中 \(d_{ik}\) 和 \(d_{kj}\) 本身不会被改变——更新它们需要 \(d_{ik}+d_{kk}<d_{ik}\)，而 \(d_{kk}=0\)。

原书图 25.4 给出图 25.1 上的 \(D^{(0)},\dots,D^{(5)}\)，\(D^{(5)}\) 就是上面的 \(L^{(4)}\)。

两种动态规划的对比：矩阵乘法版本按"路径边数"分层，每层要在 \(n\) 个中转点里取最小，一层 \(\Theta(n^3)\)、共 \(\lg n\) 层；Floyd-Warshall 按"允许的中间顶点集合"分层，每层只需判断"经过 \(k\) 是否更好"，一层 \(\Theta(n^2)\)、共 \(n\) 层。后者的子问题定义更巧妙，省掉了一个 \(\lg n\) 因子。

25.3.2 构造最短路径

有三种办法：

1. 先算出 \(D\)，再在 \(O(n^3)\) 内推出 \(\Pi\)（习题 25.1-6）。
2. 与 \(D\) 同步计算 \(\Pi^{(k)}\)，\(\pi_{ij}^{(k)}\) 是中间顶点限于 \(\{1..k\}\) 的最短路径上 \(j\) 的前驱：

走经过 \(k\) 的路径时，\(j\) 的前驱沿用 \(k\to j\) 那一段中 \(j\) 的前驱。

3. 记录最大编号的中间顶点 \(\phi_{ij}^{(k)}\)，像矩阵链乘法的 \(s\) 表一样递归重建（习题 25.2-7）。

25.3.3 负权环检测

Floyd-Warshall 在有负权环时不会出错退出，但会在对角线上留下痕迹：\(d_{ii}<0\) 当且仅当顶点 \(i\) 处在某个负权环上（习题 25.2-6）。\(d_{ii}\) 是从 \(i\) 出发回到 \(i\) 的最短"环"的权重，正常情况下为 0。

在汇率图里，这意味着一次 \(\Theta(n^3)\) 的计算就能同时标出所有处在某个套利环上的资产，比第 24 章只返回"一个"负环的信息更全面。

金融直觉：\(d_{ii}\) 是"拿 1 单位资产 \(i\) 出发，允许任意兑换，最后换回资产 \(i\) 的最低成本"。\(e^{-d_{ii}}\) 就是能换回多少单位：正常市场里最好的做法是"什么都不做"，\(d_{ii}=0\)，换回 1 单位；一旦 \(d_{ii}<0\)，说明有一条路能换回多于 1 单位。需要提醒的是，有负权环时表中其他格子的数值已经失去"最短距离"的含义（可以无限下降），只能用来判断"有没有、在哪里"，不能当成可执行的兑换率。

25.3.4 传递闭包

有向图的传递闭包（transitive closure）\(G^*=(V,E^*)\)，\(E^*=\{(i,j):G\text{ 中有 }i\text{ 到 }j\text{ 的路径}\}\)，回答"谁能到达谁"。

• 方法一：所有边权设为 1，跑 Floyd-Warshall；\(d_{ij}<n\) 表示有路径。
• 方法二：把 \(\min\) 和 \(+\) 换成逻辑或 \(\vee\) 和逻辑与 \(\wedge\)：

同样 \(\Theta(n^3)\)，但单比特运算更快、更省空间，还可以位并行（一次处理一整行）。原书图 25.5 的例子：4 个顶点，边 \(2\to3,2\to4,3\to2,4\to1,4\to3\)，最终 \(T^{(4)}\) 的第 1 行为 \((1,0,0,0)\)，其余三行全为 1。稀疏图上也可以从每个顶点做一次 BFS/DFS，\(O(VE)\)（习题 25.2-8）。

背后的统一框架是闭半环（closed semiring）：\((\min,+)\) 给出最短路径，\((\vee,\wedge)\) 给出可达性，\((\max,\min)\) 给出瓶颈路径（最大可通过容量的路径）。Floyd-Warshall 是这些问题共同的算法骨架。

25.4 Johnson 算法：稀疏图上的重新赋权

25.4.1 思路

若所有边权非负，对每个顶点跑一次 Dijkstra（斐波那契堆）就是 \(O(V^2\lg V+VE)\)，稀疏图上比 \(\Theta(V^3)\) 好得多。有负权边怎么办？Johnson 的办法是重新赋权（reweighting）：构造新的权函数 \(\hat w\)，满足

1. 对所有 \(u,v\)，\(p\) 是 \(w\) 下的最短路径 \(\iff\) \(p\) 是 \(\hat w\) 下的最短路径；
2. 对所有边，\(\hat w(u,v)\ge0\)。

25.4.2 势函数不改变最短路径

引理 25.1（重新赋权不改变最短路径）　对任意函数 \(h:V\to\mathbb R\)，定义

则对任意路径 \(p=\langle v_0,\dots,v_k\rangle\)，\(w(p)=\delta(v_0,v_k)\iff\hat w(p)=\hat\delta(v_0,v_k)\)；并且 \(G\) 在 \(w\) 下有负权环 \(\iff\) 在 \(\hat w\) 下有负权环。

证明：望远镜求和，

\(h(v_0)-h(v_k)\) 只取决于起点和终点，与走哪条路无关，所以同一对端点之间各路径的大小关系不变。对环 \(v_0=v_k\)，\(\hat w(c)=w(c)\)，负权环性质不变。\(\square\)

推导拆解：用原书例子的数字验证 (25.10)。路径 \(3\to2\to4\to1\)，原权重 \(4+1+2=7\)；\(h=(0,-1,-5,0,-4)\)。重新赋权后三条边是 \(\hat w(3,2)=0\)，\(\hat w(2,4)=0\)，\(\hat w(4,1)=2\)，合计 2。按 (25.10)：\(w(p)+h(3)-h(1)=7+(-5)-0=2\)，一致。 为什么不改变"哪条路最短"：从 3 到 1 的任何一条路径，重新赋权后都统一加上同一个数 \(h(3)-h(1)=-5\)。所有候选方案同加一个常数，排名不变。这和比较几只债券的收益率时统一减去同一个基准利率、排序不变是一个道理。

常见误区：不能简单地把所有边权都加上同一个常数 \(-w^*\)（\(w^*\) 为最小边权）来消除负权。边数不同的路径会被加上不同的总量，最短路径会改变（习题 25.3-4）。势函数之所以可行，是因为它的修正量只取决于两个端点。

推导拆解：最小的反例。从 \(a\) 到 \(c\) 有两条路：直接 \(a\to c\) 权 \(-0.5\)；或 \(a\to b\to c\)，权 \(-1\) 和 \(0\)，合计 \(-1\)，两步路径最短。最小边权 \(w^*=-1\)，给每条边加 1：直接路径变成 \(0.5\)，两步路径变成 \(0+1=1\)，最短路径反转了。原因是两步路径被加了两次常数，"步数多"的路径受到额外惩罚。势函数的修正只看起点和终点，不看中间走了几步。

25.4.3 怎样找到让所有边非负的 \(h\)

新建图 \(G'\)：加一个新顶点 \(s\)，到每个顶点连一条权 0 的边。\(s\) 没有入边，所以除了以 \(s\) 为起点的路径外，任何最短路径都不经过它；\(G'\) 无负权环 \(\iff\) \(G\) 无负权环。

令 \(h(v)=\delta(s,v)\)（用 Bellman-Ford 求）。由三角不等式 \(h(v)\le h(u)+w(u,v)\)，即

这和第 24 章差分约束的定理 24.9 是同一件事：\(h\) 是差分约束系统 \(h(v)-h(u)\le w(u,v)\) 的一个可行解。

白话解释：为什么 \(h(v)=\delta(s,v)\) 恰好让所有新权重非负？三角不等式说"到 \(v\) 的最短距离，不会比'先到 \(u\) 再走一步 \((u,v)\)'更长"，即 \(h(v)\le h(u)+w(u,v)\)，移项就是 \(\hat w(u,v)\ge0\)。直观地讲，新权重衡量"走这条边比最优走法多花了多少"，最优走法本身多花 0，任何边都不会比最优还省。树上的边（最短路径用到的边）\(\hat w=0\)，正好对应原书例子里那些为 0 的 \(\hat w\)。

例（原书图 25.6，用图 25.1 的图）　\(h(1..5)=(0,-1,-5,0,-4)\)。重新赋权后：\(\hat w(1,2)=4\)，\(\hat w(1,3)=13\)，\(\hat w(1,5)=0\)，\(\hat w(2,4)=0\)，\(\hat w(2,5)=10\)，\(\hat w(3,2)=0\)，\(\hat w(4,1)=2\)，\(\hat w(4,3)=0\)，\(\hat w(5,4)=2\)，全部非负。再以每个顶点为源跑 Dijkstra，最后用

还原原权重下的距离。

25.4.4 算法与复杂度

def johnson(G, w):
    # 1. 加超级源 s，到每个顶点一条 0 权边
    Gp = G.copy(); s = Gp.add_vertex()
    for v in G.V:
        Gp.add_edge(s, v, 0)
    # 2. Bellman-Ford 求 h，并检测负环           O(VE)
    if not bellman_ford(Gp, w, s):
        raise ValueError("the input graph contains a negative-weight cycle")
    h = {v: v.d for v in Gp.V}
    # 3. 重新赋权                                O(E)
    w_hat = {(u, v): w(u, v) + h[u] - h[v] for (u, v) in Gp.E}
    # 4. 每个顶点跑一次 Dijkstra                 |V| × O(V lg V + E)
    D = {}
    for u in G.V:
        delta_hat = dijkstra(G, w_hat, u)
        for v in G.V:
            D[u, v] = delta_hat[v] + h[v] - h[u]   # 还原原权重
    return D

斐波那契堆实现 \(O(V^2\lg V+VE)\)；二叉堆 \(O(VE\lg V)\)，稀疏图上仍渐近快于 Floyd-Warshall。

超级源点必不可少：如果随便拿一个已有顶点当 \(s\)，可能有顶点从它不可达，\(h\) 为 \(\infty\)（习题 25.3-6）；只有图强连通时才可以省掉。若原图所有边权都非负，\(h\equiv0\)，\(\hat w=w\)（习题 25.3-3）。

25.4.5 金融解读：势函数就是价格

把第 24 章的汇率图带进来：\(w(i,j)=-\ln R[i,j]\)。若 \(h(i)\) 取货币 \(i\) 的对数价值，约化权重 \(\hat w(i,j)=-\ln R[i,j]+h(i)-h(j)\) 就是"按公允价值计，这笔兑换亏了多少"——也就是点差和费用。无套利时所有 \(\hat w\ge0\)；Johnson 算法的第 2 步，实质上就是从报价中反推出一组自洽的价格。势函数重新赋权在最小费用流（订单路由、资金调拨）中也常用：用它把负成本边变成非负，从而可以反复使用 Dijkstra。第 29 章会看到，这组 \(h\) 正是对偶线性规划的最优解。

25.5 思考题与注记

• 25-1 动态图的传递闭包：边逐条插入，用布尔矩阵维护传递闭包。每插入一条边可以 \(O(V^2)\) 更新；有的插入必须 \(\Omega(V^2)\)；但只对"新变成可达"的顶点对做传播（每对至多由 0 变 1 一次），任意 \(n\) 次插入的总时间可做到 \(O(V^3)\)——又是一个摊还论证。
• 25-2 ε-稠密图：\(|E|=\Theta(V^{1+\epsilon})\)。用 \(d\) 叉堆（INSERT、DECREASE-KEY \(O(\log_dn)\)，EXTRACT-MIN \(O(d\log_dn)\)），取 \(d=V^\epsilon\)，非负权单源最短路径可做到 \(O(E)\)，全源 \(O(VE)\)；结合 Johnson 重新赋权，有负权无负环时全源也是 \(O(VE)\)。这是工程中替代斐波那契堆的实用选择。

注记：Floyd-Warshall 来自 Floyd，基于 Warshall 关于布尔矩阵传递闭包的定理。后续改进有 Fredman 的 \(O(V^3(\lg\lg V/\lg V)^{1/3})\)、Han 的 \(O(V^3(\lg\lg V/\lg V)^{5/4})\)；在无向无权图上可借助快速矩阵乘法做到 \(O(V^\omega p(V))\)（\(\omega<2.376\)，\(p\) 为多对数因子）。Aho、Hopcroft、Ullman 定义的闭半环给出了有向图路径问题的统一代数框架。

25.6 量化实战

25.6.1 先核对原书例题

下面用 numpy 实现三个算法：\((\min,+)\) 矩阵乘法用广播一次算完；Floyd-Warshall 保留 \(k\) 的外层循环、把 \(i,j\) 两层向量化；Johnson 用 Bellman-Ford 求 \(h\)、重新赋权后跑 Dijkstra。

import numpy as np
import heapq

INF = np.inf
# 图 25.1（顶点 1..5 存为下标 0..4）
W = np.full((5, 5), INF); np.fill_diagonal(W, 0)
for u, v, w in [(1,2,3),(1,3,8),(1,5,-4),(2,4,1),(2,5,7),(3,2,4),(4,1,2),(4,3,-5),(5,4,6)]:
    W[u-1, v-1] = w

def minplus(A, B):
    """(min,+) 矩阵乘法：C[i,j] = min_k A[i,k] + B[k,j]，用广播一次算完。"""
    return np.min(A[:, :, None] + B[None, :, :], axis=1)

def faster_apsp(W):                       # 重复平方，Θ(n^3 lg n)
    L, m, n = W.copy(), 1, len(W)
    while m < n - 1:
        L = minplus(L, L); m *= 2
    return L

def floyd_warshall(W):                    # Θ(n^3)，k 循环串行，i、j 两层向量化
    D = W.copy(); n = len(W)
    P = np.where(np.isfinite(W) & ~np.eye(n, dtype=bool), np.arange(n)[:, None], -1)
    for k in range(n):
        via = D[:, [k]] + D[[k], :]       # d_ik + d_kj
        better = via < D
        D = np.where(better, via, D)
        P = np.where(better, P[[k], :], P)  # 式 (25.7)：前驱取 pi_kj
    return D, P

def path(P, i, j):
    if i == j: return [i]
    if P[i, j] < 0: return None
    return path(P, i, P[i, j]) + [j]

L = faster_apsp(W)
D, P = floyd_warshall(W)
print("重复平方 L:\n", L.astype(int))
print("与 Floyd-Warshall 一致:", np.array_equal(L, D))
print("顶点 3 到 1 的最短路径:", [int(v) + 1 for v in path(P, 2, 0)], "权重", int(D[2, 0]))

# Johnson：超级源点 + Bellman-Ford 求 h，重新赋权后每个顶点跑 Dijkstra
n = len(W)
h = np.zeros(n)                            # 超级源点到各点的 0 权边 => 初值全 0
for _ in range(n):
    for u in range(n):
        for v in range(n):
            if u != v and np.isfinite(W[u, v]) and h[u] + W[u, v] < h[v]:
                h[v] = h[u] + W[u, v]
print("h =", h.astype(int))
What = W + h[:, None] - h[None, :]
print("重新赋权后的边:", {(u+1, v+1): int(What[u, v]) for u in range(n) for v in range(n)
                    if u != v and np.isfinite(W[u, v])})

def dijkstra_dense(Wh, s):
    d = np.full(len(Wh), INF); d[s] = 0; pq = [(0.0, s)]; done = set()
    while pq:
        du, u = heapq.heappop(pq)
        if u in done: continue
        done.add(u)
        for v in np.nonzero(np.isfinite(Wh[u]))[0]:
            if du + Wh[u, v] < d[v]:
                d[v] = du + Wh[u, v]; heapq.heappush(pq, (d[v], v))
    return d

DJ = np.array([dijkstra_dense(What, u) for u in range(n)]) - h[:, None] + h[None, :]
print("Johnson 与 Floyd-Warshall 一致:", np.allclose(DJ, D))

输出：

重复平方 L:
 [[ 0  1 -3  2 -4]
 [ 3  0 -4  1 -1]
 [ 7  4  0  5  3]
 [ 2 -1 -5  0 -2]
 [ 8  5  1  6  0]]
与 Floyd-Warshall 一致: True
顶点 3 到 1 的最短路径: [3, 2, 4, 1] 权重 7
h = [ 0 -1 -5  0 -4]
重新赋权后的边: {(1, 2): 4, (1, 3): 13, (1, 5): 0, (2, 4): 0, (2, 5): 10, (3, 2): 0, (4, 1): 2, (4, 3): 0, (5, 4): 2}
Johnson 与 Floyd-Warshall 一致: True

\(L^{(4)}\)、\(h\) 和全部 \(\hat w\) 都与原书一致。路径 \(3\to2\to4\to1\) 的权重 \(4+1+2=7\)，与 \(D\) 的第 3 行第 1 列相符。

向量化的要点：\(k\) 循环必须串行，因为第 \(k\) 轮要用第 \(k-1\) 轮的结果；而固定 \(k\) 后，所有 \((i,j)\) 的更新互相独立（\(d_{ik}\)、\(d_{kj}\) 本轮不变），可以一次用广播完成。这正是第 27 章习题 27.2-6 多线程 Floyd-Warshall 的结构：工作量 \(\Theta(n^3)\)，跨度 \(\Theta(n\lg n)\)。

25.6.2 多资产最优兑换表与套利定位

在数字资产市场，币对很多但盘口稀疏：几乎每个资产都对 USD、USDT、BTC、ETH 有盘口，山寨币之间的直接盘口很少且点差很宽。下面模拟 60 种资产，用 Floyd-Warshall 一次算出任意两种资产之间的最优可得兑换率 \(e^{-d_{ij}}\)，再看注入错价后对角线的变化。

import time
import numpy as np
from scipy.sparse.csgraph import floyd_warshall as sp_fw, NegativeCycleError

rng = np.random.default_rng(7)
n = 60                                     # 60 种数字资产
names = ["USD", "USDT", "BTC", "ETH"] + [f"A{k:02d}" for k in range(n - 4)]
logv = np.r_[0.0, 0.0, np.log(60000), np.log(3000), rng.normal(0, 2, n - 4)]
# 只有部分币对有市场：每个资产都对 USD/USDT/BTC/ETH 有盘口，另随机 3% 的山寨对山寨盘口
has = np.zeros((n, n), bool); has[:, :4] = has[:4, :] = True
has |= rng.random((n, n)) < 0.03; has = has | has.T; np.fill_diagonal(has, False)
hs = np.where(np.arange(n)[:, None] < 4, 2e-4, 15e-4) + np.where(np.arange(n)[None, :] < 4, 2e-4, 15e-4)
hs = hs * rng.uniform(0.5, 1.5, (n, n)); hs = (hs + hs.T) / 2       # 半点差，主流对窄、山寨对宽
R = np.where(has, np.exp(logv[:, None] - logv[None, :]) * (1 - hs), 0.0)

def fw_numpy(W):
    D = W.copy()
    for k in range(len(W)):
        D = np.minimum(D, D[:, [k]] + D[[k], :])
    return D

W = np.where(has, -np.log(np.where(has, R, 1.0)), np.inf); np.fill_diagonal(W, 0)
t0 = time.perf_counter(); D = fw_numpy(W); t1 = time.perf_counter()
print(f"numpy Floyd-Warshall n={n}: {1e3*(t1-t0):.1f} ms，对角线最小值 {D.diagonal().min():.2e}（无负环）")
D_sp = sp_fw(W, directed=True)
print("与 scipy 结果一致:", np.allclose(D, D_sp))

best = np.exp(-D)                          # 任意两资产之间的最优可得兑换率
direct = np.where(has, R, np.nan)
gain_bp = (best / direct - 1) * 1e4        # 多跳路由相对直接兑换的改善
m = has & (gain_bp > 0.01)
print(f"有直接盘口的 {has.sum()} 个有向币对中，{m.sum()} 个走多跳更划算，"
      f"中位改善 {np.median(gain_bp[m]):.1f} bp，最大 {np.nanmax(gain_bp):.1f} bp")

# 注入一个错价，看对角线
R2 = R.copy()
a, b = names.index("A07"), names.index("BTC")
R2[a, b] *= 1.004; R2[b, a] /= 1.004       # A07/BTC 盘口整体偏离 40bp
W2 = np.where(has, -np.log(np.where(has, R2, 1.0)), np.inf); np.fill_diagonal(W2, 0)
D2 = fw_numpy(W2)
neg = np.nonzero(D2.diagonal() < -1e-12)[0]
print("对角线为负（处在某个套利环上）的资产:", [names[i] for i in neg])
try:
    sp_fw(W2, directed=True)
except NegativeCycleError as e:
    print("scipy 报告:", type(e).__name__)

# 规模测试：纯 Python 三重循环 vs numpy 按 k 向量化
def fw_python(W):
    n = len(W); D = [list(r) for r in W]
    for k in range(n):
        Dk = D[k]
        for i in range(n):
            dik = D[i][k]; Di = D[i]
            for j in range(n):
                if dik + Dk[j] < Di[j]:
                    Di[j] = dik + Dk[j]
    return D
for size in (100, 200, 400):
    Wb = np.where(rng.random((size, size)) < 0.1, rng.uniform(1, 10, (size, size)), np.inf)
    np.fill_diagonal(Wb, 0)
    t0 = time.perf_counter(); fw_python(Wb.tolist()); t1 = time.perf_counter(); fw_numpy(Wb); t2 = time.perf_counter()
    print(f"n={size}: 纯 Python {t1-t0:.2f} s，numpy {1e3*(t2-t1):.1f} ms")

输出（计时随机器而异）：

numpy Floyd-Warshall n=60: 0.3 ms，对角线最小值 0.00e+00（无负环）
与 scipy 结果一致: True
有直接盘口的 628 个有向币对中，268 个走多跳更划算，中位改善 3.4 bp，最大 15.5 bp
对角线为负（处在某个套利环上）的资产: ['USD', 'USDT', 'BTC', 'A07']
scipy 报告: NegativeCycleError
n=100: 纯 Python 0.01 s，numpy 0.9 ms
n=200: 纯 Python 0.10 s，numpy 5.1 ms
n=400: 纯 Python 0.96 s，numpy 35.5 ms

几点解读：

• 最优兑换表：628 个有直接盘口的有向币对里，有 268 个"绕道"更划算，最多能省 15.5 bp。这在山寨币之间尤其明显——直接盘口点差宽，而经 USDT 或 BTC 两跳的点差之和更窄。这张表就是一个智能路由器的核心数据结构。
• 对角线定位：A07/BTC 的盘口偏离 40 bp 后，对角线为负的资产恰好是 USD、USDT、BTC、A07——它们都在某个包含错价边的套利环上（例如 A07→BTC→USD→A07）。注意 ETH 不在其中：环 A07→BTC→ETH→A07 的权重为 +2.4 bp（沿途点差略大于 40 bp 的错价），而 A07→BTC→USD→A07 为 −2.6 bp。套利是否存在，取决于错价与整条环上点差之和的比较。scipy.sparse.csgraph.floyd_warshall 在有负环时直接抛出 NegativeCycleError，不告诉你在哪；自己实现时读对角线信息更丰富。
• 向量化的收益：同样 \(\Theta(n^3)\) 的算法，按 \(k\) 向量化后比纯 Python 快约 20–30 倍。对几百种资产，numpy 版本在几十毫秒内完成，足以在每次行情快照后重算。

25.6.3 传递闭包：因子流水线的依赖分析

量化研究平台里，因子之间、数据表之间存在大量依赖。某个原始数据源被修订后，哪些下游结果需要重算？这就是传递闭包问题。

import numpy as np
# 日终因子流水线的依赖关系：边 (a, b) 表示 b 的计算要用到 a
nodes = ["行情原始", "成交原始", "复权价", "日收益", "波动率20", "动量60", "流动性", "综合打分", "风险模型"]
idx = {v: i for i, v in enumerate(nodes)}
deps = [("行情原始", "复权价"), ("复权价", "日收益"), ("日收益", "波动率20"), ("日收益", "动量60"),
        ("成交原始", "流动性"), ("复权价", "流动性"), ("动量60", "综合打分"), ("流动性", "综合打分"),
        ("波动率20", "风险模型"), ("日收益", "风险模型")]
n = len(nodes)
T = np.eye(n, dtype=bool)
for a, b in deps:
    T[idx[a], idx[b]] = True
for k in range(n):                           # 式 (25.8)：t_ij |= t_ik & t_kj，按 k 向量化
    T |= T[:, [k]] & T[[k], :]
for src in ["成交原始", "行情原始"]:
    print(f"{src} 修订后需重算:", [nodes[j] for j in np.nonzero(T[idx[src]])[0] if nodes[j] != src])
print("综合打分 依赖的全部上游:", [nodes[i] for i in np.nonzero(T[:, idx["综合打分"]])[0] if nodes[i] != "综合打分"])

输出：

成交原始 修订后需重算: ['流动性', '综合打分']
行情原始 修订后需重算: ['复权价', '日收益', '波动率20', '动量60', '流动性', '综合打分', '风险模型']
综合打分 依赖的全部上游: ['行情原始', '成交原始', '复权价', '日收益', '动量60', '流动性']

行读出"下游影响面"，列读出"上游血缘"。真实平台节点数可达数千，此时更适合用每个节点做一次 DFS（\(O(VE)\)，习题 25.2-8），或像思考题 25-1 那样增量维护。

本章小结

全源最短路径有三条路线。第一，把"再走一条边"的递推看作 \((\min,+)\) 半环上的矩阵乘法，朴素 \(\Theta(n^4)\)，重复平方 \(\Theta(n^3\lg n)\)。第二，Floyd-Warshall 按"允许使用的中间顶点集合 \(\{1..k\}\)"做动态规划，\(\Theta(n^3)\) 时间、\(\Theta(n^2)\) 空间，代码最短；把运算换成布尔或/与就得到传递闭包。第三，Johnson 算法用势函数 \(h(v)=\delta(s,v)\) 重新赋权，把负权变非负，再跑 \(|V|\) 次 Dijkstra，\(O(V^2\lg V+VE)\)，适合稀疏图。重新赋权的核心等式 \(\hat w(p)=w(p)+h(v_0)-h(v_k)\) 说明它不改变最短路径，也不改变任何环的权重。负权环检测：Floyd-Warshall 看对角线是否为负，Johnson 由 Bellman-Ford 报告。在汇率图上，Floyd-Warshall 一次给出最优兑换表并标出所有处于套利环上的资产，势函数 \(h\) 则对应一组自洽的价格。

练习

基础

1. 在原书图 25.2（6 个顶点）上手算 FASTER-ALL-PAIRS-SHORTEST-PATHS 和 Floyd-Warshall，写出每步的矩阵。（原书 25.1-1、25.2-1。）
2. 为什么要求 \(w_{ii}=0\)？如果 \(w_{ii}\) 取 \(\infty\)，(25.2) 中的 \(l_{ij}^{(m)}\) 表示什么？（原书 25.1-2。提示：变成"恰好 \(m\) 条边"。）
3. 证明 \((\min,+)\) 矩阵乘法满足结合律。（原书 25.1-4。）
4. 证明原地更新的 Floyd-Warshall 正确。（原书 25.2-4。提示：第 \(k\) 轮中第 \(k\) 行、第 \(k\) 列不变。）
5. 说明如何用 Floyd-Warshall 的输出检测负权环。（原书 25.2-6。）
6. 给出"把所有边权加上 \(-w^*\)"会改变最短路径的反例。（原书 25.3-4。提示：一条两边的路径与一条单边路径。）

进阶

7. 修改 FASTER-ALL-PAIRS-SHORTEST-PATHS 使其能检测负权环；再设计一个算法求边数最少的负权环的边数。（原书 25.1-9、25.1-10。）
8. 证明：若 \(c\) 是零权环，则 Johnson 重新赋权后 \(c\) 上每条边的 \(\hat w\) 都为 0。（原书 25.3-5。）
9. 瓶颈路径：边上的数是通道容量，路径的容量是沿途容量的最小值。改写 Floyd-Warshall 求任意两点间的最大容量路径（"\((\max,\min)\) 半环"），并用它回答"从交易所 A 到交易所 B 一次最多能转多少资金"。
10. 在 25.6.2 的代码中，把"最优兑换率"改成考虑每一跳固定手续费 \(f\)（按比例）。注意多跳路径要付多次费用。比较 \(f=0,5,10\) bp 时"走多跳更划算"的币对数量。

原书推荐习题：25.1-4、25.1-8、25.1-9、25.2-3、25.2-4、25.2-6、25.3-4、25.3-6，思考题 25-2。

原书对照

页码换算：原书页码 = PDF 页码 − 21。