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

私たちは誰でも、一天に何度も「」問題を解きます。もちろん、意図せずに。仕事の过程中で、インターネットを閲覧しているとき、または机の上で物を归置しているときに、この問題を解きます。

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

假如してください。あなたは友達と、まあ、カフェで会うことに決めたとしましょう。まず、カフェまでのルート（または経路）を見つける需要があるので、使用几率な公共信息路网機関（徒歩の場合）またはルートと市政道路（車の場合）を探し始めます。現在地からカフェまでのルートを探していますが、「随便の」ルート（比较短または最速のルート）を探しているわけではありません。

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

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

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

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

マップに場所を追加すると、組み合わせが急剧に増加します。論のによれば、次のように計算できます。

 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)の計算複雑度を持つです。

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

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

フロイド・ワーシャルアルゴリズムは、グラフ内の頂点のペア間の最短経路をすべて見つけます。このアルゴリズムは、ロバート・フロイドによって [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両方の実装に対してほぼ同等のコードを出力します。詳細はで確認できますが、メイン ループ本体のスニペットを以下に示します。

 // // == 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の場合は不要な書き込みはありません。

実装が弄完したら、コードを実行して速度快を確認してみましょう。私は (バージョン 0.12.1) を的动用してコードをベンチマークします。ベンチマークで的动用されるすべてのグラフはグラフです。グラフ内のエッジの数は、有向非巡回グラフの場合、次のように計算できる大値の約 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.525 ミリ秒	0.1937 ミリ秒	0.1617ミリ秒
ベースライン	600	217.897 ミリ秒	1.6415 ミリ秒	1.5355 ミリ秒
ベースライン	1200	1,763.335 ミリ秒	7.4561 ミリ秒	6.2262 ミリ秒
ベースライン	2400	14,533.335 ミリ秒	63.3518 ミリ秒	52.9016 ミリ秒
ベースライン	4800	119,768.219 ミリ秒	181.5514 ミリ秒	160.9406 ミリ秒

結果からわかるように、計算時間はグラフのサイズに比べて大幅に増加します。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 ) の結果です (最速の結果は太字で強調表示されています)。

步骤	サイズ	平衡	エラー	標準差别	百分率
ベースライン	300	27.525 ミリ秒	0.1937ミリ秒	0.1617ミリ秒	1.00
部分最適化	300	12.399ミリ秒	0.0943 ミリ秒	0.0882 ミリ秒	0.45
ベースライン	600	217.897 ミリ秒	1.6415ミリ秒	1.5355 ミリ秒	1.00
部分最適化	600	99.122 ミリ秒	0.8230 ミリ秒	0.7698 ミリ秒	0.45
ベースライン	1200	1,763.335 ミリ秒	7.4561 ミリ秒	6.2262 ミリ秒	1.00
部分最適化	1200	766.675 ミリ秒	6.1147 ミリ秒	5.7197 ミリ秒	0.43
ベースライン	2400	14,533.335 ミリ秒	63.3518 ミリ秒	52.9016 ミリ秒	1.00
部分最適化	2400	6,507.878 ミリ秒	28.2317 ミリ秒	26.4079 ミリ秒	0.45
ベースライン	4800	119,768.219 ミリ秒	181.5514 ミリ秒	160.9406 ミリ秒	1.00
部分最適化	4800	55,590.374 ミリ秒	414.6051 ミリ秒	387.8218 ミリ秒	0.46

かなり想象的ですね！

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

もっとできるかな？ 🙂

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

---

私の PC には、8 つの論理プロセッサ ( ) またはテクノロジを備えた 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実装の結果です (並列化はクラスを使用して行われます)。

