前回の投稿では、フロイド・ワーシャル アルゴリズムを使用して全ペア最短経路問題を解決する方法について説明しました。また、並列処理とベクトル化を使用してアルゴリズムのパフォーマンスを向上させるいくつかの方法についても検討しました。
しかし、フロイド・ワーシャル アルゴリズムの結果について考えると、すぐに興味深い点に気付くでしょう。グラフ内の頂点間の距離 (最短経路の値) は分かっていますが、「ルート」は分かりません。つまり、どの頂点が各最短経路に寄与しているかは分かりません。また、得られた結果からこれらの「ルート」を再構築することもできません。
この投稿では、皆さんと一緒にフロイド・ワーシャル アルゴリズムを拡張してみましょう。今回は、距離を計算して「ルート」を再構築できることを確認します。
コードとベンチマークに進む前に、Floyd-Warshall アルゴリズムの仕組みに関する理論を再確認しましょう。
アルゴリズムのテキストによる説明は次のとおりです。
まず、サイズN
のグラフ ( G
) をサイズN x N
の行列 ( W
) として表します。ここで、各セルW[i,j]
には、頂点i
から頂点j
へのエッジの重みが含まれます (図 1 を参照)。頂点間にエッジがない場合は、セルは特別なNO_EDGE
値に設定されます (図 1 では暗いセルとして表示)。
これからは、 W[i,j]
には頂点i
とj
の間の距離が含まれると言います。
次に、頂点k
を取り、すべての頂点ペアW[i,j]
を反復処理して、距離i ⭢ k ⭢ j
がi
とj
間の現在既知の距離よりも小さいかどうかを確認します。
最小の値はW[i,j]
に格納され、 G
のすべての頂点がk
として使用されるまで、次のk
に対してステップ3が繰り返されます。
疑似コードは次のとおりです。
algorithm FloydWarshall(W) do for k = 0 to N - 1 do for i = 0 to N - 1 do for j = 0 to N - 1 do W[i,j] = min(W[i,j], W[i,k] + W[k,j]) end for end for end for end algorithm
完了すると、行列W
のすべてのセルW[i,j]
には、頂点i
とj
間にパスがない場合はNO_EDGE
値が含まれます。
冒頭で述べたように、この情報だけでは、これらの頂点間のルートを再構築することは不可能です。では、これらのルートを再構築するにはどうすればよいでしょうか?
まあ…当たり前のように聞こえるかもしれませんが…もっとデータを収集する必要があります!
フロイド・ワーシャルアルゴリズムの本質は、中間頂点k
を通るi
からj
までの潜在的に新しいパスの距離を持つ既知の距離W[i,j]
の区画です。
この詳細にこれほど注意を向けるのは、それがルート情報をどのように保存できるかという鍵となるからです。先ほどの 5 つの頂点のグラフ (図 2 を参照) を取り、これに対して Floyd-Warshall アルゴリズムを実行してみましょう。
最初は、すべてのグラフ エッジがわかっているので、次のパスが得られます: 0 ⭢ 1 :2
0 ⭢ 4 :10
1 ⭢ 3 :1
2 ⭢ 4 :1
3 ⭢ 2 :1
および3 ⭢ 4 :3
。
また、頂点0
につながるエッジがないことにもわかっているので、 k = 0
のアルゴリズムを実行しても意味がありません。頂点0
から頂点1
につながるエッジ ( 0 ⭢ 1
) が 1 つあることも明らかです。これにより、すべてのi != 0
(ここではi
「from」です) の実行もまったく意味がなくなります。また、頂点1
には2
と4
とのエッジがあるため、 2
でも4
でもないj
(ここではj
は「to」です) に対してアルゴリズムを実行することも意味がありません。
そのため、最初のステップでは、 k = 1
、 i = 0
、 j = 2,4
のアルゴリズムを実行します。
ステップ | パス | コメント |
---|---|---|
1 | 0 ⭢ 1 ⭢ 2 | パスが見つかりました。距離 = 3 (何もありませんでした) |
2 | 0 ⭢ 1 ⭢ 4 | パスが見つかりました。距離 = 8 (以前は 10)。 |
2 つのパスが見つかりました。新しいパス ( 0 ⭢ 1 ⭢ 2
) とショートカット ( 0 ⭢ 1 ⭢ 4
) です。どちらも頂点1
を通ります。この情報 ( 1
通って2
と4
に到達したという事実) を今すぐどこかに保存しないと、その情報は永久に失われます (これは私たちが望んでいることとはまったく逆です)。
では、どうすればよいでしょうか。行列W
新しい値で更新し (図 3a を参照)、行列W
と同じサイズだがNO_EDGE
値で初期化された新しい行列R
のセルR[0,2]
とR[0,4]
にk
の値 ( k = 1
) を格納します (図 3b を参照)。
今は、行列R
が正確に何であるかに焦点を当てないでください。そのまま続けて、次のk = 2
のアルゴリズムを実行してみましょう。
ここでは、 k = 1,
の場合と同じ分析 (どの組み合わせを実行するのが理にかなっているかを理解するため) を行いますが、今回はグラフG
の代わりに行列W
使用します。行列W
、特に列 #2 と行 #2 に注目してください (図 4)。
ここで、頂点0
と頂点1
2
の 2 つのパスと、頂点2
から頂点3
と頂点4
(行 #2) への 2 つのパスがあることがわかります。これを知ると、 k = 2
、 i = 0,1
、 j = 3,4
の組み合わせに対してのみアルゴリズムを実行することが理にかなっています。
ステップ | パス | コメント |
---|---|---|
1 | 0 ⭢ 2 ⭢ 3 | パスが見つかりました。距離 = 4 (何もありませんでした) |
2 | 0 ⭢ 2 ⭢ 4 | パスが見つかりました。距離 = 6 (以前は 8) |
3 | 1 ⭢ 2 ⭢ 3 | パスが見つかりました。距離 = 2 (何もありませんでした) |
4 | 1 ⭢ 2 ⭢ 4 | パスが見つかりました。距離 = 4 (以前は 6) |
すでに行ったように、セルW[0,3]
、 W[0,4]
、 W[1,3]
、 W[1,4]
新しい距離で更新し、セルR[0,3]
、 R[0,4]
、 R[1,3]
、 R[1,4]
をk = 2
で更新します(図5を参照)。
処理すべき残りの部分はk = 3
だけです (グラフ内に頂点4
から他の頂点につながるエッジがないため)。
もう一度、行列W
(図 6) を見てみましょう。
行列W
によれば、頂点0
(列 # 2
) から頂点3
へのパスは 3 つあり、頂点3
から頂点4
1
行 #3) へのパスは 1 つあります。したがって、処理するパスは次のようになります。
ステップ | パス | コメント |
---|---|---|
1 | 0 ⭢ 3 ⭢ 4 | パスが見つかりました。距離 = 5 (以前は 6) |
2 | 1 ⭢ 3 ⭢ 4 | パスが見つかりました。距離 = 3 (以前は 4) |
3 | 2 ⭢ 3 ⭢ 4 | パスが見つかりました。距離 = 2 (以前は 3) |
これはアルゴリズムの最後の反復でした。残っているのは、セルW[0,4]
、 W[1,4]
、 W[2,4]
新しい距離で更新し、セルR[0,4]
、 R[1,4]
、 R[2,4]
を3
に設定することだけです。
アルゴリズムの最後には次のようになります (図 7)。
前回の投稿からわかるように、行列W
にはグラフG
の最短経路のペアがすべて含まれています。しかし、行列R
内には何が格納されているのでしょうか?
新しい最短経路が見つかるたびに、行列R
の対応するセルをk
の現在の値で更新していました。最初はこのアクションは不思議に思えるかもしれませんが、意味は非常に単純です。つまり、目的の頂点に到達した頂点を保存していたのです: i ⭢ k ⭢ j
(ここではk
から直接j
に到達しています)。
この瞬間は非常に重要です。なぜなら、どこから来たかがわかれば、頂点j
(「to」) から頂点i
(「from」) まで (ロブスターのように) 逆方向に移動してルートを再構築できるからです。
以下は、 i
からj
へのルートを再構築するアルゴリズムのテキストによる説明です。
X
を準備します。R[i,j]
から値z
読み取ります。z
がNO_EDGE
の場合、 i
からj
へのルートが見つかったので、ステップ 7 に進みます。X
の先頭にz
追加します。R[i,z]
の値をz
に読み込みます。i
追加します。j
X
に追加します。X
i
からj
へのルートが含まれるようになりました。上記のアルゴリズムは既存のパスに対してのみ機能することに注意してください。また、パフォーマンスの観点からも完璧ではなく、最適化できることは確かです。ただし、この投稿の範囲内で、理解しやすく、紙に再現しやすいように具体的に説明しています。
- 著者注
疑似コードではさらにシンプルになります。
algorithm RebuildRoute(i, j, R) x = array() z = R[i,j] while (z ne NO_EDGE) do x.prepend(z) z = R[i,z] end while x.prepend(i) x.append(j) return x end algorithm
グラフG
で試して、頂点0
から頂点4
までのルートを再構築してみましょう (図 8 を参照)。
上記の図のテキストによる説明は次のとおりです。
まず、 R[0,4]
から値を読み取ることから始めます。結果は3
になります。アルゴリズムによれば、これはルートが頂点3
から頂点4
に行くことを意味します (青色で表示)。
R[0,4]
の値はNO_EDGE
ではないため、先に進んでR[0,3]
の値を読み取ります。結果は2
(緑色で表示)になります。
ここでも、 R[0,3]
の値はNO_EDGE
ではないため、先に進んでR[0,2]
の値を読み取ります。結果は1
になります(赤で表示)。
最後に、 R[0,1]
の値を読み取ります。これはNO_EDGE
値になります。つまり、ルートに寄与する頂点は0
と4
以外にはありません。したがって、結果のルートは0 ⭢ 1 ⭢ 2 ⭢ 3 ⭢ 4
となり、グラフを見ると、これは確かに頂点0
から頂点4
への最短経路に対応しています。
思慮深い読者はこう尋ねるかもしれません。
行列R
から読み取ったすべてのデータが同じパスに属していることをどのように確認できますか?
- 思慮深い読者
非常に良い質問です。行列W
の最短経路値を更新するときに、行列R
新しい値で更新するので、確実です。したがって、最終的には、行列R
のすべての行に最短経路に関連するデータが含まれます。それ以上でもそれ以下でもありません。
さて、理論は終わりました。次は実装の時間です。
前回の投稿では、Floyd-Warshall アルゴリズムのオリジナル バージョンを実装するだけでなく、スパース グラフ、並列処理、ベクトル化のサポートなど、さまざまな最適化も統合し、最終的にこれらすべてを組み合わせました。
今日はこれをすべて繰り返す必要はありません。ただし、ルート記憶をアルゴリズムのオリジナル バージョンとベクトル化バージョン (スパース グラフをサポート) に統合する方法を紹介したいと思います。
信じ難いかもしれませんが、ルート記憶をアルゴリズムのオリジナルバージョンに統合するには、次の操作を行うだけです。
関数シグネチャを拡張して、 R
マトリックスを別のパラメーターとして含めます ( int[] routes
。
最短経路が変更されるたびに、 routes
に k を保存します (行: 2 および 14)。
以上です。コードはたった 1 行半です。
public void BaselineWithRoutes( int[] matrix, int[] routes, int sz) { for (var k = 0; k < sz; ++k) { for (var i = 0; i < sz; ++i) { for (var j = 0; j < sz; ++j) { var distance = matrix[i * sz + k] + matrix[k * sz + j]; if (matrix[i * sz + j] > distance) { matrix[i * sz + j] = distance; routes[i * sz + j] = k; } } } } }
ベクトル化されたバージョン(スパース グラフのサポート付き)への統合には、完了までにもう少し労力がかかりますが、概念的な手順は同じです。
関数シグネチャを拡張して、 R
マトリックスを別のパラメーターとして含めます ( int[] routes
。
各反復で、 k
値の新しいベクトルを初期化します (行: 6)。
最短経路が変更されるたびに、 k
値のベクトルをroutes
に保存します (行: 31-32)。
元のアルゴリズムを更新したのと同じ方法で、アルゴリズムのベクトル化されていない部分を更新します (行: 41)。
public void SpartialVectorOptimisationsWithRoutes( int[] matrix, int[] routes, int sz) { for (var k = 0; k < sz; ++k) { var k_vec = new Vector<int>(k); for (var i = 0; i < sz; ++i) { if (matrix[i * sz + k] == Constants.NO_EDGE) { continue; } var ik_vec = new Vector<int>(matrix[i * sz + k]); var j = 0; for (; j < sz - Vector<int>.Count; j += Vector<int>.Count) { var ij_vec = new Vector<int>(matrix, i * sz + j); var ikj_vec = new Vector<int>(matrix, k * sz + j) + ik_vec; var lt_vec = Vector.LessThan(ij_vec, ikj_vec); if (lt_vec == new Vector<int>(-1)) { continue; } var r_vec = Vector.ConditionalSelect(lt_vec, ij_vec, ikj_vec); r_vec.CopyTo(matrix, i * sz + j); var ro_vec = new Vector<int>(routes, i * sz + j); var rk_vec = Vector.ConditionalSelect(lt_vec, ro_vec, k_vec); rk_vec.CopyTo(routes, i * sz + j); } for (; j < sz; ++j) { var distance = matrix[i * sz + k] + matrix[k * sz + j]; if (matrix[i * sz + j] > distance) { matrix[i * sz + j] = distance; routes[i * sz + j] = k; } } } } }
ちょっとした注意点 – Vector.ConditionalSelect
操作は、入力ベクトルの対応する 2 つの値のうち小さい方に等しい値を持つ新しいベクトルを返します。つまり、ベクトルlt_vec
の値が-1
に等しい場合、操作はij_vec
から値を選択し、それ以外の場合はk_vec
から値を選択します。
- 著者注
ルート情報を収集するのは、十分に合理的に思えます。なぜなら、既存のアルゴリズムに簡単に統合できるからです。ただし、デフォルトで統合することがどれだけ実用的かを見てみましょう。
ここに、ありとなしの 2 セットのベンチマークがあります (どちらもアルゴリズムのオリジナル バージョンとベクトル化されたバージョンを実行します)。すべてのベンチマークは、Intel Core i5-6200U CPU 2.30GHz (Skylake) プロセッサを搭載し、Windows 10 を実行しているマシンで、 BenchmarkDotNet v0.13.1 (.NET 6.0.4 x64) によって実行されました。
すべての実験は、 前回の投稿で使用した定義済みのグラフ セット (300、600、1200、2400、および 4800 頂点) で実行されました。
ソース コードと実験グラフは、GitHub のリポジトリで入手できます。実験グラフは、 Data/Sparse-Graphs.zip
ディレクトリにあります。この記事のすべてのベンチマークは、 APSP02.csファイルに実装されています。
以下はベンチマーク結果です。Baseline およびBaselineWithRoutes
メソッドはアルゴリズムのオリジナル バージョンに対応し、 SpartialVectorOptimisations
およびSpartialVectorOptimisationsWithRoutes
メソッドBaseline
アルゴリズムのベクトル化バージョン (スパース グラフのサポート付き) に対応します。
方法 | サイズ | 平均(ミリ秒) | エラー (ミリ秒) | 標準偏差 (ミリ秒) |
---|---|---|---|---|
ベースライン | 300 | 40.233 | 0.0572 | 0.0477 |
ルート付きベースライン | 300 | 40.349 | 0.1284 | 0.1201 |
部分ベクトル最適化 | 300 | 4.472 | 0.0168 | 0.0140 |
ルート付き部分ベクトル最適化 | 300 | 4.517 | 0.0218 | 0.0193 |
ベースライン | 600 | 324.637 | 5.6161 | 4.6897 |
ルート付きベースライン | 600 | 321.173 | 1.4339 | 1.2711 |
部分ベクトル最適化 | 600 | 32.133 | 0.2151 | 0.1679 |
ルート付き部分ベクトル最適化 | 600 | 34.646 | 0.1286 | 0.1073 |
ベースライン | 1200 | 2,656.024 | 6.9640 | 5.8153 |
ルート付きベースライン | 1200 | 2,657.883 | 8.8814 | 7.4164 |
部分ベクトル最適化 | 1200 | 361.435 | 2.5877 | 2.2940 |
ルート付き部分ベクトル最適化 | 1200 | 381.625 | 3.6975 | 3.2777 |
ベースライン | 2400 | 21,059.519 | 38.2291 | 33.8891 |
ルート付きベースライン | 2400 | 20,954.852 | 56.4719 | 50.0609 |
部分ベクトル最適化 | 2400 | 3,029.488 | 12.1528 | 11.3677 |
ルート付き部分ベクトル最適化 | 2400 | 3,079.006 | 8.6125 | 7.1918 |
ベースライン | 4800 | 164,869.803 | 547.6675 | 427.5828 |
ルート付きベースライン | 4800 | 184,305. 980 | 210.9535 | 164.6986 |
部分ベクトル最適化 | 4800 | 21,882.379 | 200.2820 | 177.5448 |
ルート付き部分ベクトル最適化 | 4800 | 21,004.612 | 79.8752 | 70.8073 |
ベンチマーク結果の解釈は簡単ではありません。詳しく調べると、ルート収集によって実際にアルゴリズムの実行速度が速くなった場合と、まったく大きな影響がなかった場合の組み合わせが見つかります。
これは非常に混乱しているように見えます (もしそうなら、デニス・バクヴァロフとアンドレイ・アキンシンが出演するこの 番組を聞いて、計測にどのようなトリッキーな要素が影響するかをよりよく理解することを強くお勧めします)。ここでの私の最善の見解は、「混乱」する動作は優れた CPU キャッシュ パフォーマンスによって発生するということです (グラフがキャッシュを飽和させるほど大きくないため)。部分的に、この理論は、最大のグラフで大幅な劣化が見られる「太字」のサンプルに基づいています。ただし、私はそれを検証していません。
一般的に、ベンチマークでは、極めて高性能なシナリオや巨大なグラフについて話しているのではない限り、デフォルトで両方のアルゴリズムにルート記憶を統合しても問題ないことを示しています (ただし、1 つの行列W
とR
ではなく 2 つの行列を割り当てる必要があるため、メモリ消費量が 2 倍になることに注意してください)。
残っているのはルート再構築アルゴリズムの実装だけです。
C# でのルート再構築アルゴリズムの実装は簡単で、これまでに示した疑似コードをほぼ完全に反映しています (動的配列を表すためにLinkedList<T>
を使用します)。
public static IEnumerable<int> RebuildWithLinkedList( int[] routes, int sz, int i, int j) { var x = new LinkedList<int>(); var z = routes[i * sz + j]; while (z != Constants.NO_EDGE) { x.AddFirst(z); z = routes[i * sz + z]; } x.AddFirst(i); x.AddLast(j); return x; }
上記のアルゴリズムは完璧ではなく、間違いなく改善できます。最も簡単な改善方法は、 sz
サイズの配列を事前に割り当てて、逆の順序で埋めることです。
public static IEnumerable<int> RebuildWithArray( int[] routes, int sz, int i, int j) { var x = new int[sz]; var y = sz - 1; // Fill array from the end x[y--] = j; var z = routes[i * sz + j]; while (z != Constants.NO_EDGE) { x[y--] = z; z = routes[i * sz + z]; } x[y] = i; // Create an array segment from 'y' to the end of the array // // This is required because 'sz' (and therefore length of the array x) // equals to the number of vertices in the graph // // So, in cases when route doesn't span through // all vertices - we need return a segment of the array return new ArraySegment<int>(x, y, sz - y); }
この「最適化」により割り当て回数は 1 回に減りますが、必ずしもアルゴリズムが「高速化」されたり、「割り当てるメモリが少なくなる」わけではありません。ここでのポイントは、 i
からj
へのルートを順序付ける必要がある場合、常に行列R
から取得したデータを「反転」する必要があり、それを回避する方法はないということです。
しかし、データを逆の順序で表示できれば、LINQ を利用して不要な割り当てを回避できます。
public static IEnumerable<int> RebuildWithReverseYield( int[] routes, int sz, int i, int j) { yield return j; var z = routes[i * sz + j]; while (z != Constants.NO_EDGE) { yield return z; z = routes[i * sz + z]; } yield return i; }
このコードを「変更」または「改善」する方法には、さらに多くのバリエーションがあります。ここで重要なのは、どのような「変更」にも、メモリや CPU サイクルの面でマイナス面があることを覚えておくことです。
ぜひ自由に実験してみてください。可能性はほぼ無限です。
すべてのルート アルゴリズムの実装は、GitHub のRoutes.csファイルで確認できます。
この投稿では、フロイド・ワーシャル アルゴリズムの背後にある理論を深く掘り下げ、最短経路「ルート」を記憶できるように拡張しました。次に、 前回の投稿の C# 実装 (オリジナルとベクトル化) を更新しました。最後に、アルゴリズムのいくつかのバージョンを実装して、「ルート」を再構築しました。
これまで何度も繰り返してきましたが、同時に、この中で何か新しい興味深いものを見つけられたことを願っています。これは、全ペア最短経路問題の終わりではありません。これは始まりに過ぎません。
楽しんで読んでいただけたなら幸いです。また次回お会いしましょう!