C# でフロイド・ワーシャル アルゴリズムを使用して全ペア最短経路問題を解く

私たちは誰でも、一日に何度も「最短経路」問題を解きます。もちろん、意図せずに。仕事の途中で、インターネットを閲覧しているとき、または机の上で物を整理しているときに、この問題を解きます。

ちょっと…やりすぎのように聞こえますか？調べてみましょう。

想像してください。あなたは友達と、まあ、カフェで会うことに決めたとしましょう。まず、カフェまでのルート（または経路）を見つける必要があるので、利用可能な公共交通機関（徒歩の場合）またはルートと道路（車の場合）を探し始めます。現在地からカフェまでのルートを探していますが、「任意の」ルート（最短または最速のルート）を探しているわけではありません。

もう一つの例を挙げます。これは前の例ほど明白ではありません。経路探索中に、脇道を通る「近道」を取ることにしました。なぜなら、それは「近道」であり、この道を通る方が「速い」からです。しかし、この「近道」の方が速いことをどうやって知ったのでしょうか。個人的な経験に基づいて「最短経路」問題を解き、脇道を通るルートを選択しました。

どちらの例でも、「最短」の経路は、ある場所から別の場所へ移動するために必要な距離または時間で決定されます。移動の例は、グラフ理論、特に「最短経路」問題の非常に自然な応用です。しかし、机の上に物を並べる例はどうでしょうか。この場合、「最短」は、たとえば 1 枚の紙を手に入れるために実行する必要があるアクションの数または複雑さを表します。机を開け、フォルダーを開き、1 枚の紙を取り出すか、机から直接 1 枚の紙を取り出すかです。

上記の例はすべて、グラフ内の 2 つの頂点間の最短経路 (2 つの場所間のルートまたはパス、ある場所から別の場所へ 1 枚の紙を取得するアクションの数または複雑さ) を見つける問題を表しています。このクラスの最短経路問題はSSSP (Single Source Shortest Path)と呼ばれ、これらの問題を解決するための基本アルゴリズムは、計算複雑度がO(n^2)の Dijkstra アルゴリズムです。

しかし、すべての頂点間の最短経路をすべて見つける必要がある場合もあります。次の例を考えてみましょう。自宅、職場、劇場間の定期的な移動のマップを作成しているとします。このシナリオでは、 work ⭢ home 、 home ⭢ work 、 work ⭢ theatre 、 theatre ⭢ work 、 home ⭢ theatre 、 theatre ⭢ homeの 6 つのルートが作成されます (逆のルートは、一方通行などの理由で異なる場合があります)。

マップに場所を追加すると、組み合わせが大幅に増加します。組み合わせ論のn の順列の公式によれば、次のように計算できます。

 A(k, n) = n! / (n - m)! // where n - is a number of elements, // k - is a number of elements in permutation (in our case k = 2)

これにより、4 つの場所に 12 のルート、10 の場所には 90 のルートが得られます (これは印象的です)。注: これは、場所間の中間点を考慮していない場合です。つまり、自宅から職場に行くには、4 つの道路を渡り、川に沿って歩き、橋を渡る必要があります。いくつかのルートに共通の中間点があると想像してください...そうですね...結果として、多数の頂点を持つ非常に大きなグラフが得られ、各頂点は場所または重要な中間点のいずれかを表します。グラフ内のすべての頂点ペア間の最短経路をすべて見つける必要がある問題のクラスは、 APSP (全ペア最短経路)と呼ばれ、これらの問題を解決するための基本アルゴリズムは、 O(n^3)の計算複雑度を持つFloyd-Warshall アルゴリズムです。

これが今日実装するアルゴリズムです 🙂

フロイド・ワーシャルアルゴリズム

フロイド・ワーシャルアルゴリズムは、グラフ内の頂点のペア間の最短経路をすべて見つけます。このアルゴリズムは、ロバート・フロイドによって [1] で公開されました (詳細については、「参考文献」セクションを参照してください)。同じ年に、ピーター・インガーマンは [2] で、3 つのネストされたforループの形式でアルゴリズムの最新の実装を説明しました。

 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 // where W - is a weight matrix of N x N size, // min() - is a function which returns lesser of it's arguments

