閉路がある (=木ではない) グラフの最短路も、幅優先探索を使い高速で求められる...!?

このAtCoderの過去問を解いた際、ダイクストラ法を使いACはできたが、Python (PyPy) だとどうしても実行時間制限 (2000ms) ギリギリになってしまう。
過去の提出を探したところ、1000ms超えが多い中、789msという高速な解答を見つけた。

木に対する幅優先探索による最短経路計算のようなことをしているが、グラフが木でなく閉路がある場合、時間計算量は $O(V^2)$ であるように見える (木であるなら $O(V) = O(E)$ ) 。
しかし実際には高速で答えが出ている。
計算量の推定が間違っている可能性も考え、以下の3つの方法を比較してみた。

dijkstra0 : 上記リンクの、幅優先探索のような方法 (一部改変, $O(V^2)$ のはず)
dijkstra1 : 優先度付きキューを使ったダイクストラ法 ( $O((E+V)\log{V})$ )
dijkstra2 : 最も基本的なダイクストラ法 ( $O(V^2)$ )

参考: Wikipediaの記事の擬似コード

時間計算量だけを考えると、多くの場合で dijkstra1 が最も速そうだ。

問題設定

連結な無向グラフ (頂点数 $V$ 、辺数 $E$ 、多重辺なし、自己ループなし) を用いる。
最初に下図のような、始点から順に全ての点を1回ずつ通る一本道を作っておくことで、連結を担保する (道中の頂点の順序はランダム) 。

その他の辺は、まだ辺のない頂点間にランダムに張る。
辺の重みは、各辺について $1$ ～ $10^7$ からランダムに決める。

$V \in \{ 10^3,~ 10^4,~ 10^5,~ 10^6,~ 10^7 \} $
$E \in \{ 10^3,~ 10^4,~ 10^5,~ 3\cdot 10^5,~ 10^6,~ 10^7 \} $
(辺数に $3\cdot 10^5$ を入れたのは競プロの制約を想定したため)

$V-1 \leqq E \leqq \frac{V(V-1)}{2}$ を満たす $(V, E)$ について、時間を測定した。

コード

worstcase.py

def dijkstra0(a, in_connect_list, v, money_kind): 
    out_shortest_list = [10**15 for i in range(v)]   # 頂点aからの距離
    out_shortest_list[a] = 0
    searching_list = [a]
    
    while searching_list != []:
        new_search_list = []
        for i in searching_list:
            for j in in_connect_list[i]:
                c = j[0]
                p = j[money_kind] 
                if out_shortest_list[c] > p + out_shortest_list[i]:
                    out_shortest_list[c] = p + out_shortest_list[i]
                    new_search_list.append(c)
        searching_list = list(set(new_search_list))  # 元コードから変更
    return out_shortest_list


def dijkstra1(s, M):
    D = [10**15] * v  # 頂点sからの距離
    D[s] = 0
    V = [0] * v  # その頂点のD[i]が最短距離と確定したら1
    Q = []  # 優先度付きキュー
    for v in range(v):
        heapq.heappush(Q, [D[v], v])

    le = len(Q)
    while le > 0:
        q = heapq.heappop(Q)
        u = q[1]
        du = q[0]
        if V[u] == 0:
            V[u] = 1
            le -= 1
            for i in range(len(M[u])):
                v = M[u][i][0]
                luv = M[u][i][1]
                if V[v] == 0:
                    alt = du + luv
                    if D[v] > alt:
                        D[v] = alt
                        heapq.heappush(Q, [alt, v])
    return D


def dijkstra2(s, M):
    D = [10**15] * v  # 頂点sからの距離
    D[s] = 0
    Q = [i for i in range(v)]
    
    while len(Q) > 0:
        mi_ind = -1
        mi = 10**15
        for i in Q:
            if D[i] < mi:
                mi = D[i]
                mi_ind = i
        Q.pop(Q.index(mi_ind))
        for i in range(len(M[mi_ind])):
            v = M[mi_ind][i][0]
            luv = M[mi_ind][i][1]
            D[v] = min(D[v], D[mi_ind] + luv)
    return D


from random import randint, random, sample
import heapq
import time

V = [10**3, 10**4, 10**5, 10**6, 10**7]
E = [10**3, 10**4, 10**5, 3 * 10**5, 10**6, 10**7] 

for v in V:
    for e in E:
        if v-1 <= e <= v*(v-1)//2:  # 連結かつ多重辺が存在しないための条件
            T = []
            for _ in range(10):  # 10回測定した平均をとる
                
                N_shuf = [0] + sample([i for i in range(1, v)], v-1)
                M = [[] for i in range(v)]

                A = [min(N_shuf[i],N_shuf[i+1]) * v
                     + max(N_shuf[i],N_shuf[i+1]) for i in range(v-1)]
                while len(A) < e:
                    e_add = e - len(A)
                    ADD = [0] * e_add
                    i = 0
                    while i < e_add:
                        v0 = randint(0, v-1)
                        v1 = randint(0, v-1)
                        if v0 != v1:
                            ADD[i] = min(v0,v1) * v + max(v0,v1)
                            i += 1
                    A += ADD
                    A = list(set(A))  # 多重辺を除去

                for ii in A:
                    i = ii // v
                    j = ii % v
                    M[i].append([j, randint(1, 10**7)])
                    M[j].append([i, randint(1, 10**7)])

                t0 = time.time()
                
                D = dijkstra0(0, M, v, 1)
                # D = dijkstra1(0, M)
                # D = dijkstra2(0, M)

                t1 = time.time() - t0
                T.append(t1)
            print(v, e, sum(T) / 10)