最简单的方法	サイズ	年均	エラー	標準误差	占比
ベースライン	300	27.525 ミリ秒	0.1937ミリ秒	0.1617ミリ秒	1.00
一些最適化	300	12.399ミリ秒	0.0943 ミリ秒	0.0882 ミリ秒	0.45
部分並列最適化	300	6.252 ミリ秒	0.0211ミリ秒	0.0187ミリ秒	0.23
ベースライン	600	217.897 ミリ秒	1.6415 ミリ秒	1.5355 ミリ秒	1.00
地方最適化	600	99.122 ミリ秒	0.8230 ミリ秒	0.7698 ミリ秒	0.45
部分並列最適化	600	35.825 ミリ秒	0.1003 ミリ秒	0.0837 ミリ秒	0.16
ベースライン	1200	1,763.335 ミリ秒	7.4561 ミリ秒	6.2262 ミリ秒	1.00
环节最適化	1200	766.675 ミリ秒	6.1147 ミリ秒	5.7197 ミリ秒	0.43
部分並列最適化	1200	248.801 ミリ秒	0.6040 ミリ秒	0.5043 ミリ秒	0.14
ベースライン	2400	14,533.335 ミリ秒	63.3518 ミリ秒	52.9016 ミリ秒	1.00
位置最適化	2400	6,507.878 ミリ秒	28.2317 ミリ秒	26.4079 ミリ秒	0.45
部分並列最適化	2400	2,076.403 ミリ秒	20.8320 ミリ秒	19.4863 ミリ秒	0.14
ベースライン	4800	119,768.219 ミリ秒	181.5514 ミリ秒	160.9406 ミリ秒	1.00
个部分最適化	4800	55,590.374 ミリ秒	414.6051 ミリ秒	387.8218 ミリ秒	0.46
部分並列最適化	4800	15,614.506 ミリ秒	115.6996 ミリ秒	102.5647 ミリ秒	0.13

結果から、 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 でベクトル化を実装するには、 と型を选用できます ( 名前空間の型も选用できますが、現在の実装では少し宽度なため、ここでは选用しませんが、自由度に試してみてください)。

 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 番目の非ベクトル化ループが残りの反復を計算します。

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

方式	サイズ	均	エラー	標準较差	比列
ベースライン	300	27.525 ミリ秒	0.1937 ミリ秒	0.1617ミリ秒	1.00
部位最適化	300	12.399ミリ秒	0.0943 ミリ秒	0.0882 ミリ秒	0.45
局部並列最適化	300	6.252 ミリ秒	0.0211ミリ秒	0.0187ミリ秒	0.23
部分ベクトル最適化	300	3.056ミリ秒	0.0301ミリ秒	0.0281 ミリ秒	0.11
ベースライン	600	217.897 ミリ秒	1.6415ミリ秒	1.5355 ミリ秒	1.00
这部分最適化	600	99.122 ミリ秒	0.8230 ミリ秒	0.7698 ミリ秒	0.45
地方並列最適化	600	35.825 ミリ秒	0.1003 ミリ秒	0.0837 ミリ秒	0.16
部分ベクトル最適化	600	24.378 ミリ秒	0.4320 ミリ秒	0.4041 ミリ秒	0.11
ベースライン	1200	1,763.335 ミリ秒	7.4561 ミリ秒	6.2262 ミリ秒	1.00
一部分最適化	1200	766.675 ミリ秒	6.1147 ミリ秒	5.7197 ミリ秒	0.43
这部分並列最適化	1200	248.801 ミリ秒	0.6040 ミリ秒	0.5043 ミリ秒	0.14
部分ベクトル最適化	1200	185.628 ミリ秒	2.1240 ミリ秒	1.9868 ミリ秒	0.11
ベースライン	2400	14,533.335 ミリ秒	63.3518 ミリ秒	52.9016 ミリ秒	1.00
部门最適化	2400	6,507.878 ミリ秒	28.2317 ミリ秒	26.4079 ミリ秒	0.45
部分並列最適化	2400	2,076.403 ミリ秒	20.8320 ミリ秒	19.4863 ミリ秒	0.14
部位ベクトル最適化	2400	2,568.676 ミリ秒	31.7359 ミリ秒	29.6858 ミリ秒	0.18
ベースライン	4800	119,768.219 ミリ秒	181.5514 ミリ秒	160.9406 ミリ秒	1.00
这部分最適化	4800	55,590.374 ミリ秒	414.6051 ミリ秒	387.8218 ミリ秒	0.46
部分並列最適化	4800	15,614.506 ミリ秒	115.6996 ミリ秒	102.5647 ミリ秒	0.13
区域ベクトル最適化	4800	18,257.991 ミリ秒	84.5978 ミリ秒	79.1329 ミリ秒	0.15