行列形式で表現されたグラフを扱ったことがない場合は、上記のアルゴリズムが何を行うのか理解するのが難しいかもしれません。そこで、理解を深めるために、グラフを行列形式で表現する方法と、そのような表現が最短経路問題を解決するのになぜ役立つのかを見てみましょう。

下の図は、5 つの頂点を持つ有向の重み付きグラフを示しています。左側のグラフは、円と矢印で構成された視覚的な形式で表示されています。各円は頂点を表し、矢印は方向のあるエッジを表します。円内の数字は頂点番号に対応し、エッジの上の数字はエッジの重みに対応します。右側では、同じグラフが重みマトリックスの形式で表示されています。重みマトリックスは隣接マトリックスの一種で、各マトリックス セルには「重み」、つまり頂点i (行) と頂点j (列) の間の距離が含まれます。重みマトリックスには、頂点間の「パス」に関する情報 ( iからjに移動するときに経由する頂点のリスト) は含まれず、これらの頂点間の重み (距離) のみが含まれます。

重みマトリックスでは、セルの値が隣接する頂点間の重みに等しいことがわかります。そのため、頂点0 (行0 ) からのパスを調べると、 0から1へのパスが 1 つしかないことがわかります。ただし、視覚的に表現すると、頂点0から頂点2と3 (頂点1経由) へのパスが明確にわかります。この理由は単純です。初期状態では、重みマトリックスには隣接する頂点間の距離のみが含まれます。ただし、この情報だけで残りを見つけるのに十分です。

どのように動作するか見てみましょう。セルW[0, 1]に注目してください。その値は、頂点0から頂点1への重みが1のパスがあることを示しています (簡単に言うと、 0 ⭢ 1: 1と記述できます)。この知識があれば、行1のすべてのセル (頂点1からのすべてのパスのすべての重みを含む) をスキャンし、この情報を行0にバックポートして、重みを1 ( W[0, 1]の値) 増やすことができます。

同じ手順を使用して、頂点0から他の頂点までのパスを見つけることができます。検索中に、同じ頂点につながるパスが複数あることがあり、さらに重要なことに、これらのパスの重みが異なる可能性があります。このような状況の例を下の図に示します。頂点0から頂点2までの検索で、重みの小さい頂点3へのパスがもう 1 つ見つかりました。

パスは 2 つあります。元のパス0 ⭢ 3: 4と、先ほど発見した新しいパス0 ⭢ 2 ⭢ 3: 3 (重み行列にはパスが含まれていないため、元のパスにどの頂点が含まれているかはわかりません)。ここで、最短のパスを保持し、セルW[0, 3]に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

k上のfor最も外側のサイクルでは、グラフ内のすべての頂点が反復処理され、各反復で、変数kパスを検索している頂点を表します。 i上のfor最も内側のサイクルでも、グラフ内のすべての頂点が反復処理され、各反復で、 iパスを検索している頂点を表します。最後に、 j上のfor最も内側のサイクルでは、グラフ内のすべての頂点が反復処理され、各反復で、 jパスを検索している頂点を表します。これらを組み合わせると、次のようになります。各反復kで、すべての頂点iから頂点kを経由してすべての頂点jに至るパスを検索します。サイクル内では、パスi ⭢ j ( W[i, j]で表される) とパスi ⭢ k ⭢ j ( W[I, k]とW[k, j]の合計で表される) を比較し、最短のパスをW[i, j]に書き戻します。

さて、仕組みを理解したら、アルゴリズムを実装する時です。

基本的な実装

実装に入る前に、いくつかの技術的な点を明確にする必要があります。

1. すべての実装は、線形配列の形式で表される重み行列で動作します。
2. すべての実装では整数演算を使用します。頂点間のパスが存在しないことは、特別な定数重みNO_EDGE = (int.MaxValue / 2) - 1で表されます。

さて、その理由を考えてみましょう。

1 に関して。マトリックスについて話すとき、セルを「行」と「列」で定義します。このため、マトリックスを「正方形」または「長方形」の形で想像するのが自然なようです (図 4a)。

