241

さまざまなサイズの正方行列でいくつかの実験を行った後、パターンが浮かび上がりました。常に、size の行列の転置は、 size の行列の2^n転置よりも遅くなります2^n+1。の値が小さい場合n、違いは大きくありません。

ただし、値が 512 を超えると大きな違いが生じます。(少なくとも私にとっては)

免責事項:要素の二重スワップのために、関数が実際に行列を転置しないことは知っていますが、違いはありません。

コードに従います:

#define SAMPLES 1000
#define MATSIZE 512

#include <time.h>
#include <iostream>
int mat[MATSIZE][MATSIZE];

void transpose()
{
   for ( int i = 0 ; i < MATSIZE ; i++ )
   for ( int j = 0 ; j < MATSIZE ; j++ )
   {
       int aux = mat[i][j];
       mat[i][j] = mat[j][i];
       mat[j][i] = aux;
   }
}

int main()
{
   //initialize matrix
   for ( int i = 0 ; i < MATSIZE ; i++ )
   for ( int j = 0 ; j < MATSIZE ; j++ )
       mat[i][j] = i+j;

   int t = clock();
   for ( int i = 0 ; i < SAMPLES ; i++ )
       transpose();
   int elapsed = clock() - t;

   std::cout << "Average for a matrix of " << MATSIZE << ": " << elapsed / SAMPLES;
}

変更MATSIZEすると、サイズを変更できます (当たり前!)。ideone に 2 つのバージョンを投稿しました。

私の環境(MSVS 2010、完全最適化)では、違いは似ています:

  • サイズ 512 - 平均2.19 ミリ秒
  • サイズ 513 - 平均0.57 ミリ秒

なぜこうなった?

4

3 に答える 3

208

説明は、C++ でのソフトウェアの最適化の Agner Fog から来ており、データへのアクセス方法とキャッシュへの保存方法に還元されます。

用語と詳細な情報については、キャッシングに関するウィキ エントリを参照してください。ここで絞り込みます。

キャッシュは、セットラインで編成されます。一度に 1 つのセットのみが使用され、そのセットに含まれるすべての行を使用できます。ラインがミラーリングできるメモリの倍のライン数がキャッシュサイズになります。

特定のメモリ アドレスについて、次の式を使用して、どのセットがそれをミラーリングする必要があるかを計算できます。

set = ( address / lineSize ) % numberOfsets

この種の式は、各メモリアドレスが読み取られる可能性が高いため、セット全体に均一な分布を与えるのが理想的です (理想的には と言いました)。

オーバーラップが発生する可能性があることは明らかです。キャッシュ ミスの場合、メモリがキャッシュに読み込まれ、古い値が置き換えられます。各セットには多数の行があり、そのうち最も使用頻度の低い行が新しく読み取られたメモリで上書きされることに注意してください。

アグナーの例にいくらか従おうとします:

各セットに 4 行があり、それぞれが 64 バイトを保持しているとします。0x2710最初に set に入るaddress を読み取ろうとします28。そして、アドレス0x2F000x37000x3F00およびも読み取ろうとします0x4700。これらはすべて同じセットに属します。を読み取る前0x4700は、セット内のすべての行が占有されていたはずです。そのメモリを読み取ると、セット内の既存の行、つまり最初に保持されていた行が追い出され0x2710ます。問題は、(この例では)0x800離れたアドレスを読み取るという事実にあります。これが重要なストライドです (この例でも)。

クリティカルストライドも計算できます。

criticalStride = numberOfSets * lineSize

間隔を空けて配置された、または複数離れた変数criticalStrideが、同じキャッシュ ラインを求めて競合します。

これが理論の部分です。次に、説明です (Agner も、間違いを避けるために注意深くフォローしています)。

64x64 のマトリックスを想定します (効果はキャッシュによって異なることに注意してください)。キャッシュは 8kb で、セットあたり 4 ライン * ライン サイズ 64 バイトです。各行は、行列 (64 ビットint) の 8 つの要素を保持できます。

重要なストライドは 2048 バイトで、これは行列の 4 行 (メモリ内で連続している) に対応します。

行 28 を処理していると仮定します。この行の要素を取得し、列 28 の要素と交換しようとしています。行の最初の 8 つの要素はキャッシュ ラインを構成しますが、8 つの異なる行になります。列 28 のキャッシュ行。重要なストライドは 4 行離れていることを思い出してください (列内の 4 つの連続する要素)。

列内の要素 16 に到達すると (セットごとに 4 つのキャッシュ ラインと 4 つの行が離れている = 問題)、ex-0 要素がキャッシュから削除されます。列の最後に到達すると、以前のすべてのキャッシュ ラインが失われ、次の要素へのアクセス時にリロードが必要になります (ライン全体が上書きされます)。

クリティカル ストライドの倍数ではないサイズを使用すると、この完璧な災害シナリオが台無しになります。垂直方向にクリティカル ストライド離れた要素を扱う必要がなくなるため、キャッシュのリロード回数が大幅に減少します。

もう 1 つの免責事項- 説明に頭を悩ませていて、それを釘付けにしたいと思っていますが、間違っている可能性があります。とにかく、 Mystialからの応答 (または確認) を待っています。:)

