STOC 2025 的这篇论文《Breaking the Sorting Barrier for Directed Single-Source Shortest Paths》将单源最短路径问题的时间复杂度做到了 O ( m log 2 3 n ) O(m\log^{\frac{2}{3}}n) O ( m log 3 2 n ) ,在稀疏图上优于 Dijkstra 的 O ( m + n log n ) O(m + n\log n) O ( m + n log n ) 。
1. 假设
论文首先做了两个假设:
图上每个点的入度和出度均为常数,比如 2 2 2 。这样就有 m = Θ ( n ) m=\Theta(n) m = Θ ( n ) 。这样的假设可以保证后续一个点向外松弛时时间复杂度不会炸(最多松弛两次)
图上每个点的最短路长度都是不同的。这样的假设可以保证每个点的最短路前驱最多只会有一个,最短路关系一定形成的是最短路树而不是 DAG。
在最后我们会解释这两个假设是不失一般性的。
2. 算法过程
算法有两个参数,k k k 和 t t t ,并且设 k = ⌊ log 1 3 n ⌋ k = \left \lfloor \log^{\frac{1}{3}}n \right \rfloor k = ⌊ log 3 1 n ⌋ ,t = ⌊ log 2 3 n ⌋ t = \left \lfloor \log^{\frac{2}{3}}n \right \rfloor t = ⌊ log 3 2 n ⌋ 。后续我们会看到,这样的两个参数可以得到最后 O ( m log 2 3 n ) O(m\log^{\frac{2}{3}}n) O ( m log 3 2 n ) 的时间复杂度。
论文中为了叙述方便,还引入了如下符号:d ( x ) d(x) d ( x ) 表示 x x x 最终的最短路距离,d ^ ( x ) \hat{d}(x) d ^ ( x ) 表示 x x x 当前距离。
为了方便后面叙述,我们称一个节点 complete ,就表示这个节点的最短路已经找到了。一个 complete 的点即满足 d ( x ) = d ^ ( x ) d(x) = \hat{d}(x) d ( x ) = d ^ ( x ) 。
算法的主体部分是 BMSSP 算法,BMSSP 又用到了两个部分:Finding Pivots 和一个神秘的数据结构。对于后者,我们在这里将其称之为分块优先队列。
2.1 BMSSP
BMSSP 算法是整篇论文的主体。BMSSP 算法接收三个参数 l l l ,S S S 和 B B B 。
l l l 表示当前层数,BMSSP 的过程有点类似于分治,有一个递归过程,而在不断递归的时候,BMSSP 的层数会一直降低,直到 l = 0 l=0 l = 0 ,我们就说到达了一个 BaseCase。在 BaseCase 时,会走另一个算法,我们一会再说。
S S S 表示一个集合。
B B B 表示一个上界。
BMSSP 仅考虑最短路长度小于 B B B ,且没有 complete 的点。我们将这种点称为目标节点。这些节点组成集合 U ~ \widetilde{U} U 。
同时 BMSSP 还对上述参数做了如下要求:
S S S 不能太大:∣ S ∣ ≤ 2 l t |S| \le 2^{lt} ∣ S ∣ ≤ 2 lt 。
你需要保证 U ~ \widetilde{U} U 中所有的点的最短路,都经过 S S S 中某个已经 complete 的点(一个点的最短路肯定是由某个已经找到最短路的点推过来的)。这也就意味着 S S S 中必须有至少一个 complete 的点(但不意味着 S S S 中的点都是 complete 的,S S S 里面混杂着 complete 和没 complete 的点)(同时后续算法会保证,B B B 一定大于 S S S 中最大的 d ^ \hat{d} d ^ ,这也就意味着 S S S 中没被 complete 的点的最短路一定是经过已经 complete 的点)。
当执行完一次 BMSSP 算法之后,会产生两种情况:
成功执行:BMSSP 会把 U ~ \widetilde{U} U 的点全部 complete,(由于第二个要求,这件事是肯定能办到的)
部分执行:BMSSP 为了控制时间复杂度,只 complete 掉 U ~ \widetilde{U} U 内的一部分点。
然后算法会返回一个集合 U U U 和一个边界 B ′ B' B ′ 。表示说所有最短路长度小于 B ′ B' B ′ ,且没有 complete 的点,都已经被 complete 掉,且放入集合 U U U 中了。
对于成功执行情况,就会有 U = U ~ U=\widetilde{U} U = U ,B ′ = B B'=B B ′ = B 。对于部分执行情况,就会有 U ⊂ U ~ U\subset \widetilde{U} U ⊂ U ,B ′ < B B'<B B ′ < B 。返回这两个值可以让调用方得知 BMSSP 完成的进度。
然后我们要求点 x x x 的单源最短路的话,只需要令 l = ⌈ log n t ⌉ l=\left\lceil \frac{\log n}{t} \right\rceil l = ⌈ t l o g n ⌉ (为什么是这样后面会解释),S = { x } S=\left\{x\right\} S = { x } ,B = ∞ B=\infty B = ∞ ,跑 BMSSP 就可以了。
然后我们来看 BMSSP 是怎么做的:
既然我们已经在 BMSSP 中的第二个要求中说,我们的目标节点最短路一定是经过 S S S 的,那么朴素的想法是从 S S S 出发往外松弛。论文则提出,我们可以先搞一个 Finding Pivots 操作,先把 S S S 的规模缩一下,缩成 P P P ,然后再从 P P P 往外松弛。同时 Finding Pivots 操作在缩小规模的时候还顺便 complete 掉一些点,这些点被放在 W W W 中(但 W W W 里面还混杂着一部分没有 complete 的点,需要后面再筛出来)。
朴素的 Dijkstra 会使用堆,把 S S S 都放入堆中,然后每次把 d ^ \hat{d} d ^ 最小的点拿出来往外松弛(d ^ \hat{d} d ^ 最小的点一定是已经 complete 的,或者说有 d ^ = d \hat{d} = d d ^ = d )。堆的这个排序是 Dijkstra 时间复杂度的瓶颈,所以论文为了避开这点,选择使用一个分块的优先队列来维护 d ^ \hat{d} d ^ 。分块的大小为 M M M 。可以往这个分块优先队列里插入,或是批量插入(批量插入的时间复杂度会小一点,但同时需要满足批量插入的数据比分块优先队列中最小的数据还要小),或是取出前 M M M 个最小的数据。
于是 BMSSP 就用分块优先队列,以类似朴素的 Dijkstra 的方法进行求解。首先把 P P P 都放入这个分块优先队列中,然后每次循环,把分块优先队列中前 M M M 小的点拉出来,记为 S i S_i S i (S i S_i S i 内可能是无序的)。拉出来的时候数据结构还会返回一个上界 B i B_i B i ,B i B_i B i 为分块优先队列中剩余的最小的 d ^ \hat{d} d ^ (也意味着 B i B_i B i 比 S i S_i S i 中最大的 d ^ \hat{d} d ^ 还大)。
由于 S i S_i S i 内的点不是全部 complete 的,所以不能直接往外松弛。BMSSP 就会先拿着 S i S_i S i 和 B i B_i B i ,去调用 l − 1 l-1 l − 1 层的 BMSSP。首先,由于 ∣ S i ∣ ≤ M = 2 ( l − 1 ) t |S_i|\le M = 2^{(l-1)t} ∣ S i ∣ ≤ M = 2 ( l − 1 ) t (当分块优先队列里的数据不足 M M M 个的时候,∣ S i ∣ < M |S_i|< M ∣ S i ∣ < M ),所以满足下一层 BMSSP 的第一个要求。事实上,对于第二个要求也是成立的,即所有最短路长度小于 B i B_i B i 的点的最短路事实上就是经过 S i S_i S i 的(后面正确性部分会细说,这个对应正确性部分的 Pull Minimum 引理)。调用完之后就拿到了子 BMSSP 返回的 U i U_i U i 和 B i ′ B_i' B i ′ 。U i U_i U i 是子 BMSSP complete 掉的点的集合,我们把他直接并入当前 BMSSP 的 U U U 中。
现在 U i U_i U i 已经是 complete 的了,于是从 U i U_i U i 往外松弛。对于边 ( u , v ) (u,v) ( u , v ) ,如果松弛成功的话,根据朴素的 Dijkstra,我们需要把 v v v 放入堆中。我们这里也是一样,只不过我们需要对 d ^ ( v ) \hat{d}(v) d ^ ( v ) 进行讨论:
d ^ ( v ) ≥ B \hat{d}(v)\ge B d ^ ( v ) ≥ B ,这个我们就不管了,因为当前 BMSSP 只负责最短路小于 B B B 的部分。
B i ≤ d ^ ( v ) < B B_i\le \hat{d}(v) < B B i ≤ d ^ ( v ) < B ,我们把 v v v 直接插入分块优先队列中。
B i ′ ≤ d ^ ( v ) < B i B_i'\le \hat{d}(v) <B_i B i ′ ≤ d ^ ( v ) < B i ,我们也把 v v v 插入分块优先队列中。但是注意,我们之前说,分块优先队列中剩余数据的最小值是 B i B_i B i ,这意味着 d ^ ( v ) \hat{d}(v) d ^ ( v ) 比最小值还要小,可以使用批量插入来加速。
d ^ ( v ) < B i ′ \hat{d}(v) < B_i' d ^ ( v ) < B i ′ ,不可能,因为既然子 BMSSP 返回了 B i ′ B_i' B i ′ ,就说明最短路距离小于 B i ′ B_i' B i ′ 的点都被 complete 了。一个已经 complete 的点是不可能再被有效松弛的。
同时,由于子 BMSSP 部分执行的原因,可能 S i S_i S i 有部分点他本身的 d ^ \hat{d} d ^ 就大于了 B i ′ B_i' B i ′ ,没派上用场,所以需要把这部分点也重新放回分块优先队列中。这部分点的 d ^ \hat{d} d ^ 也是小于 B i B_i B i 的,所以可以批量插入。需要批量插入的点都放入集合 K K K 中,在最后统一批量插入。
每次循环开始还会检查当前已经 complete 的点数量,即集合 U U U 的大小,是否超过了 k 2 l t k2^{lt} k 2 lt ,超过的话就说明当前 BMSSP 的工作量有点大了,需要立刻终止。因此,为了保证最顶层的 BMSSP 一定是成功执行的,我们必须让顶层的 l l l 满足 k 2 l t > n k2^{lt}>n k 2 lt > n 。l = ⌈ log n t ⌉ l=\left\lceil \frac{\log n}{t} \right\rceil l = ⌈ t l o g n ⌉ 是满足的,所以我们最顶层的 BMSSP 取了这个值。(这个值偏大,但是不影响时间复杂度)。
在 BMSSP 的最后,会把最后一次子 BMSSP 的返回值 B i ′ B_i' B i ′ 作为返回值 B ′ B' B ′ 。然后这里有一点:W W W 集合中,满足 d ^ \hat{d} d ^ 小于 B ′ B' B ′ 的点就是先前 W W W 中已经 complete 的点(为什么是这样会在后面将 Finding Pivots 的时候说)。把这些点都筛出来放入集合 U U U ,就组成了当前 BMSSP 所有 complete 掉的点。
你可以这么直观地理解 BMSSP:我们的目标是把最短路距离位于 [ 0 , B ) [0, B) [ 0 , B ) 的点全部 complete 掉。前 i − 1 i-1 i − 1 次循环,子 BMSSP 把这个进度推进到了 B i − 1 ′ B_{i-1}' B i − 1 ′ 。第 i i i 次循环,我们希望子 BMSSP 把进度推进到 B i B_i B i ,但是子 BMSSP 可能由于部分执行的原因,只把进度推到了 B i ′ B_i' B i ′ 。当前 BMSSP 也可能由于工作量太大,强行终止,使得最后的进度也没到 B B B ,只到了最后一次子 BMSSP 推进到的进度。
2.2 BaseCase
注意,上面 BMSSP 的分块优先队列 D D D 保证了每次取出的 S i S_i S i 的大小不超过 M M M 。当 l = 1 l=1 l = 1 时,M = 1 M=1 M = 1 ,所以再往 l = 0 l=0 l = 0 这一层传的 S i S_i S i 只有一个元素了。同时根据 BMSSP 的第二个要求,可以得知此时跑一个 Dijkstra 是没问题的:
在 BaseCase 中同样的需要满足 ∣ U ∣ < k 2 l t |U|<k2^{lt} ∣ U ∣ < k 2 lt ,以及只关心最短路小于 B B B 的点。
2.3 分块优先队列
然后我们先讲点数据结构。
对于分块优先队列的直观的感觉已经在上面说了,就是顾名思义,分块的优先队列,把数据一块一块地维护,块内的数据是乱序的,块间是有序的(也就是前一个块的最大值小于后一个块的最小值)。
更具体地说,这个数据结构维护了若干个 key-value,有两个给定的参数 M M M 和 B B B ,数据结构中最多会维护 N N N 个 key-value。然后想要有如下的操作:
Insert :插入一对 key-value,如果 key 已经存在于数据结构中,那么就更新它的 value。所有的 value 值不超过一个边界值 B B B 。
Batch Prepend :如果有若干个 key 的 value 值比数据结构中已有的最小 value 还要小,那么这堆 key 可以批量插入。
Pull :从这个数据结构中取出 value 最小的 M M M 个 key 并从数据结构中删除。如果已有的 key 不足 M M M 个,则把所有的 key 都取出来。同时该操作还返回一个上界 x x x ,表示数据结构中剩余的元素中的最小 value 的值(如果最后数据结构里不剩东西了,就把 x x x 设置为 B B B )。
论文中使用链表实现。
分块优先队列维护两个链表 D 0 D_0 D 0 和 D 1 D_1 D 1 。D 0 D_0 D 0 用来放 Batch Prepend 操作添加的 key,D 1 D_1 D 1 用来放 Insert 操作添加的 key。D 0 D_0 D 0 和 D 1 D_1 D 1 是把块作为元素连接起来的链表。同一块内再维护一个链表,把块内的数据连接起来。(也就是链表套链表)。每个块都有一个属性:上界,块中所有的元素的 value 都不会超过这个上界。同时,对每个 key 维护其在链表中的指针(这里假设 key 不太大,可以按 node_map[key] = ptr 这样开桶维护,并且维护的开销也是 O ( 1 ) O(1) O ( 1 ) 的)。
初始化时,D 0 D_0 D 0 直接为空链表,D 1 D_1 D 1 里面放入一个空的块,块的上界为 B B B 。
插入时,首先看这个 key 之前是否已经被插入过了,如果是的话先直接删除(由于我们之前已经对 key 在链表中的指针进行了维护,所以这里无论是查找还是链表内删除,都是 O ( 1 ) O(1) O ( 1 ) 的)。然后找到上界最小且大于要插入的 value 的块进行插入。同一个块内的 key 的 value 是无序的,所以直接插,O ( 1 ) O(1) O ( 1 ) 。而找到要往哪个块插,需要我们先使用平衡树对每一个块的上界进行维护。后面我们会看到块的数量是 O ( max { 1 , N M } ) O(\max\left \{ 1, \frac{N}{M} \right \}) O ( max { 1 , M N } ) 的,于是查找要插的块的时间复杂度是 O ( max { 1 , log N M } ) O(\max\left\{ 1, \log \frac{N}{M} \right\}) O ( max { 1 , log M N } ) 。
值得注意的是,当把一个元素从块内删除之后,不需要更新块的上界。如果块内元素全被删干净时,这个块没了,此时不得将块的上界从平衡树中删掉,那也还是只需要 O ( max { 1 , log N M } ) O(\max\left\{ 1, \log \frac{N}{M} \right\}) O ( max { 1 , log M N } ) 的时间复杂度,不影响整个插入操作的时间复杂度。
将元素插入到一个块之后,有可能这个块的元素数量大于了 M M M ,于是需要对该块进行分割。分割的时候按 value 的中位数,分成两个 M 2 \frac{M}{2} 2 M 的块。之所以这么做是为了确保每个块内的元素的数量位于 M 2 \frac{M}{2} 2 M 和 M M M 之间,即块内元素的数量是 Θ ( M ) \Theta(M) Θ ( M ) 的,这样就能确保块的数量是 O ( max { 1 , N M } ) O(\max\left \{ 1, \frac{N}{M} \right \}) O ( max { 1 , M N } ) 的。
寻找中位数使用线性寻找中位数算法:BFPRT,所以这个分割操作的时间复杂度是 O ( M ) O(M) O ( M ) 的。
对于 Batch Prepend 操作,设批量插入 L L L 个元素。如果 L ≤ M L\le M L ≤ M ,那么直接在 D 0 D_{0} D 0 的最前面开一个块,把这 L L L 个元素丢进去,时间复杂度大概是 O ( L ) O(L) O ( L ) 的。反之,我们在 D 0 D_0 D 0 最前面创建 O ( L M ) O(\frac{L}{M}) O ( M L ) 个块,然后对于这 L L L 个元素不断地跑 BFPRT 算法按中位数切割(也就是先找中位数,把 L L L 劈成两部分,然后对于每个部分再去找中位数,这样总共会劈成 4 4 4 部分......这样一直分下去,直到劈成了 O ( L M ) O(\frac{L}{M}) O ( M L ) 部分)。这个不断切割的过程就像是一个递归树,层数为 O ( log L M ) O(\log \frac{L}{M}) O ( log M L ) ,同一层递归节点都要跑线性的 BFPRT 算法,所以总时间复杂度为 O ( L log L M ) O(L\log \frac{L}{M}) O ( L log M L ) 。
对于 Pull 操作,我们先从 D 0 D_0 D 0 中从前往后收集块,直到收集到的元素数量大于 M M M 或者没有更多块了。然后再从 D 1 D_1 D 1 中进行类似的操作。把从两个链表中收集到的元素全部拿到一个集合 S S S 中。如果 ∣ S ∣ ≤ M |S| \le M ∣ S ∣ ≤ M ,那么非常好处理,这里不多说。反之,由于每一个块内元素大小不超过 M M M ,所以最后 S S S 内的元素数量不会超过 M M M 太多,最多也就是 4 M 4M 4 M ,即 S = Θ ( M ) S=\Theta(M) S = Θ ( M ) 。然后我们可以直接在 S S S 上面跑 BFPRT 算法,在 Θ ( ∣ S ∣ ) = Θ ( M ) \Theta(|S|)=\Theta(M) Θ ( ∣ S ∣ ) = Θ ( M ) 的时间来把前 M M M 小的元素找到,然后把这 M M M 个元素逐个从删掉。逐个 O ( 1 ) O(1) O ( 1 ) 删除也只会有 O ( M ) O(M) O ( M ) 的时间。可能存在的问题是,删除时如果某个块变成空块了,就需要再付 O ( max { 1 , log N M } ) O(\max\left\{ 1, \log \frac{N}{M} \right\}) O ( max { 1 , log M N } ) 的块删除时间。这里我的理解是,你可以根据类似势能分析的思想,把这里删除块的时间复杂度算在插入块上。
总之,我们有:
操作
时间复杂度
Insert
O ( max { 1 , log N M } ) O(\max\left\{ 1, \log \frac{N}{M} \right\}) O ( max { 1 , log M N } )
Batch Prepend
O ( L log L M ) O(L\log \frac{L}{M}) O ( L log M L )
Pull
O ( M ) O(M) O ( M )
2.4 Finding Pivots
我们上面提到过说,Finding Pivots 会把 S S S 集合缩成 P P P ,后续只需要在 P P P 上做就是对的。我们来看 Finding Pivots 是怎么做的:
Finding Pivots 操作会进行类似于 Bellman ford 的 k k k 轮松弛。和 Bellman ford 不同的是,Bellman ford 每轮松弛会松弛整个图的边,而这里的松弛只涉及到起点位于 W i − 1 W_{i-1} W i − 1 内的边,W i − 1 W_{i-1} W i − 1 表示上一轮松弛中成功松弛到的点。直观上理解,每轮松弛都相当于一圈一圈往外延展。这里的 k k k 轮松弛还有种做 k k k 层 BFS 的味道。
最终会得到两个东西,一个东西是 W W W 集合,也就是这 k k k 轮松弛下来能够到达的所有点,或者说是从 S S S 内的点出发 k k k 步以内能到达的所有点。一个是 P P P 集合,表示后续的 BMSSP 算法只需要对 P P P 内连出去的边进行松弛即可,相当于把 S S S 缩小成了 P P P 。我们把 P P P 中的点称为 pivot。
关于 P P P 集合是怎么来的,以及为什么只需要关心 P P P 集合就是对的,接下来慢慢解释。
回忆起我们之前关于目标节点 U ~ \widetilde{U} U 的说法。上图中第 9 行的判断,只有当前距离小于 B B B 的 v v v 才加入到 W W W 并且在下一轮松弛。当前距离大于 B B B 的 v v v ,最短路肯定更是大于 B B B 了,那不归我们管,我们不放入 W W W 中。
U ~ \widetilde{U} U 必然可以被分成两类:
从 S S S 出发,走 k k k 步就能 complete,这种点已经被我们加入到了 W W W 中。
从 S S S 出发,走 k k k 步不能 complete,这种点可能需要再多走几步才能 complete。不过,至少我们可以知道,这样的点,从 S S S 出发的最短路径上,前 k k k 个点必然是已经在 W W W 中且已经 complete 的了(即对应上面那一类)。
于是 Finding Pivots 在最后,根据第 15 行的关系式,构建了一个”最短路树“(严谨点说应该是最短路森林)(最开始我们提出的第二个假设保证只会构成树而不会构成 DAG)。同时,这里的最短路树其实是假的,因为我们只做了 k k k 轮松弛,无法保证被松弛的到的点已经是 complete 的(甚至 S S S 中还有没 complete 的点)。也就是说,已经 complete 的点(也就是第一类的点)肯定在这个最短路树上,但是这个最短路树上还有没有 complete 的点,但这些点并不影响算法的正确性和复杂度。
然后,对于每个最短路树,如果其内部节点数量大于 k k k 的话,就将其加入 P P P 集合作为 pivot。
那么这里,对于走 k k k 步不能 complete 的点,必然有:其最短路经过 pivot。
于是同时也解释了为什么只需要看 P P P 集合就是对的。原来在 S S S 中但没在 P P P 中的点,其树的大小小于 k k k ,也就是从这个点出发,走不了 k k k 步就没节点可以继续走了。走 k k k 步不能 complete 的点必然不会经过这些点。
然后值得注意的是 W W W 。虽然走 k k k 步就能 complete 的点一定在 W W W 中,但是 W W W 中还有没 complete 的点,所以后面 BMSSP 算法还要再把这些点从 W W W 中筛出来。W W W 集合的产生既是顺手的事,也是必要的,W W W 集合用来处理走 k k k 步就能 complete 的点的最短路,P P P 集合和后面的 BMSSP 算法用来处理走 k k k 步不能 complete 的点的最短路。
然后回忆起 BMSSP 在最后说,如果 W W W 内的点 x x x 满足 d ^ ( x ) < B ′ \hat{d}(x) < B' d ^ ( x ) < B ′ ,那么他就是 complete 的。我们有如下证明:假如说 x x x 不是 complete 的,根据上面所说,x x x 的最短路必然经过某个 pivot。而 BMSSP 在若干次循环之后,一定做到了把最短路小于 B ′ B' B ′ 且经过 P P P 的点 complete ,于是这里发生矛盾。
Finding Pivots 操作需要保证给到后续 BMSSP 算法的集合 P P P 一定满足 P ≤ ∣ U ~ ∣ k P\le \frac{|\widetilde{U}|}{k} P ≤ k ∣ U ∣ 以保证后面 BMSSP 的时间复杂度。Finding Pivots 操作在最后有 ∣ W ∣ ≤ ∣ U ~ ∣ |W| \le |\widetilde{U} | ∣ W ∣ ≤ ∣ U ∣ ,同时由于我们 pivots 选择的规则,有 k ∣ P ∣ ≤ ∣ W ∣ k|P|\le |W| k ∣ P ∣ ≤ ∣ W ∣ ,两个合起来,就有 ∣ P ∣ ≤ ∣ U ~ ∣ k |P|\le \frac{|\widetilde{U} |}{k} ∣ P ∣ ≤ k ∣ U ∣ 。
同时注意到,如果在某轮松弛之后,出现了 ∣ W ∣ > k ∣ S ∣ |W|>k|S| ∣ W ∣ > k ∣ S ∣ ,而 ∣ W ∣ ≤ ∣ U ~ ∣ |W|\le |\widetilde{U} | ∣ W ∣ ≤ ∣ U ∣ ,于是我们有 k ∣ S ∣ < ∣ U ~ ∣ k|S| < |\widetilde{U}| k ∣ S ∣ < ∣ U ∣ ,∣ S ∣ < ∣ U ~ ∣ k |S|<\frac{|\widetilde{U} |}{k} ∣ S ∣ < k ∣ U ∣ ,那么此时就说明 S S S 其实本身规模就挺小的,直接把 S S S 作为 P P P 也可以接受,于是就直接把 S S S 返回了。在这种情况下返回的 W W W 我认为只是顺手的事了,貌似不太必要。
所以,∣ W ∣ = O ( min { k ∣ S ∣ , ∣ U ~ ∣ } ) |W|=O(\min\left\{ k|S|, |\widetilde{U}| \right\}) ∣ W ∣ = O ( min { k ∣ S ∣ , ∣ U ∣ } ) ,每轮松弛花费 ∣ W ∣ |W| ∣ W ∣ 的时间,如果没有提前返回,最后也是花 ∣ W ∣ |W| ∣ W ∣ 的时间生成集合 P P P 。最后总的时间复杂度主要在多轮松弛上,为 O ( k ∣ W ∣ ) = O ( min { k 2 ∣ S ∣ , k ∣ U ~ ∣ } ) O(k|W|)=O(\min\left\{ k^2|S|, k|\widetilde{U}| \right\}) O ( k ∣ W ∣ ) = O ( min { k 2 ∣ S ∣ , k ∣ U ∣ } ) 。
3. 正确性
首先是一堆数学符号:
T ( u ) T(u) T ( u ) 表示以 u u u 为根的最短路树。T ( u ) T(u) T ( u ) 内的点的最短路经过 u u u 。
T ( S ) = ⋃ v ∈ S T ( v ) T(S)=\bigcup_{v\in S}T(v) T ( S ) = ⋃ v ∈ S T ( v )
S ∗ = { v ∈ S : v is complete } S^*=\left \{ v\in S: v \text{ is complete} \right \} S ∗ = { v ∈ S : v is complete }
显然有 T ( S ∗ ) ⊆ T ( S ) T(S^*)\subseteq T(S) T ( S ∗ ) ⊆ T ( S ) 。
T < B ( S ) = { v ∈ T ( S ) : d ( v ) < B } T_{<B}(S) = \left \{ v \in T(S): d(v) < B \right \} T < B ( S ) = { v ∈ T ( S ) : d ( v ) < B } ,即最短路经过 S S S 且最短距离小于 B B B 的点。T < B ( S ) T_{<B}(S) T < B ( S ) 和上文提到的 U ~ \widetilde{U} U 是一致的。
T [ b , B ) ( S ) = { v ∈ T ( S ) : d ( v ) ∈ [ b , B } } T_{\left [ b, B \right )}(S) = \left \{ v \in T(S): d(v) \in \left [ b, B \right\} \right\} T [ b , B ) ( S ) = { v ∈ T ( S ) : d ( v ) ∈ [ b , B } }
Pull Minimum 引理
这个引理针对的是数据结构 D D D (分块优先队列)的 Pull 操作。回忆 Pull 操作,它会把 D D D 中的一个前缀 S i S_i S i 取出来,并返回一个上界 B i B_i B i 。也就是说,原本 D D D 中当前距离小于 B i B_i B i 的点都在 S i S_i S i 中,当前距离大于等于 B i B_i B i 的点都在 D ∖ S i D \setminus S_i D ∖ S i 中。
Pull Minimum 引理宣称,如果所有最短路距离小于 B B B 且还没 complete 的点都在 T ( D ∗ ) T(D^*) T ( D ∗ ) 里,那么
(a)所有最短路距离小于 B i B_i B i (这里的 B B B 是当前 BMSSP 的参数 B B B ,所以必然有 B i < B B_i < B B i < B )且还没 complete 的点都在 T ( S i ∗ ) T(S_i^*) T ( S i ∗ ) 中。
这意味着,在前提条件成立的情况下,如果我们把 D D D 以 B i B_i B i 分成两部分,那么所有最短路距离小于 B i B_i B i 且还没 complete 的点就归 S i S_i S i 管了,D D D 剩下的部分就不管了。也就解释了为什么 BMSSP 中 S i S_i S i 和 B i B_i B i 是满足子 BMSSP 第二个条件的。
(b)如果 B i ′ < B i B_i' < B_i B i ′ < B i ,那么 T < B i ′ ( D ) = T < B i ′ ( S i ) T_{<B_i'}(D)=T_{<B_i'}(S_i) T < B i ′ ( D ) = T < B i ′ ( S i )
这意味着,即使后续下层 BMSSP 是部分执行的,仅处理完了 S i S_i S i 中最短路距离小于 B i ′ B_i' B i ′ 的部分,但实际上是处理完了整个 D D D 中最短路距离小于 B i ′ B_i' B i ′ 的部分。
证明:
所有最短路距离小于 B i B_i B i 且还没 complete 的点 v v v ,根据引理的前提,其最短路必然是经过 D D D 中某个已经 complete 的点 u u u 。那么就有 d ^ ( u ) = d ( u ) ≤ d ( v ) < B i \hat{d}(u)=d(u)\le d(v)<B_i d ^ ( u ) = d ( u ) ≤ d ( v ) < B i ,于是 u u u 肯定是位于 S i S_i S i 中的,于是 v v v 肯定是位于 T ( S i ∗ ) T(S_i^*) T ( S i ∗ ) 中的。(a)得证。
对于 (b),肯定有 T < B i ′ ( S i ) ⊆ T < B i ′ ( D ) T_{<B_i'}(S_i)\subseteq T_{<B_i'}(D) T < B i ′ ( S i ) ⊆ T < B i ′ ( D ) 。然后对于 v ∈ T < B i ′ ( D ) v\in T_{<B_i'}(D) v ∈ T < B i ′ ( D ) 有 d ( v ) < B i ′ < B i < B d(v) < B_i' < B_i < B d ( v ) < B i ′ < B i < B ,那么根据 (a)有 v ∈ T ( S i ∗ ) v \in T(S_i^*) v ∈ T ( S i ∗ ) 。设 v v v 最短路经过 S i ∗ S_i^* S i ∗ 中的点 x x x ,那么就有 d ^ ( x ) = d ( x ) ≤ d ( v ) < B i ′ \hat{d}(x) = d(x) \le d(v) < B_i' d ^ ( x ) = d ( x ) ≤ d ( v ) < B i ′ ,于是 v v v 位于 T < B i ′ ( S i ) T_{<B_i'}(S_i) T < B i ′ ( S i ) 中,即 T < B i ′ ( D ) ⊆ T < B i ′ ( S i ) T_{<B_i'}(D) \subseteq T_{<B_i'}(S_i) T < B i ′ ( D ) ⊆ T < B i ′ ( S i ) ,于是(b)得证。
BMSSP 正确性
然后我们证明 BMSSP 的正确性,也就是证明:给定 l l l 、B B B 和集合 S S S ,且 S S S 满足两点大前提:
∣ S ∣ ≤ 2 l t |S| \le 2^{lt} ∣ S ∣ ≤ 2 lt
所有 d ( x ) < B d(x) < B d ( x ) < B 且还没 complete 的 x x x ,其最短路必须经过 S S S 中某个已经 complete 的点,也就是 T < B ( S ∗ ) = T < B ( S ) T_{<B}(S^*) = T_{<B}(S) T < B ( S ∗ ) = T < B ( S )
那么 BMSSP 一定会给出一个 U U U 和 B ′ B' B ′ ,满足经过 S S S 且最短路距离小于 B ′ B' B ′ 的点均在 U U U 中且已经 complete。
证明:
使用双重数学归纳法。
对于 l = 0 l=0 l = 0 ,S S S 中必然只有一个点 x x x (根据前提 1)。然后 BMSSP 的 BaseCase 会从 x x x 开始做传统 Dijkstra,根据 Dijkstra 的正确性和前提 2,上述结论显然是成立的。
现在假设 l − 1 l-1 l − 1 成立,我们要证明 l l l 也成立。
在第 l l l 层,我们会进行若干次循环迭代(BMSSP 算法第 8 8 8 行)。我们令 D i D_i D i 为第 i i i 次循环前 D D D 的样子。然后我们宣称下面两个命题在第 i i i 次循环前必然是成立的:
所有最短路距离小于 B B B 且还没 complete 的点,都位于 T [ B i − 1 ′ , B ) ( P ) T_{\left [ B_{i-1}',B \right)}(P) T [ B i − 1 ′ , B ) ( P ) 中。P P P 是第 l l l 层最开始 Finding Pivots 操作得到的集合 P P P 。
这意味着...其实挺显然的,你可以想象 BMSSP 的这个循环主要是在松弛,然后松弛就意味着延申最短路树。而这个延申最开始就是从 P P P 出发的。
T [ B i − 1 ′ , B ) ( P ) = T < B ( D i ) = T < B ( D i ∗ ) T_{\left [ B_{i-1}', B \right )}(P) = T_{<B}(D_i) = T_{<B}(D_i^*) T [ B i − 1 ′ , B ) ( P ) = T < B ( D i ) = T < B ( D i ∗ )
这意味着,这些还没 complete 的点可以直接在 D i D_i D i 上甚至是 D i ∗ D_i^* D i ∗ 上处理。
为了证明上述两个命题成立,我们再用数学归纳法。
首先对于 i = 1 i=1 i = 1 ,我们已经在 Finding Pivots 操作中说过,所有最短路距离小于 B B B 且还没 complete 的点(已经 complete 的位于 W W W 集合了),最短路一定经过 P P P 。由于 T < B 0 ′ ( P ) T_{<B_0'}(P) T < B 0 ′ ( P ) 是空集,所以 T < B ( P ) = T [ B 0 ′ , B ) ( P ) T_{<B}(P)=T_{\left[ B_0', B \right)}(P) T < B ( P ) = T [ B 0 ′ , B ) ( P ) ,命题一成立。
同时,必然有 T < B ( P ) = T < B ( P ∗ ) T_{<B}(P)=T_{<B}(P^*) T < B ( P ) = T < B ( P ∗ ) (也就是 P P P 中所有没被 complete 的点,最短路经过 P ∗ P^* P ∗ 。假如说不经过 P ∗ P^* P ∗ ,那么根据 BMSSP 的第二个要求,必然是经过 ( S ∖ P ) ∗ (S\setminus P)^* ( S ∖ P ) ∗ 。而根据 Finding Pivots 操作,最短路经过 S ∖ P S\setminus P S ∖ P 的点必然已经 complete 掉并放入 W W W 中了,否则被经过的点就不在 S ∖ P S\setminus P S ∖ P 里了)。再加上初始时 P = D i P = D_i P = D i ,命题二也成立。
假设对于 i i i 成立。由于上述命题是在 i i i 循环之前成立的,所以我们把 i i i 循环的过程走一遍就会得到 i + 1 i+1 i + 1 循环前的情况,也就能判断对于 i + 1 i+1 i + 1 是否成立。
首先第 i i i 个循环先从 D i D_i D i 中取出了 S i S_i S i 。然后去调用了 l − 1 l-1 l − 1 层的 BMSSP。由于我们前面假设了 l − 1 l-1 l − 1 层 BMSSP 成立,那么我们完全相信这个 BMSSP 的结果是对的。不过我们还需要确保满足 l − 1 l-1 l − 1 层 BMSSP 的两个前提要求:
由于 D D D 的 M = 2 ( l − 1 ) t M=2^{(l-1)t} M = 2 ( l − 1 ) t ,所以 ∣ S i ∣ < 2 ( l − 1 ) t |S_i| < 2^{(l-1)t} ∣ S i ∣ < 2 ( l − 1 ) t ,因此满足第一个前提
由于 i i i 成立,根据命题二,所有最短路距离小于 B B B 且还没 complete 的点一定位于 T ( D i ∗ ) T(D_i^*) T ( D i ∗ ) 中,那么对于 B i < B B_i<B B i < B 肯定也是成立的,于是满足 Pull Minimum 引理的前提。再根据 Pull Minimum 引理的(a),最短路距离小于 B i B_i B i 且还没 complete 的点一定位于 T ( S i ∗ ) T(S_i^*) T ( S i ∗ ) 中,因此满足第二个前提。
l − 1 l-1 l − 1 层的 BMSSP 返回了一个 U i U_i U i 和 B i ′ B_i' B i ′ 。原来最短路距离位于 [ B i − 1 ′ , B i ′ ) [B_{i-1}',B_i') [ B i − 1 ′ , B i ′ ) 的还没 complete 的点,即 T [ B i − 1 ′ , B i ′ ) ( P ) T_{[B_{i-1}',B_i')} (P) T [ B i − 1 ′ , B i ′ ) ( P ) ,如今已经全部放入 U i U_i U i 且 complete 了。由于 i i i 成立,根据命题一,原来剩下的还要处理的点在 T [ B i − 1 ′ , B ) ( P ) T_{\left [ B_{i-1}',B \right)}(P) T [ B i − 1 ′ , B ) ( P ) 中,现在 T [ B i − 1 ′ , B i ′ ) ( P ) T_{[B_{i-1}',B_i')} (P) T [ B i − 1 ′ , B i ′ ) ( P ) 解决了,于是还剩下 T [ B i ′ , B ) ( P ) T_{[B_i', B)} (P) T [ B i ′ , B ) ( P ) 没搞,于是命题一对于 i + 1 i+1 i + 1 成立。
由于 i i i 成立,根据命题二,T [ B i − 1 ′ , B ) ( P ) = T < B ( D i ) = T < B ( D i ∗ ) T_{\left [ B_{i-1}', B \right )}(P) = T_{<B}(D_i) = T_{<B}(D_i^*) T [ B i − 1 ′ , B ) ( P ) = T < B ( D i ) = T < B ( D i ∗ ) 。首先,根据这个算法步骤我们可以看到,D i D_i D i 中最小的元素必然也是大于 B i − 1 ′ B_{i-1}' B i − 1 ′ 的,那么经过 D i D_i D i 的距离小于 B i − 1 ′ B_{i-1}' B i − 1 ′ 的最短路是不存在的,于是实际上有:T [ B i − 1 ′ , B ) ( P ) = T [ B i − 1 ′ , B ) ( D i ) = T [ B i − 1 ′ , B ) ( D i ∗ ) T_{\left [ B_{i-1}', B \right )}(P) = T_{\left [ B_{i-1}', B \right )}(D_i) = T_{\left [ B_{i-1}', B \right )}(D_i^*) T [ B i − 1 ′ , B ) ( P ) = T [ B i − 1 ′ , B ) ( D i ) = T [ B i − 1 ′ , B ) ( D i ∗ ) 。三个集合相等,对应的区间也都是一样的,于是我们对这三个集合同时取距离位于 [ B i , B ) [B_i,B) [ B i , B ) 的这部分,就有了 T [ B i ′ , B ) ( P ) = T [ B i ′ , B ) ( D i ) = T [ B i ′ , B ) ( D i ∗ ) T_{\left [ B_{i}', B \right )}(P) = T_{\left [ B_{i}', B \right )}(D_i) = T_{\left [ B_{i}', B \right )}(D_i^*) T [ B i ′ , B ) ( P ) = T [ B i ′ , B ) ( D i ) = T [ B i ′ , B ) ( D i ∗ ) 。
要想证明命题二,我们可以试试能不能证明下两个子集关系:
T [ B i ′ , B ) ( D i ∗ ) ⊆ T < B ( D i + 1 ∗ ) T_{\left [ B_{i}', B \right )}(D_i^*) \subseteq T_{<B}(D_{i+1}^*) T [ B i ′ , B ) ( D i ∗ ) ⊆ T < B ( D i + 1 ∗ )
T < B ( D i + 1 ) ⊆ T [ B i ′ , B ) ( D i ) T_{<B}(D_{i+1})\subseteq T_{\left [ B_{i}', B \right )}(D_i) T < B ( D i + 1 ) ⊆ T [ B i ′ , B ) ( D i )
如果可以的话,我们就能有:
T [ B i ′ , B ) ( P ) ⊆ T [ B i ′ , B ) ( D i ∗ ) ⊆ T < B ( D i + 1 ∗ ) ⊆ T < B ( D i + 1 ) ⊆ T [ B i ′ , B ) ( D i ) ⊆ T [ B i ′ , B ) ( P ) T_{\left [ B_{i}', B \right )}(P) \subseteq T_{\left [ B_{i}', B \right )}(D_i^*) \subseteq T_{<B}(D_{i+1}^*) \subseteq T_{<B}(D_{i+1}) \subseteq T_{\left [ B_{i}', B \right )}(D_i) \subseteq T_{\left [ B_{i}', B \right )}(P)
T [ B i ′ , B ) ( P ) ⊆ T [ B i ′ , B ) ( D i ∗ ) ⊆ T < B ( D i + 1 ∗ ) ⊆ T < B ( D i + 1 ) ⊆ T [ B i ′ , B ) ( D i ) ⊆ T [ B i ′ , B ) ( P )
由于首位相等,于是中间的全部集合都是相等的,于是命题二对于 i + 1 i+1 i + 1 也成立。
D D D 中的元素大多数都是上轮循环继承到下轮循环。少数发生变动的元素是由于 Pull 操作和后面松弛完了之后再插入进去导致的。
对于子集关系一,取一个节点 y ∈ D i ∗ ∖ D i + 1 ∗ y\in D_i^* \setminus D_{i+1}^* y ∈ D i ∗ ∖ D i + 1 ∗ 。这个节点在上一轮不在下一轮,说明他被 Pull 出去了。同时他还没由于 BatchPrepend 加回到这一轮,说明 d ^ ( y ) < B i ′ \hat{d}(y) < B_i' d ^ ( y ) < B i ′ ,因此 y ∈ U i y \in U_i y ∈ U i 。现在任取 x ∈ T [ B i ′ , B ) ( y ) x\in T_{\left [ B_{i}', B \right )}(y) x ∈ T [ B i ′ , B ) ( y ) ,由于 d ( x ) ≥ B i ′ d(x)\ge B_i' d ( x ) ≥ B i ′ ,所以 x x x 肯定不在 U i U_i U i 中。考虑从 y y y 到 x x x 的最短路径,y ∈ U i y\in U_i y ∈ U i 而 x ∉ U i x\notin U_i x ∈ / U i ,说明路径上必然存在一条边 ( u , v ) (u,v) ( u , v ) 使得 u ∈ U i u\in U_i u ∈ U i 而 v ∉ U i v\notin U_i v ∈ / U i 。我们后面对 U i U_i U i 内所有点的出边都进行了松弛,且 ( u , v ) (u,v) ( u , v ) 在最短路径上,因此 ( u , v ) (u,v) ( u , v ) 必然得到了有效松弛。而 v v v 必然是位于了 D i + 1 ∗ D_{i+1}^* D i + 1 ∗ 中的,而 x x x 是 v v v 的下游,因此 x ∈ T ( v ) ⊆ T ( D i + 1 ∗ ) x\in T(v) \subseteq T(D_{i+1}^*) x ∈ T ( v ) ⊆ T ( D i + 1 ∗ ) ,得证。
对于子集关系二,还是取一个节点 y ∈ D i + 1 ∖ D i y\in D_{i+1}\setminus D_i y ∈ D i + 1 ∖ D i 。这个点不在上一轮而在下一轮,这是怎么来的?是因为我们在处理 U i U_i U i 的时候松弛了边 ( u , y ) (u,y) ( u , y ) ,并且这个松弛是有效松弛,所以就被加进 D D D 了。我们取最后一次让 y y y 有效更新的边 ( u , y ) (u, y) ( u , y ) ,然后会分成如下几种情况:
如果 y y y 是 complete 的,由于 ( u , y ) (u, y) ( u , y ) 是 y y y 的最后一次有效松弛,那么一定有 y ∈ T ( u ) y\in T(u) y ∈ T ( u ) ,即 y y y 的真实最短路是从 u u u 来,且 d ( y ) ≥ B i ′ d(y) \ge B_i' d ( y ) ≥ B i ′ 。又因为 u ∈ U i = T < B i ′ ( D i ) u\in U_i = T_{<B_i'}(D_i) u ∈ U i = T < B i ′ ( D i ) ,所以 y ∈ T [ B i ′ , B ) ( D i ) y \in T_{\left [ B_{i}', B \right )}(D_i) y ∈ T [ B i ′ , B ) ( D i )
如果 y y y 不是 complete 的,根据已经证明的 i + 1 i+1 i + 1 轮成立的命题一,这样的 y ∈ T [ B i ′ , B ) ( P ) y \in T_{[B_i', B)}(P) y ∈ T [ B i ′ , B ) ( P ) ,然后由对于 i i i 轮成立的命题二,y ∈ T [ B i ′ , B ) ( D i ) y\in T_{[B_i', B)}(D_i) y ∈ T [ B i ′ , B ) ( D i )
无论是哪种情况,都有 y ∈ T [ B i ′ , B ) ( D i ) y \in T_{\left [ B_{i}', B \right )}(D_i) y ∈ T [ B i ′ , B ) ( D i ) 。D i + 1 D_{i+1} D i + 1 出发的某个路径可能是从 y y y 出发的,但是 y y y 又是从 D i D_i D i 出发的,所以 T < B ( D i + 1 ) ⊆ T [ B i ′ , B ) ( D i ) T_{<B}(D_{i+1})\subseteq T_{\left [ B_{i}', B \right )}(D_i) T < B ( D i + 1 ) ⊆ T [ B i ′ , B ) ( D i ) 得证。于是关于循环的两个命题得证。
然后我们可以看到,假如说 l l l 层 BMSSP 进行了 q q q 次循环,那么这 q q q 次循环实际会收集 q q q 个 complete 的点的集合:U 1 , U 2 , … , U q U_1, U_2, \dots, U_q U 1 , U 2 , … , U q ,其中 U i = T [ B i − 1 ′ , B i ′ ) ( P ) U_i = T_{\left [ B_{i-1}',B_i' \right)}(P) U i = T [ B i − 1 ′ , B i ′ ) ( P ) (根据循环的命题一)。我们看到我们其实是把 P P P 出去的距离小于 B ′ B' B ′ 的最短路,不重不漏地分成了 q q q 个部分,并且全部 complete 掉了。
现在还差一步。对于 S ∖ P S\setminus P S ∖ P 出去的最短路,这些点已经在 Finding Pivots 操作中就 complete 了,BMSSP 在最后也对其进行了收集。因此 BMSSP 是正确的。
4. 时间复杂度
BMSSP 的过程构成了一个递归树。我们设递归树上的某个节点 x x x 的参数是 l x l_x l x 、B x B_x B x 和 S x S_x S x ,目标节点为 U ~ x \widetilde{U}_x U x ,Finding Pivots 操作获得的集合为 P x P_x P x ,BMSSP 返回 B x ′ B_x' B x ′ 和 U x U_x U x 。
首先,在递归树的同一深度上,所有节点 x x x 返回的集合 U x U_x U x 是互不相交的,这一点已经在 BMSSP 正确性的部分中提到过了。因此同一层的 ∣ U x ∣ |U_x| ∣ U x ∣ 加起来应该是小于等于 n n n 的。
由于我们之前说顶层 BMSSP 取 l = ⌈ log n t ⌉ l=\left\lceil \frac{\log n}{t} \right\rceil l = ⌈ t l o g n ⌉ ,所以整个递归树的层数是 log n t = O ( log 1 3 n ) \frac{\log n}{t}=O(\log ^{\frac{1}{3}} n) t l o g n = O ( log 3 1 n ) 的。
∣ P ∣ |P| ∣ P ∣ 规模
对于 Finding Pivots 操作返回的 P P P ,有 ∣ P ∣ ≤ ∣ U ~ ∣ k |P| \le \frac{|\widetilde{U}|}{k} ∣ P ∣ ≤ k ∣ U ∣ 。进一步,我们还有:
如果 x x x 是成功执行的,就有 U = U ~ U = \widetilde{U} U = U ,那么 ∣ P ∣ ≤ ∣ U ∣ k |P| \le \frac{|U|}{k} ∣ P ∣ ≤ k ∣ U ∣
如果 x x x 是部分执行的,∣ U ∣ ≥ k 2 l t ≥ k ∣ S ∣ |U| \ge k2^{lt} \ge k|S| ∣ U ∣ ≥ k 2 lt ≥ k ∣ S ∣ ,因此 ∣ P ∣ ≤ ∣ S ∣ ≤ ∣ U ∣ k |P| \le |S| \le \frac{|U|}{k} ∣ P ∣ ≤ ∣ S ∣ ≤ k ∣ U ∣
总之就是 ∣ P ∣ ≤ ∣ U ∣ k |P| \le \frac{|U|}{k} ∣ P ∣ ≤ k ∣ U ∣ 。
∣ U i ∣ |U_i| ∣ U i ∣ 规模
完成完 BMSSP 之后,∣ U ∣ ≤ 4 k 2 l t |U| \le 4k2^{lt} ∣ U ∣ ≤ 4 k 2 lt 。如果 BMSSP 部分执行的话,我们还有 ∣ U ∣ ≥ k 2 l t |U| \ge k2^{lt} ∣ U ∣ ≥ k 2 lt 。
对于后者显然成立,下面证明前者:
对 BMSSP 的层数 l l l 使用数学归纳法。l = 0 l=0 l = 0 时,朴素的 Dijkstra 每次循环都只往 ∣ U ∣ |U| ∣ U ∣ 中加一个点,加到 ∣ U ∣ ≥ k 2 l t |U| \ge k2^{lt} ∣ U ∣ ≥ k 2 lt 就立刻停止。所以 ∣ U ∣ ≤ 4 k 2 l t |U| \le 4k2^{lt} ∣ U ∣ ≤ 4 k 2 lt 是显然成立的。
假设对于 l − 1 l-1 l − 1 成立。
设 BMSSP 一共进行了 q q q 次循环,那么在第 q q q 次循环之前,必然有 ∣ U before_q ∣ < k 2 l t |U_{\text{before\_q}}|<k2^{lt} ∣ U before_q ∣ < k 2 lt 。
第 q q q 次循环,会调用下层 BMSSP,得到返回结果 U q U_q U q ,并把 U q U_q U q 加进 U U U 里面。由于对于 l − 1 l-1 l − 1 成立,所以 ∣ U q ∣ ≤ 4 k 2 ( l − 1 ) t |U_q|\le 4k2^{(l-1)t} ∣ U q ∣ ≤ 4 k 2 ( l − 1 ) t 。
循环结束之后,还会对 W W W 集合进行过滤得到 W filtered W_{\text{filtered}} W filtered ,并把 W filtered W_{\text{filtered}} W filtered 加入到 U U U 里面。根据 Finding Pivots 的性质,∣ W ∣ ≤ k ∣ S ∣ |W|\le k|S| ∣ W ∣ ≤ k ∣ S ∣ ,又因为 BMSSP 的前提条件 ∣ S ∣ ≤ 2 l t |S| \le 2^{lt} ∣ S ∣ ≤ 2 lt ,因此有 ∣ W filtered ∣ ≤ ∣ W ∣ ≤ k 2 l t |W_{\text{filtered}}| \le |W| \le k2^{lt} ∣ W filtered ∣ ≤ ∣ W ∣ ≤ k 2 lt (这里就体现出当时 Finding Pivots 控制 W W W 规模的重要性了)
最终返回大小
∣ U ∣ ≤ ∣ U before_q ∣ + ∣ U q ∣ + ∣ W filtered ∣ ≤ k 2 l t + 4 k 2 ( l − 1 ) t + k 2 l t = 2 k 2 l t + 4 k ⋅ 2 l t 2 t \begin{align}
|U|
&\le |U_{\text{before\_q}}| + |U_q| + |W_{\text{filtered}}| \\
&\le k2^{lt} + 4k2^{(l-1)t} + k2^{lt} \\
&= 2k2^{lt} + 4k\cdot \frac{2^{lt}}{2^t}
\end{align}
∣ U ∣ ≤ ∣ U before_q ∣ + ∣ U q ∣ + ∣ W filtered ∣ ≤ k 2 lt + 4 k 2 ( l − 1 ) t + k 2 lt = 2 k 2 lt + 4 k ⋅ 2 t 2 lt
由于 t = ⌊ log 2 3 n ⌋ t = \left \lfloor \log^{\frac{2}{3}}n \right \rfloor t = ⌊ log 3 2 n ⌋ ,对于足够大的 n n n ,显然有 2 t ≥ 4 2^t\ge 4 2 t ≥ 4 ,因此 ∣ U ∣ ≤ 2 k 2 l t + k 2 l t = 3 k 2 l t ≤ 4 k 2 l t |U| \le 2k2^{lt} + k2^{lt} = 3k2^{lt} \le 4k2^{lt} ∣ U ∣ ≤ 2 k 2 lt + k 2 lt = 3 k 2 lt ≤ 4 k 2 lt ,证毕。
分块优先队列的时间复杂度
在第 l l l 层的 BMSSP 中,有 M = 2 ( l − 1 ) t M=2^{(l-1)t} M = 2 ( l − 1 ) t 。然后观察算法过程,D D D 中元素一部分来自 P P P ,另一部分来自松弛操作加入的点。根据出度为常数的假设,松弛操作最多会进行 O ( ∣ U ∣ ) O(|U|) O ( ∣ U ∣ ) 次,而 ∣ U ∣ |U| ∣ U ∣ 受到 k 2 l t k2^{lt} k 2 lt 的限制。同时由于 ∣ P ∣ ≤ U k |P| \le \frac{U}{k} ∣ P ∣ ≤ k U ,因此 D D D 最多会被塞入 N = O ( k 2 l t ) N = O(k2^{lt}) N = O ( k 2 lt ) 个点。
Insert:时间复杂度为 O ( max { 1 , log N M } ) O(\max\left\{ 1, \log \frac{N}{M} \right\}) O ( max { 1 , log M N } ) 。带入 N N N 和 M M M ,有 O ( log k 2 l t 2 ( l − 1 ) t ) = O ( log ( k 2 t ) ) = O ( t + log k ) O(\log \frac{k2^{lt}}{2^{(l-1)t}})=O(\log (k2^t)) = O(t+\log k) O ( log 2 ( l − 1 ) t k 2 lt ) = O ( log ( k 2 t )) = O ( t + log k ) ,由于 log k = log ( log 1 3 n ) ≪ t \log k = \log (\log ^{\frac{1}{3}}n) \ll t log k = log ( log 3 1 n ) ≪ t ,所以可以认为单次 Insert 操作时间复杂度为 O ( t ) O(t) O ( t ) 。
BatchPrepend:时间复杂度为 O ( L log L M ) O(L\log \frac{L}{M}) O ( L log M L ) 。而批量前插的元素集合 K K K 的大小是 O ( ∣ U i ∣ ) = O ( k 2 l t ) = O ( k M ) O(|U_i|) = O(k2^{lt}) = O(kM) O ( ∣ U i ∣ ) = O ( k 2 lt ) = O ( k M ) ,因此均摊到每个元素上的时间复杂度为 O ( log L M ) = O ( log k M M ) = O ( log k ) O(\log \frac{L}{M}) = O(\log \frac{kM}{M}) = O(\log k) O ( log M L ) = O ( log M k M ) = O ( log k ) 。
Pull:拉取 M M M 个的时间复杂度为 O ( M ) O(M) O ( M ) ,那么均摊到每个元素上,时间复杂度为 O ( 1 ) O(1) O ( 1 ) 。
Lemma 3.11
min { k ∣ S ∣ , ∣ U ~ ∣ } ≤ ∣ U ∣ \min\left\{ k|S|, |\widetilde{U}| \right\} \le |U| min { k ∣ S ∣ , ∣ U ∣ } ≤ ∣ U ∣
其中 U ~ \widetilde{U} U 是我们的 BMSSP 理论上想要处理的所有目标(具体定义在 Finding Pivots 那里),U U U 表示 BMSSP 实际处理的节点数量。如果 BMSSP 成功执行了,那么 U = U ~ U=\widetilde{U} U = U 。如果部分执行,则 U ⊂ U ~ U\subset \widetilde{U} U ⊂ U 。
证明:
如果成功执行了,那么 U = U ~ U=\widetilde{U} U = U ,上式显然成立。
如果没成功执行,必然有 ∣ U ∣ ≥ k 2 l t |U| \ge k2^{lt} ∣ U ∣ ≥ k 2 lt 。根据 BMSSP 的第一个前提,∣ S ∣ ≤ 2 l t |S| \le 2^{lt} ∣ S ∣ ≤ 2 lt ,因此 k ∣ S ∣ ≤ k 2 l t ≤ ∣ U ∣ k|S| \le k2^{lt} \le |U| k ∣ S ∣ ≤ k 2 lt ≤ ∣ U ∣ ,证毕。
然后整个算法的过程的时间开销主要分成如下几个部分:
(1)全图所有 FindPivots 开销
递归树某个节点 x x x 的时间复杂度为 O ( min { k 2 ∣ S x ∣ , k ∣ U ~ x ∣ } ) = O ( k ∣ U x ∣ ) O(\min\left\{ k^2|S_x|, k|\widetilde{U}_x| \right\}) = O(k|U_x|) O ( min { k 2 ∣ S x ∣ , k ∣ U x ∣ } ) = O ( k ∣ U x ∣ )
同一层的所有节点开销加起来,为 O ( k n ) O(kn) O ( kn )
递归树一共 O ( log 1 3 n ) O(\log ^{\frac{1}{3}} n) O ( log 3 1 n ) 层,于是该部分总耗时为 O ( n log 2 3 n ) O(n\log^{\frac{2}{3}}n) O ( n log 3 2 n ) 。
(2)Pivots 插入数据结构 D D D 中开销
由于 ∣ P ∣ ≤ ∣ U ∣ k |P| \le \frac{|U|}{k} ∣ P ∣ ≤ k ∣ U ∣ 。因此把 P P P 插入到数据结构 D D D 的开销是 O ( ∣ U x ∣ k ⋅ t ) O(\frac{|U_x|}{k} \cdot t) O ( k ∣ U x ∣ ⋅ t ) 。
同一层全部节点加起来的总开销为 O ( n t k ) O(\frac{nt}{k}) O ( k n t ) 。
整个递归树在该部分的开销为 O ( n t k ) × O ( log 1 3 n ) = O ( n log 2 3 n ) O(\frac{nt}{k}) \times O(\log ^{\frac{1}{3}} n) = O(n\log^{\frac{2}{3}}n) O ( k n t ) × O ( log 3 1 n ) = O ( n log 3 2 n ) 。
(3)边松弛导致直接插入开销
如果在某次松弛中,一条边 ( u , v ) (u,v) ( u , v ) 让 v v v 直接触发了 Insert 操作,说明 v v v 经过松弛后 B i ≤ d ^ ( v ) = d ( u ) + w u v ≤ B B_i \le \hat{d}(v) = d(u) + w_{uv} \le B B i ≤ d ^ ( v ) = d ( u ) + w uv ≤ B 。事实上,每条边,在整个算法的过程中,只会触发一次直接 Insert 操作。
证明:
首先由于 BMSSP 得到的 U i U_i U i 都是不交的,所以在当前 BMSSP 的过程中,u u u 只会被扫到一次,因此一条边也只会被松弛一次,也最多只会触发一次 Insert 操作。
假如 u u u 所在的 U i U_i U i 是由某个子 BMSSP 得到的,该子 BMSSP 的参数 B B B 是 B i B_i B i 。在该子 BMSSP 中,假如说 u u u 又作为某个 U i U_i U i 内的点往外松弛的话,遇到 ( u , v ) (u,v) ( u , v ) 这条边之后,由于 B i ≤ d ( u ) + w u v B_i \le d(u) + w_{uv} B i ≤ d ( u ) + w uv ,所以根本不会触发松弛,也就不会触发 Insert 操作。
u u u 在调用当前 BMSSP 的父节点里往外松弛的时候,如果遇到了 ( u , v ) (u,v) ( u , v ) 这条边,由于 d ( u ) + w u v ≤ B d(u) + w_{uv} \le B d ( u ) + w uv ≤ B ,这里的 B B B 是当前 BMSSP 的上限,却是父 BMSSP 想要触发 Insert 的下限,因此不会触发 Insert 操作。
边的数量是 Θ ( n ) \Theta(n) Θ ( n ) 的,所以总的时间复杂度为 O ( n t ) = O ( n log 2 3 n ) O(nt) = O(n\log^{\frac{2}{3}}n) O ( n t ) = O ( n log 3 2 n ) 。
(4)BatchPrepend 开销
每次 BatchPrepend,单个元素的均摊耗时为 O ( log k ) O(\log k) O ( log k ) 。最坏情况下,整个一层得到的所有 U i U_i U i 内的点松弛出去的边都会导致终点被加入到集合 K K K 从而被 BatchPrepend 到 D D D 中。由于 U i U_i U i 无交,以及出度为常数的假设,所以被 BatchPrepend 的点的数量是最多 O ( n ) O(n) O ( n ) 的。一共是 O ( log 1 3 n ) O(\log ^{\frac{1}{3}} n) O ( log 3 1 n ) 层,因此总耗费时间为 O ( n ⋅ log ( log 1 3 n ) log 1 3 n ) O(n\cdot \log(\log^{\frac{1}{3}} n) \log ^{\frac{1}{3}} n) O ( n ⋅ log ( log 3 1 n ) log 3 1 n ) ,比 O ( n log 2 3 n ) O(n\log ^{\frac{2}{3}}n) O ( n log 3 2 n ) 小,可以忽略不计。
于是我们就得到了 O ( n log 2 3 n ) O(n\log^{\frac{2}{3}} n) O ( n log 3 2 n ) 的时间复杂度。
5. 不失一般性
上面所有的分析都是在论文的两个假设上进行的。下面我们说明论文的这两个假设不失一般性:
常度数图
对于任意的图 G G G ,都可以进行如下的转换,得到图 G ′ G' G ′ 。图 G ′ G' G ′ 满足入度和出度最多只有 2 2 2 ,且 G G G 中的最短路在 G ′ G' G ′ 中仍得到保留:
对于所有的节点 v v v ,假如说他有个相邻节点 w w w (无论这个节点是从 v v v 出去的还是指向 v v v 的),我们都在 G ′ G' G ′ 中对应一个节点 x v w x_{vw} x v w 。所有关于 v v v 的这种节点都用边权为 0 0 0 的变强连通起来(比如构成一个环)。
然后对于 G G G 中的所有边 ( u , v ) (u,v) ( u , v ) ,都加一个从 x u v x_{uv} x uv 指向 x v u x_{vu} x vu 的有向边,边权为 w u v w_{uv} w uv 。
比如:
转为:
我们可以看到,红圈内的点相当于原图中的一个点(相当于我们把原图中的点拆开了),而同一个圈内的节点是强连通,且互相之间的距离为 0 0 0 ,因此最短路得到了保留。而且从这个构造过程中也能发现,点入度和出度最多只有 2 2 2 ,且 n ′ = Θ ( m ′ ) = Θ ( m ) n'=\Theta(m')=\Theta(m) n ′ = Θ ( m ′ ) = Θ ( m ) (n ′ n' n ′ 和 m ′ m' m ′ 为 G ′ G' G ′ 的点和边的数量,m m m 为 G G G 的边的数量)
每个点的最短路均不同
如果两个最短路的长度是相同的,我们还可以再去比较跳数(即经过的点数),跳数相同则继续比最短路终点编号,终点编号相同再去比较倒数第二个点的节点编号,还相同则继续比较倒数第三个节点的编号....这样一直比下去肯定能决出胜负,否则的话这两个最短路就是完全一致的。
这样的比较看起来是 O ( n ) O(n) O ( n ) 的。但其实论文中的最短路比较一共分成两种情况:
当对边 ( u , v ) (u,v) ( u , v ) 进行松弛时,此时 v v v 原来的最短路可能是这样 s → ⋯ → x → v s\rightarrow \dots \rightarrow x \rightarrow v s → ⋯ → x → v ,现在要尝试的最短路是 s → ⋯ → u → v s\rightarrow \dots \rightarrow u \rightarrow v s → ⋯ → u → v 。由于我们不考虑重边,因此 x x x 和 u u u 必然是不同的,也就是说即使两条最短路的距离和跳数都相同,终点也相同(显然一定是),倒数第二个点也不同。
当比较两个不同终点的最短路,比如 d ^ ( u ) \hat{d}(u) d ^ ( u ) 和 d ^ ( v ) \hat{d}(v) d ^ ( v ) (u ≠ v u\neq v u = v )时,即使距离和跳数都相同,也能在终点编号这里决出胜负。
这种情况一般是 BaseCase 中的堆以及分块优先队列中使用。
不过可能出现这种情况:要插入的 key 已经在分块优先队列中存在了,此时就需要终点相同的最短路的比较。而算法中都是成功松弛之后才会把 key 插入数据结构,那么根据第一种情况,这个 key 所对应的最短路的倒数第二个节点必然是不一样的。
所以我们可以把终点为 u u u 的最短路视为一个四元组 ( d ( u ) , h o p ( u ) , u , p r e d ( u ) ) (d(u),hop(u),u,pred(u)) ( d ( u ) , h o p ( u ) , u , p re d ( u )) ,其中 d ( u ) d(u) d ( u ) 为最短路,h o p ( u ) hop(u) h o p ( u ) 为最短路跳数,p r e d ( u ) pred(u) p re d ( u ) 为 u u u 最短路的前驱,然后比较的时候直接按这个四元组比较就可以了。维护 h o p ( u ) hop(u) h o p ( u ) 和 p r e d ( u ) pred(u) p re d ( u ) 就跟维护 d ( u ) d(u) d ( u ) 一样。可以 O ( 1 ) O(1) O ( 1 ) 做。
6. 复现
代码在 https://github.com/caiwen666/hnu-algo 的 src/algorithms/bmssp 中:
const_graph.rs 实现了一般图到常度数图的转化。
path_dist.rs 实现了上面所说的路径四元组的比较。整个算法的距离都是用 PathDist 这个类型来表示的。
block_ds.rs 实现了分块优先队列。
mod.rs 实现了算法的主体。
tests/bmssp.rs 提供了对 bmssp 的测试,benches/bmssp.rs 提供了对 bmssp 的性能测试。
复现的时候需要注意:
分块优先队列的 D 1 D_1 D 1 需要保证时刻有一个上界是 B B B 的块。这也就意味着,上界为 B B B 的块内的点如果被删空了,那么是不能删掉这个块的,还需要保留。
论文中分块优先队列的块分裂使用 BFPRT 来找中位数,但显然直接手搓一个 BFPRT 的话常数肯定爆炸了。好在 Rust 标准库里的 select_nth_unstable_by_key 也能平均做到 O ( n ) O(n) O ( n ) 。
同理,Batch Prepend 递归搞 BFPRT 的话常数估计也爆了,所以我这里的实现是直接使用 sort_unstable_by_key 。
其他的地方基本就是照着论文上的伪代码搓。https://github.com/caiwen666/hnu-algo/pull/1 里展示了复现的过程。
复现的第一版代码跑了 3.8s,大概比 Dijkstra 慢了 30 倍。。。显然常数太大了。跑了个火焰图,发现 PathDist 的比较过程是热点,于是把 path_dist.rs 那里的四元组压成了一个 u128 和一个 u32,这样从原来的比较四次优化成了比较两次。同时还稍微搞了下内存复用。不过只是优化到了 3.6s。
然后又跑火焰图,发现计算哈希的开销比较大。最开始的代码用了不少 HashMap 和 HashSet,而 Rust 默认使用的 SipHash 为了保证安全性,计算速度比较慢。我的 HashMap 和 HashSet 的 key 都是节点编号,不太需要高级的哈希算法,于是换成了速度更快但是安全性比较差的 FxHash,把时间优化到了 2.8s。
后来实在没活整了,就直接让 Claude Opus 4.6 给我优化了一下,把时间优化到了 1.6s。我在 review 的时候顺便学习了一下,AI 主要干了这么几件事:
将 PathDist 重新布局
原来是 dis(64) + hop(32) + end(32) 合成一个 u128,然后 pred 单独一个 u32
现在变成了 dis(64) 单独拿出来,hop(22) + end(21) + pred(21) 合成一个 u64
在底层 u128 可能还是分高 64 位和低 64 位两次比较的,而如今两个 u64 应该是稍微好点的。同时由于把 hop、end 和 pred 缩成了 21 位,所以现在算法只能跑在点数在 2e6 以内的图上。
将常度数图变为 CSR 表示
这样的话所有的边在内存中都是连续的,可能可以提高缓存的命中率。
Finding Pivots 在建最短路森林的时候,原来是使用 HashMap<usize, Vec<usize>>,现在是使用链式前向星,减少了内存分配。
把 HashMap 和 HashSet 能砍的全都砍了:
BMSSP 中的 u_set 变为 Vec,然后使用公共的内存池,通过在桶中发放时间戳来判断一个点是否已经在集合中了。
但是由于 BMSSP 是递归的,所以公共内存池可能又被子 BMSSP 占用。AI 利用了子 BMSSP 的 U 是没有交集的这个性质,确保了子 BMSSP 不会破坏父 BMSSP 的 u_set 标记。
Finding Pivots 中的 w_set 也是变为 Vec,然后使用公共的桶来判断一个点是否已经在集合里了。
在生成 p_set 的时候应该没有重复节点,所以 p_set 直接改用 Vec。
BlockDs 里面的 key 到 node_id 的映射,在原来是需要每个 BlockDs 实例里面都建 HashMap。AI 的思路是复用一个桶。这个桶建在 BMSSP 实例里面。
首先第一个问题,这个复用的桶需要在 BlockDs 销毁的时候将其还原。AI 的思路是在 block_ds 里面维护一个 Vec,存储所有插入过的 key,还原的时候就根据这个 Vec 复位。
第二个问题是,父 BMSSP 在调用子 BMSSP 的时候,父 BMSSP 会建一个 BlockDs,子 BMSSP 也会建一个 BlockDs,这就需要我们需要准备多个可以复用的桶。这中间可能遇到 Rust 的 mut 引用检查的阻碍,需要 unsafe。
还有其他的能内存复用的地方就复用。
预留容量,比如 BaseCase 中的 u_set 大小大概是 k + 1 k+1 k + 1 ,那么就让 u_set 的 Vec 预留 k + 2 k+2 k + 2 的空间。
使用占用内存更小的类型。这样主要是减小缓存行的压力
比如点的编号基本都从 usize 变为 u32 了
除此以外,AI 还把要用 Option 的地方直接用 u32::MAX 来代替 None。u32 是 4 字节,但是 Option 需要给加 Tag,对齐完了之后,Option<u32> 就是 8 字节了。
使用洛谷 Cyaron 造了四组数据:
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
from cyaron import *
_n = ati([0 , 10 , 1E3 , 5E5 , 5E5 ])
_m = ati([0 , 30 , 2E4 , 1E6 , 1E6 ])
for i in range (1 , 5 ):
test_data = IO(file_prefix="ssp" , data_id=i)
n = _n[i]
m = _m[i]
s = 1
if i == 4 :
graph = Graph.hack_spfa(n, self_loop=False , directed=True , repeated_edges=False , weight_limit=int (1E9 ))
m = int (1.5 * n)
else :
graph = Graph.graph(n, m, self_loop=False , directed=True , repeated_edges=False , weight_limit=int (1E9 ))
test_data.input_writeln(n, m, s)
test_data.input_writeln(graph)
test_data.output_gen("ssp_std.exe" )
(本来第四组数据想卡 spfa,结果发现洛谷的 Cyaron 的 hack_spfa 是假的,相关 issue 已经挂了一万年了)
在腾讯云 C6.LARGE8 实例上跑了下 benchmark,对比了一下 Dijkstra,Spfa 和 BMSSP:
Case1
Case2
Case3
Case4
Dijkstra
268.89 ns
161.11 µs
106.43 ms
8.4268 ms
Spfa
176.81 ns
141.64 µs
43.333 ms
5.0149 ms
BMSSP
19.683 µs
19.349 ms
1.3813 s
265.75 ms
可以看到,BMSSP 即使理论上比 Dijkstra 的时间复杂度低,但是由于算法过程比较复杂,常数太大,实测还是比 Dijkstra 慢十好几倍。
实际上理论上,不考虑任何常数的影响,BMSSP 也就比 Dijkstra 快 log 1 3 n \log^{\frac{1}{3}} n log 3 1 n 倍,当 n = 10 6 n=10^6 n = 1 0 6 时,理论上也就快 2.7 2.7 2.7 倍,但是 BMSSP 的代码实现的常数肯定比 Dijkstra 高了不止 3 3 3 倍,得到这样的结果是不意外的。但是复现这个论文是算法课作业,于是就图一乐 。