しかし、これは必ずしも、想像力を働かせるために、行列を配列の配列の形式で表現しなければならないという意味ではありません (図 4b)。その代わりに、各セルのインデックスが次の式を使用して計算される線形配列の形式で行列を表現することができます (図 4c)。

 i = row * row_size + col; // where row - cell row index, // col - cell column index, // row_size - number of cells in a row.

重み行列の線形配列は、フロイド・ワーシャル アルゴリズムを効果的に実行するための前提条件です。その理由は単純ではなく、詳細な説明には別の投稿、またはいくつかの投稿が必要です。ただし、現時点では、このような表現によってデータの局所性が大幅に向上し、実際にアルゴリズムのパフォーマンスに大きな影響を与えることに注意してください。

2 について。アルゴリズムの擬似コードをよく見ると、2 つの頂点間のパスの存在に関連するチェックは見つかりません。代わりに、擬似コードではmin()関数が使用されています。理由は簡単です。元々、2 つの頂点間にパスがない場合、セルの値はinfinityに設定され、JavaScript を除くすべての言語では、すべての値はinfinity未満です。整数の場合、「パスなし」の値としてint.MaxValue使用したくなるかもしれません。ただし、 i ⭢ kパスとk ⭢ jパスの両方の値がint.MaxValueに等しい場合、整数オーバーフローが発生します。そのため、 int.MaxValueの半分より 1 小さい値を使用します。

確かに可能ですが、残念ながらパフォーマンスが大幅に低下します。簡単に言うと、CPU は分岐評価結果の統計を保持します。たとえば、 if文の一部がtrueまたはfalseに評価された場合などです。この統計を使用して、実際のif文が評価される前に「統計的に予測された分岐」のコードを事前に実行します (これを投機的実行と呼びます)。これにより、コードの実行がより効率的になります。ただし、CPU の予測が不正確な場合、CPU は停止して正しい分岐を計算する必要があるため、正しい予測と無条件実行に比べてパフォーマンスが大幅に低下します。

各反復kで重みマトリックスの重要な部分を更新するため、コード パターンがなく、すべての分岐がデータのみに基づいているため、CPU 分岐統計は役に立たなくなります。そのため、このようなチェックでは、分岐予測の誤りが大量に発生します。

ここで私は、（今のところは）私を信じていただき、その後、時間をかけてこのテーマを研究していただきたいと考えています。

うーん、理論的な部分は終わったようです。アルゴリズムを実装しましょう (この実装をBaselineと呼びます)。

 public void Baseline(int[] matrix, 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; } } } } }

上記のコードは、1 つの例外を除いて、前述の疑似コードとほぼ同じコピーです。Math.Min Math.Min()の代わりに、必要な場合にのみセルを更新するためにifを使用します。

理由は簡単です。執筆時点では、JIT はifとMath.Min両方の実装に対してほぼ同等のコードを出力します。詳細はsharplab.ioで確認できますが、メイン ループ本体のスニペットを以下に示します。

 // // == Assembly code of implementation of innermost loop for of j using if. // 53: movsxd r14, r14d // // compare matrix[i * sz + j] and distance (Condition) // 56: cmp [rdx+r14*4+0x10], ebx 5b: jle short 64 // // if matrix[i * sz + j] greater than distance write distance to matrix // 5d: movsxd rbp, ebp 60: mov [rdx+rbp*4+0x10], ebx 64: // end of loop // // == Assembly code of implementation of innermost loop for of j using Math.Min. // 4f: movsxd rbp, ebp 52: mov r14d, [rdx+rbp*4+0x10] // // compare matrix[i * sz + j] and distance (Condition) // 57: cmp r14d, ebx. // // if matrix[i * sz + j] less than distance write value to matrix // 5a: jle short 5e // // otherwise write distance to matrix // 5c: jmp short 61 5e: mov ebx, r14d 61: mov [rdx+rbp*4+0x10], ebx 65: // end of loop

ifとMath.Minどちらを使用しても、条件チェックが行われていることがわかります。ただし、 ifの場合は不要な書き込みはありません。

実装が完了したら、コードを実行して速度を確認してみましょう。