于 2012-07-10T13:00:21.473 に答える
80

Luchian はなぜこの動作が発生するのかを説明していますが、この問題の解決策の 1 つを示し、同時にキャッシュ忘却アルゴリズムについて少し説明するのは良い考えだと思いました。

あなたのアルゴリズムは基本的に次のことを行います。

for (int i = 0; i < N; i++) 
   for (int j = 0; j < N; j++) 
        A[j][i] = A[i][j];

これは、最新の CPU にとっては恐ろしいことです。解決策の 1 つは、キャッシュ システムの詳細を把握し、アルゴリズムを微調整してこれらの問題を回避することです。これらの詳細を知っている限り、うまく機能します..特にポータブルではありません。

それよりもうまくできるでしょうか?はい、できます: この問題に対する一般的なアプローチは、名前が示すように、特定のキャッシュ サイズに依存することを回避するキャッシュ無視アルゴリズムです [1]。

解決策は次のようになります。

void recursiveTranspose(int i0, int i1, int j0, int j1) {
    int di = i1 - i0, dj = j1 - j0;
    const int LEAFSIZE = 32; // well ok caching still affects this one here
    if (di >= dj && di > LEAFSIZE) {
        int im = (i0 + i1) / 2;
        recursiveTranspose(i0, im, j0, j1);
        recursiveTranspose(im, i1, j0, j1);
    } else if (dj > LEAFSIZE) {
        int jm = (j0 + j1) / 2;
        recursiveTranspose(i0, i1, j0, jm);
        recursiveTranspose(i0, i1, jm, j1);
    } else {
    for (int i = i0; i < i1; i++ )
        for (int j = j0; j < j1; j++ )
            mat[j][i] = mat[i][j];
    }
}

少し複雑ですが、VS2010 x64 リリースの古い e8400 で非常に興味深いことが短いテストで示されています。MATSIZE 8192

int main() {
    LARGE_INTEGER start, end, freq;
    QueryPerformanceFrequency(&freq);
    QueryPerformanceCounter(&start);
    recursiveTranspose(0, MATSIZE, 0, MATSIZE);
    QueryPerformanceCounter(&end);
    printf("recursive: %.2fms\n", (end.QuadPart - start.QuadPart) / (double(freq.QuadPart) / 1000));

    QueryPerformanceCounter(&start);
    transpose();
    QueryPerformanceCounter(&end);
    printf("iterative: %.2fms\n", (end.QuadPart - start.QuadPart) / (double(freq.QuadPart) / 1000));
    return 0;
}

results: 
recursive: 480.58ms
iterative: 3678.46ms

編集: サイズの影響について: 1 まで再帰する (再帰アルゴリズムの通常の最適化) の代わりに、反復ソリューションをリーフ ノードとして使用しているため、ある程度目立ちますが、それほど顕著ではありません。LEAFSIZE = 1 に設定すると、キャッシュの影響はありません [ 8193: 1214.06; 8192: 1171.62ms, 8191: 1351.07ms- それは誤差の範囲内です。変動は 100ms 領域にあります。この「ベンチマーク」は、完全に正確な値が必要な場合は、あまり快適ではありません])

[1] この資料の情報源: Leiserson や共同研究者と一緒にこの問題について研究した人から講義を受けられない場合は..彼らの論文が良い出発点になると思います。これらのアルゴリズムについては、まだほとんど説明されていません。CLR には、それらに関する脚注が 1 つあります。それでも、人々を驚かせる素晴らしい方法です。


編集(注:私はこの回答を投稿した人ではありません。これを追加したかっただけです):
上記のコードの完全なC++バージョンは次のとおりです。

template<class InIt, class OutIt>
void transpose(InIt const input, OutIt const output,
    size_t const rows, size_t const columns,
    size_t const r1 = 0, size_t const c1 = 0,
    size_t r2 = ~(size_t) 0, size_t c2 = ~(size_t) 0,
    size_t const leaf = 0x20)
{
    if (!~c2) { c2 = columns - c1; }
    if (!~r2) { r2 = rows - r1; }
    size_t const di = r2 - r1, dj = c2 - c1;
    if (di >= dj && di > leaf)
    {
        transpose(input, output, rows, columns, r1, c1, (r1 + r2) / 2, c2);
        transpose(input, output, rows, columns, (r1 + r2) / 2, c1, r2, c2);
    }
    else if (dj > leaf)
    {
        transpose(input, output, rows, columns, r1, c1, r2, (c1 + c2) / 2);
        transpose(input, output, rows, columns, r1, (c1 + c2) / 2, r2, c2);
    }
    else
    {
        for (ptrdiff_t i1 = (ptrdiff_t) r1, i2 = (ptrdiff_t) (i1 * columns);
            i1 < (ptrdiff_t) r2; ++i1, i2 += (ptrdiff_t) columns)
        {
            for (ptrdiff_t j1 = (ptrdiff_t) c1, j2 = (ptrdiff_t) (j1 * rows);
                j1 < (ptrdiff_t) c2; ++j1, j2 += (ptrdiff_t) rows)
            {
                output[j2 + i1] = input[i2 + j1];
            }
        }
    }
}
于 2012-07-10T13:26:41.610 に答える