結果

10回行った平均、単位は秒、有効数字2桁とした。

dijkstra0 : 幅優先探索のような方法 ( $O(V^2)$ )

$V$ \ $E$	$10^3$	$10^4$	$10^5$	$3\cdot 10^5$	$10^6$	$10^7$
$10^3$	0.00060	0.0074	0.14	0.40
$10^4$		0.0069	0.29	0.96	3.1	25
$10^5$			0.13	1.0	4.4	53
$10^6$					1.6	75
$10^7$						18

dijkstra1 : 優先度付きキューを使ったダイクストラ法 ( $O((E+V)\log{V})$ )

$V$ \ $E$	$10^3$	$10^4$	$10^5$	$3\cdot 10^5$	$10^6$	$10^7$
$10^3$	0.0016	0.0059	0.070	0.23
$10^4$		0.020	0.16	0.41	1.1	8.9
$10^5$			0.36	1.3	2.6	15
$10^6$					4.4	37
$10^7$						53

dijkstra2 : 最も基本的なダイクストラ法 ( $O(V^2)$ )

$V$ \ $E$	$10^3$	$10^4$	$10^5$	$3\cdot 10^5$	$10^6$	$10^7$
$10^3$	0.017	0.023	0.11	0.29
$10^4$		1.9	2.1	2.2	3.1	14
$10^5$			180	-	-	-
$10^6$					-	-
$10^7$						-
- : 数百秒以上であったため測定せず

考察

最も基本的なダイクストラ法は、$V$ が増えてくると明らかに遅くなる。
しかし、 dijkstra0 と優先度付きキューを使った方法 ( dijkstra1 ) には、表の範囲内では大きな差はなさそうだ。
競プロにありそうな条件 ( $V \leqq 10^5,~ E \leqq 3 \cdot 10^5$ ) ならば、どちらでも通るか、場合によっては dijkstra0 の方が速そうでもある (下の表で時間の比較を記載した) 。
$E$ を固定して $V$ を大きくすると、むしろ速くなる場合があるというのも dijkstra0 の特徴だ。
どうやら、 dijkstra0 はスパースで木に近い ( = $V$ が大きい割に $E$ が小さい) グラフを得意とするようである。木に近いとほぼ $O(V)$ になるのだろう。

下表: 所要時間の比 ( dijkstra0 / dijsktra1 )

$V$ \ $E$	$10^3$	$10^4$	$10^5$	$3\cdot 10^5$	$10^6$	$10^7$
$10^3$	0.38	1.3	2.0	1.7
$10^4$		0.35	1.8	2.4	2.8	2.8
$10^5$			0.37	0.82	1.7	3.6
$10^6$					0.36	2.0
$10^7$						0.34

最悪ケースだと…

dijkstra0 の最悪ケースとして、以下のようなグラフが考えられる。

最短距離の更新回数が、一番左上の頂点は $(V-1)$ 回、その1つ右の頂点は $(V-2)$ 回、...、一番右上の頂点は $1$ 回となり、合計 $\frac{V(V-1)}{2}$ 回の更新が行われるため、計算時間も $V^2$ にほぼ比例して増大していくだろう。
ダイクストラ法の関数部分以外のコードを以下のように書き換え、$V \in \{ 10^3,~ 5\cdot 10^3,~ 10^4,~ 5\cdot 10^4\}$ で比較してみた。

worstcase.py

from random import randint, random, sample
import heapq
import time

for v in [1000, 5000, 10000, 50000]:
        T = []
        for _ in range(10):

            N_shuf = [0] + sample([i for i in range(1, v)], v-1)
            M = [[] for i in range(v)]

            for i in range(1, v):
                M[0].append([N_shuf[i], 10**7 - i * 2])
                M[N_shuf[i]].append([0, 10**7 - i * 2])
            for i in range(1, v-1):
                M[N_shuf[i+1]].append([N_shuf[i], 1])
                M[N_shuf[i]].append([N_shuf[i+1], 1])

            t0 = time.time()
            D = dijkstra0(0, M, v, 1)
            # D = dijkstra1(0, M)
            # D = dijkstra2(0, M)
            t1 = time.time() - t0
            T.append(t1)
        print(v, sum(T) / 10)

$V$	$E$	dijkstra0	dijkstra1	dijkstra2
1000	1997	0.11	0.0024	0.017
5000	9997	3.1	0.016	0.45
10000	19997	12	0.033	1.8
50000	99997	740	0.25	48

やはり dijkstra0 では $V$ が10倍になると計算時間は100～200倍になり ( これは dijkstra2 も同様) 前述のランダムなグラフの場合とは全く違う結果となっている。
$O((E+V)\log{V})$ である dijkstra1 はやはり安定して速い。
冒頭の問題のようにテストケース内にこの最悪ケースが入っていない場合もあるだろうが、競プロでは dijkstra1 を使うのが良さそうだ。
時間制限ギリギリの場合は、ダイクストラ法以外の部分を含め、速く動く実装をすることで対処するのが良いだろう。

まとめ

木に対する幅優先探索を閉路のあるグラフに用いると、時間計算量は $O(V^2)$ だが、ダイクストラ法を工夫したもの ( $O((E+V)\log{V})$ ) より速い場合がある。
競プロで出題される程度のサイズのグラフに対しては、通ってしまうこともある。
しかし最悪ケースでは遅いので、やはりダイクストラ法を使ったほうが良い。

更新履歴

2020/04/20 記事タイトル変更 (旧題：O(V^2)でもダイクストラより速い？閉路があるグラフに幅優先探索)、一部本文修正