私はBenchmark.NET (バージョン 0.12.1) を使用してコードをベンチマークします。ベンチマークで使用されるすべてのグラフは、300、600、1200、2400、および 4800 頂点の有向非巡回グラフです。グラフ内のエッジの数は、有向非巡回グラフの場合、次のように計算できる最大値の約 80% です。

 var max = v * (v - 1)) / 2; // where v - is a number of vertexes in a graph.

ロックしましょう！

以下は私の PC (Windows 10.0.19042、Intel Core i7-7700 CPU 3.60Ghz (Kaby Lake) / 8 論理プロセッサ / 4 コア) で実行したベンチマークの結果です。

結果からわかるように、計算時間はグラフのサイズに比べて大幅に増加します。300 頂点のグラフでは 27 ミリ秒、2400 頂点のグラフでは 14.5 秒、4800 頂点のグラフでは 119 秒、つまりほぼ 2 分かかります。

アルゴリズムのコードを見ると、計算を高速化するために何かできることがあるとは想像しにくいかもしれません... なぜなら、ループが 3 つあるからです。ループは 3 つだけです。

しかし、よくあることですが、可能性は細部に隠されています 🙂

データを知る – 疎グラフ

フロイド・ワーシャル アルゴリズムは、特に密なグラフや完全なグラフの場合に、すべてのペアの最短経路問題を解決するための基本アルゴリズムです (アルゴリズムはすべての頂点ペア間の経路を検索するため)。

ただし、私たちの実験では、有向非巡回グラフを使用します。このグラフには、頂点1から頂点2へのパスがある場合、頂点2から頂点1へのパスは存在しないという素晴らしい特性があります。私たちにとって、これは、 iからkへのパスがない場合にスキップできる隣接していない頂点が多数あることを意味します (この実装をSpartialOptimisationと呼びます)。

 public void SpartialOptimisation(int[] matrix, int sz) { for (var k = 0; k < sz; ++k) { for (var i = 0; i < sz; ++i) { if (matrix[i * sz + k] == NO_EDGE) { continue; } 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; } } } } }

以下は、同じグラフ セットに対する以前の実装 ( Baseline ) と現在の実装 ( SpartialOptimisation ) の結果です (最速の結果は太字で強調表示されています)。

かなり印象的ですね！

計算時間を半分に短縮しました。もちろん、グラフの密度が高くなるほど、速度の向上は少なくなります。ただし、これは、処理する予定のデータのクラスが事前にわかっている場合に便利な最適化の 1 つです。

もっとできるかな？ 🙂

ハードウェアを理解する – 並列処理

私の PC には、8 つの論理プロセッサ ( HW ) またはハイパースレッディングテクノロジを備えた 4 つのコアを備えたInter Core i7-7700 CPU 3.60GHz (Kaby Lake)プロセッサが搭載されています。コアが複数あるということは、作業に使える「手」が増えるようなものです。では、作業のどの部分を複数のワーカー間で効率的に分割できるかを確認し、それを並列化してみましょう。

ループは、ほとんどの場合、反復が独立しており、同時に実行できるため、常に並列化の最も明白な候補です。アルゴリズムには、分析して並列化を妨げる依存関係があるかどうかを確認する必要がある 3 つのループがあります。

for of kループから始めましょう。各反復で、ループは頂点kを通るすべての頂点からすべての頂点へのパスを計算します。新しいパスと更新されたパスは重み行列に格納されます。反復を見ると、結果に影響を与えずに 0、1、2、3 または 3、1、2、0 の任意の順序で実行できることに気付くかもしれません。ただし、それらは順番に実行する必要があります。そうしないと、一部の反復では、新しいパスまたは更新されたパスが重み行列にまだ書き込まれていないため、それらを使用できなくなります。このような不一致は、間違いなく結果を台無しにします。そのため、引き続き調べる必要があります。

次の候補はfor of iループです。各反復で、ループは頂点iから頂点kへのパス (セル: W[i, k] )、頂点kから頂点jへのパス (セル: W[k, j ]) を読み取り、次にiからjへの既知のパス (セル: W[i, j] ) がi ⭢ k ⭢ jパス (合計: W[i, k] + W[k, j] ) より短いかどうかを確認します。アクセス パターンを詳しく見ると、各反復でiループがk行からデータを読み取り、 i行を更新していることに気付くかもしれません。これは基本的に、反復が独立しており、任意の順序で実行できるだけでなく、並列に実行できることを意味しています。