結果から明らかなように、ベクトル化により計算時間が大幅に短縮されました。非並列化バージョン ( 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 頂点のグラフの場合) 短縮されました。

形式	サイズ	一般	エラー	標準问题	比例
ベースライン	300	27.525 ミリ秒	0.1937ミリ秒	0.1617ミリ秒	1.00
地方最適化	300	12.399ミリ秒	0.0943 ミリ秒	0.0882 ミリ秒	0.45
区域並列最適化	300	6.252 ミリ秒	0.0211ミリ秒	0.0187ミリ秒	0.23
环节ベクトル最適化	300	3.056ミリ秒	0.0301ミリ秒	0.0281 ミリ秒	0.11
部分並列ベクトル最適化	300	3.008 ミリ秒	0.0075ミリ秒	0.0066ミリ秒	0.11
ベースライン	600	217.897 ミリ秒	1.6415 ミリ秒	1.5355 ミリ秒	1.00
那部分最適化	600	99.122 ミリ秒	0.8230 ミリ秒	0.7698 ミリ秒	0.45
地方並列最適化	600	35.825 ミリ秒	0.1003 ミリ秒	0.0837 ミリ秒	0.16
环节ベクトル最適化	600	24.378 ミリ秒	0.4320 ミリ秒	0.4041 ミリ秒	0.11
部分並列ベクトル最適化	600	13.425ミリ秒	0.0319ミリ秒	0.0283 ミリ秒	0.06
ベースライン	1200	1,763.335 ミリ秒	7.4561 ミリ秒	6.2262 ミリ秒	1.00
一些最適化	1200	766.675 ミリ秒	6.1147 ミリ秒	5.7197 ミリ秒	0.43
部件並列最適化	1200	248.801 ミリ秒	0.6040 ミリ秒	0.5043 ミリ秒	0.14
部件ベクトル最適化	1200	185.628 ミリ秒	2.1240 ミリ秒	1.9868 ミリ秒	0.11
部分並列ベクトル最適化	1200	76.770 ミリ秒	0.3021 ミリ秒	0.2522 ミリ秒	0.04
ベースライン	2400	14,533.335 ミリ秒	63.3518 ミリ秒	52.9016 ミリ秒	1.00
地方最適化	2400	6,507.878 ミリ秒	28.2317 ミリ秒	26.4079 ミリ秒	0.45
个部分並列最適化	2400	2,076.403 ミリ秒	20.8320 ミリ秒	19.4863 ミリ秒	0.14
一些ベクトル最適化	2400	2,568.676 ミリ秒	31.7359 ミリ秒	29.6858 ミリ秒	0.18
部分並列ベクトル最適化	2400	1,281.877 ミリ秒	25.1127 ミリ秒	64.8239 ミリ秒	0.09
ベースライン	4800	119,768.219 ミリ秒	181.5514 ミリ秒	160.9406 ミリ秒	1.00
区域最適化	4800	55,590.374 ミリ秒	414.6051 ミリ秒	387.8218 ミリ秒	0.46
部份並列最適化	4800	15,614.506 ミリ秒	115.6996 ミリ秒	102.5647 ミリ秒	0.13
有些ベクトル最適化	4800	18,257.991 ミリ秒	84.5978 ミリ秒	79.1329 ミリ秒	0.15
部分並列ベクトル最適化	4800	12,785.824 ミリ秒	98.6947 ミリ秒	87.4903 ミリ秒	0.11

結論

この文章发表では、全ペア较长経路問題に少し踏み込んで、それを解決するために C# で Floyd-Warshall アルゴリズムを実装しました。また、データを珍惜し、.NET とハードウェアの低レベル機能を运用するように実装を提升しました。

この网上投稿では、私たちは「オールイン」をしました。しかし、実際のアプリケーションでは、私たちだけではないということを覚えておくことが很重要です。ハードコアな並列化は、システム パフォーマンスを适度に太低させ、メリットよりもデメリットをもたらす几率性があります。双方、ベクトル化は、プラットフォームに依存しない技术で実行すれば、少し健康安全です。一台のベクトル运行命令は、对应の CPU でのみ安全使用できる場合があることに目光してください。

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

参考文献

---