これは有望そうなので、実装してみましょう (この実装をSpartialParallelOptimisationsと呼びます)。

public void SpartialParallelOptimisations(int[] matrix, int sz) { for (var k = 0; k < sz; ++k) { Parallel.For(0, sz, i => { if (matrix[i * sz + k] == NO_EDGE) { return; } 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; } } }); } }

以下は、同じグラフ セットに対するBaseline 、 SpartialOptimisation 、およびSpartialParallelOptimisations実装の結果です (並列化はParallelクラスを使用して行われます)。

結果から、 for of iループの並列化により、以前の実装 ( SpartialOptimisation ) と比較して計算時間が 2 ～ 4 倍短縮されたことがわかります。これは非常に印象的ですが、純粋な計算の並列化は利用可能なすべての計算リソースを消費し、システム内の他のアプリケーションのリソース不足につながる可能性があることを覚えておくことが重要です。並列化は注意して使用する必要があります。

ご想像のとおり、これで終わりではありません 🙂

プラットフォームを知る – ベクトル化

ベクトル化とは、単一の要素で動作するコードを、複数の要素を同時に動作するコードに変換することです。

これは複雑な問題のように聞こえるかもしれませんが、簡単な例でどのように機能するかを見てみましょう。

 var a = new [] { 5, 7, 8, 1 }; var b = new [] { 4, 2, 2, 6 }; var c = new [] { 0, 0, 0, 0 }; for (var i = 0; i < 4; ++i) c[i] = a[i] + b[i];

非常に単純化した方法では、CPU レベルでの上記のforループの反復0の実行は次のように表すことができます。

ここで何が起こっているか見てみましょう。CPU はメモリから配列aとbの値を CPU レジスタにロードします (1 つのレジスタには 1 つの値しか保持できません)。次に、CPU はこれらのレジストリに対してスカラー加算演算を実行し、その結果をメイン メモリ (配列cに直接) に書き込みます。このプロセスはループが終了する前に 4 回繰り返されます。

ベクトル化とは、特別な CPU レジスタ（ベクトル レジスタまたは SIMD (単一命令複数データ) レジスタ）と、対応する CPU 命令を使用して、複数の配列値に対して一度に操作を実行することを意味します。

ここで何が起こっているか見てみましょう。CPU はメモリから配列aとbの値を CPU レジスタにロードします (ただし、今回は 1 つのベクトル レジスタが2 つの配列値を保持できます)。次に、CPU はこれらのレジストリに対してベクトル加算演算を実行し、結果をメイン メモリ (配列cに直接) に書き込みます。一度に 2 つの値を操作するため、このプロセスは 4 回ではなく 2 回繰り返されます。

.NET でベクトル化を実装するには、 Vector 型とVector<T>型を使用できます ( System.Runtime.Intrinsics名前空間の型も使用できますが、現在の実装では少し高度なため、ここでは使用しませんが、自由に試してみてください)。

 public void SpartialVectorOptimisations(int[] matrix, int sz) { for (var k = 0; k < sz; ++k) { for (var i = 0; i < sz; ++i) { if (matrix[i * sz + k] == 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); } 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; } } } } }

ベクトル化されたコードは少し奇妙に見えるかもしれないので、ステップごとに説明していきましょう。

i ⭢ kパスのベクトルを作成します: var ik_vec = new Vector<int>(matrix[i * sz + k]) 。結果として、ベクトルがint型の値を 4 つ保持でき、 i ⭢ kパスの重みが 5 に等しい場合、4 つの 5 のベクトル ([5, 5, 5, 5]) が作成されます。次に、各反復で、頂点iから頂点j, j + 1, ..., j + Vector<int>.Countへのパスを同時に計算します。

計算プロセス全体は、次の 3 つのステップに分けられます。

1. 重み行列からの情報をベクトルij_vecとikj_vecに読み込みます。
2. ij_vecベクトルとikj_vecベクトルを比較し、最小値を新しいベクトルr_vecに選択します。
3. r_vec重み行列に書き戻します。

1と3 は非常に簡単ですが、 2 は説明が必要です。前述のように、ベクトルでは複数の値を同時に操作します。したがって、一部の操作の結果は 1 つにはなりません。つまり、比較操作Vector.LessThan(ij_vec, ikj_vec)複数の値を比較するため、 trueまたはfalse返すことができません。そのため、代わりに、ベクトルij_vecとikj_vecの対応する値の比較結果を含む新しいベクトルを返します ( ij_vecの値がikj_vecの値より小さい場合は-1 、それ以外の場合は0 )。返されたベクトル (それ自体) はあまり役に立ちませんが、 Vector.ConditionalSelect(lt_vec, ij_vec, ikj_vec)ベクトル操作を使用して、 ij_vecベクトルとikj_vecベクトルから必要な値を抽出し、新しいベクトルr_vecにすることができます。この操作は、入力ベクトルの対応する 2 つの値のうち小さい方に等しい値を持つ新しいベクトルを返します。つまり、ベクトルlt_vecの値が-1に等しい場合、操作はij_vecから値を選択し、それ以外の場合はikj_vecから値を選択します。

2 番の他に、説明が必要な部分がもう 1 つあります。2 番目の非ベクトル化ループです。これは、重み行列のサイズがVector<int>.Count値の積ではない場合に必要です。その場合、メイン ループはすべての値を処理できず (ベクトルの一部をロードできないため)、2 番目の非ベクトル化ループが残りの反復を計算します。

この実装とこれまでのすべての実装の結果は次のとおりです。

結果から明らかなように、ベクトル化により計算時間が大幅に短縮されました。非並列化バージョン ( SpartialOptimisation ) と比較して 3 倍から 4 倍短縮されました。ここで興味深いのは、ベクトル化バージョンは、より小さなグラフ (2400 頂点未満) でも並列バージョン ( SpartialParallelOptimisations ) よりも優れていることです。

そして最後に、ベクトル化と並列処理を組み合わせてみましょう。

プラットフォームとハードウェアを理解しましょう – ベクトル化と並列処理!

最後の実装は並列化とベクトル化の取り組みを組み合わせたもので、正直に言うと個性に欠けています 🙂 なぜなら、 SpartialParallelOptimisationsのParallel.For本体をSpartialVectorOptimisationsのベクトル化されたループに置き換えただけだからです。

 public void SpartialParallelVectorOptimisations(int[] matrix, int sz) { for (var k = 0; k < sz; ++k) { Parallel.For(0, sz, i => { if (matrix[i * sz + k] == NO_EDGE) { return; } 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); } 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; } } }); } }

この投稿のすべての結果は以下に記載されています。予想どおり、並列処理とベクトル化を同時に使用すると最良の結果が得られ、最初の実装と比較して計算時間が最大 25 倍 (1200 頂点のグラフの場合) 短縮されました。

結論

この投稿では、全ペア最短経路問題に少し踏み込んで、それを解決するために C# で Floyd-Warshall アルゴリズムを実装しました。また、データを尊重し、.NET とハードウェアの低レベル機能を利用するように実装を更新しました。

この投稿では、私たちは「オールイン」をしました。しかし、実際のアプリケーションでは、私たちだけではないということを覚えておくことが重要です。ハードコアな並列化は、システム パフォーマンスを大幅に低下させ、メリットよりもデメリットをもたらす可能性があります。一方、ベクトル化は、プラットフォームに依存しない方法で実行すれば、少し安全です。一部のベクトル命令は、特定の CPU でのみ使用できる場合があることに注意してください。

この記事を楽しんでいただき、たった 3 サイクルでアルゴリズムからパフォーマンスを少しだけ「絞り出す」楽しさを感じていただければ幸いです 🙂

参考文献

1. Floyd, RW アルゴリズム 97: 最短経路 / RW Floyd // Communications of the ACM. – 1962. – Vol. 5、№. 6. – P. 345-.
2. Ingerman、PZ アルゴリズム 141: パス マトリックス / PZ Ingerman // Communications of the ACM. – 1962. – Vol. 5、№. 11. – P. 556。