Mengapa transposing matriks 512x512 jauh lebih lambat daripada transposing matriks 513x513?

Setelah melakukan beberapa percobaan pada matriks persegi dengan ukuran yang berbeda, sebuah pola muncul. Selalu, transposing ukuran matriks 2^nlebih lambat daripada transposing ukuran satu2^n+1 . Untuk nilai kecil n, perbedaannya tidak besar.

Namun perbedaan besar terjadi pada nilai 512. (setidaknya untuk saya)

Penafian: Saya tahu fungsi ini tidak benar-benar mengubah matriks karena pertukaran elemen yang ganda, tetapi tidak ada bedanya.

Mengikuti kode:

#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;
}

Mengubah MATSIZEmemungkinkan kita mengubah ukurannya (ya!). Saya memposting dua versi di ideone:

• ukuran 512 - rata-rata 2,46 ms - http://ideone.com/1PV7m
• ukuran 513 - rata-rata 0,75 ms - http://ideone.com/NShpo

Di lingkungan saya (MSVS 2010, optimisasi penuh), perbedaannya serupa:

• ukuran 512 - rata-rata 2,19 ms
• ukuran 513 - rata-rata 0,57 ms

Mengapa ini terjadi?

Penjelasannya berasal dari Agner Fog dalam Mengoptimalkan perangkat lunak dalam C ++ dan mengurangi cara data diakses dan disimpan dalam cache.

Untuk syarat dan info terperinci, lihat entri wiki tentang caching , saya akan mempersempitnya di sini.

Cache diatur dalam set dan baris . Pada suatu waktu, hanya satu set digunakan, dari mana setiap baris yang dikandungnya dapat digunakan. Memori yang bisa dicerminkan oleh garis berapa kali jumlah baris memberi kita ukuran cache.

Untuk alamat memori tertentu, kita dapat menghitung set mana yang harus dicerminkan dengan rumus:

set = ( address / lineSize ) % numberOfsets

Formula semacam ini idealnya memberikan distribusi yang seragam di set, karena setiap alamat memori cenderung dibaca (saya katakan idealnya ).

Jelas bahwa tumpang tindih dapat terjadi. Jika ada cache miss, memori dibaca dalam cache dan nilai lama diganti. Ingat setiap set memiliki sejumlah baris, di mana yang paling baru digunakan ditimpa dengan memori yang baru dibaca.

Saya akan mencoba mengikuti contoh Agner:

Asumsikan setiap set memiliki 4 baris, masing-masing memegang 64 byte. Kami pertama-tama mencoba membaca alamat 0x2710yang sudah diatur 28. Dan kemudian kami juga berusaha untuk alamat membaca 0x2F00, 0x3700,0x3F00 dan 0x4700. Semua ini milik set yang sama. Sebelum membaca 0x4700, semua baris di set akan ditempati. Membaca bahwa memori mengusir garis yang ada di set, garis yang awalnya dipegang 0x2710. Masalahnya terletak pada kenyataan bahwa kita membaca alamat yang (untuk contoh ini) 0x800terpisah. Ini adalah langkah kritis (sekali lagi, untuk contoh ini).

Langkah kritis juga dapat dihitung:

criticalStride = numberOfSets * lineSize

Variabel yang berjarak criticalStrideatau beberapa bagian bersaing untuk baris cache yang sama.

Ini adalah bagian teori. Selanjutnya, penjelasannya (juga Agner, saya mengikutinya dengan cermat untuk menghindari kesalahan):

Asumsikan sebuah matriks 64x64 (ingat, efeknya bervariasi sesuai dengan cache) dengan cache 8kb, 4 baris per set * ukuran garis 64 byte. Setiap baris dapat menampung 8 elemen dalam matriks (64-bit)int ).

Langkah kritis akan menjadi 2048 byte, yang sesuai dengan 4 baris matriks (yang kontinu dalam memori).

Asumsikan kita sedang memproses baris 28. Kami mencoba untuk mengambil elemen dari baris ini dan menukar mereka dengan elemen dari kolom 28. 8 elemen pertama dari baris membuat garis cache, tetapi mereka akan masuk ke 8 berbeda baris cache di kolom 28. Ingat, langkah kritis adalah 4 baris terpisah (4 elemen berturut-turut dalam kolom).

Ketika elemen 16 tercapai di kolom (4 baris cache per set & 4 baris terpisah = masalah) elemen ex-0 akan diusir dari cache. Ketika kami mencapai akhir kolom, semua baris cache sebelumnya akan hilang dan perlu memuat ulang pada akses ke elemen berikutnya (seluruh baris ditimpa).

Memiliki ukuran yang bukan kelipatan dari langkah kritis mengacaukan skenario yang sempurna ini untuk bencana, karena kita tidak lagi berurusan dengan elemen-elemen yang langkah kritis terpisah di vertikal, sehingga jumlah cache ulang sangat berkurang.

Penafian lain - saya baru saja mendapatkan penjelasan dan berharap saya berhasil, tetapi saya mungkin salah. Bagaimanapun, saya menunggu jawaban (atau konfirmasi) dari Mysticial . :)

Luchian memberikan penjelasan tentang alasannya perilaku ini terjadi, tapi saya pikir itu akan menjadi ide bagus untuk menunjukkan satu solusi yang mungkin untuk masalah ini dan pada saat yang sama menunjukkan sedikit tentang algoritma cache cache.

Algoritma Anda pada dasarnya tidak:

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

yang hanya mengerikan untuk CPU modern. Salah satu solusinya adalah mengetahui detail tentang sistem cache Anda dan menyesuaikan algoritme untuk menghindari masalah tersebut. Bekerja dengan baik selama Anda tahu detail-detail itu .. tidak terlalu portabel.

Bisakah kita berbuat lebih baik dari itu? Ya kita bisa: Pendekatan umum untuk masalah ini adalah algoritma yang tidak menyadari cache yang sehingga seperti namanya, jangan bergantung pada ukuran cache tertentu [1]

Solusinya akan terlihat seperti ini:

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];
    }
}

Sedikit lebih rumit, tetapi tes singkat menunjukkan sesuatu yang cukup menarik pada e8400 kuno saya dengan rilis VS2010 x64, testcode for 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

Sunting: Tentang pengaruh ukuran: Ini jauh lebih sedikit diucapkan meskipun masih terlihat sampai batas tertentu, itu karena kami menggunakan solusi berulang sebagai simpul daun alih-alih berulang ke 1 (pengoptimalan biasa untuk algoritma rekursif). Jika kita menetapkan LEAFSIZE = 1, cache tidak mempengaruhi saya [ 8193: 1214.06; 8192: 1171.62ms, 8191: 1351.07ms- itu ada di dalam margin of error, fluktuasi ada di area 100ms; "tolok ukur" ini bukanlah sesuatu yang akan membuat saya terlalu nyaman jika kita menginginkan nilai yang sepenuhnya akurat])

[1] Sumber untuk hal-hal ini: Ya, jika Anda tidak bisa mendapatkan kuliah dari seseorang yang bekerja dengan Leiserson dan rekan mengenai hal ini .. Saya menganggap makalah mereka sebagai titik awal yang baik. Algoritma tersebut masih sangat jarang dijelaskan - CLR memiliki catatan kaki tunggal tentang mereka. Tetap itu cara yang bagus untuk mengejutkan orang.

Sunting (catatan: Saya bukan orang yang memposting jawaban ini; Saya hanya ingin menambahkan ini):
Berikut adalah versi C ++ lengkap dari kode di atas:

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];
            }
        }
    }
}

Sebagai ilustrasi penjelasan dalam jawaban Luchian Grigore , inilah tampilan keberadaan cache matriks untuk dua kasus matriks 64x64 dan 65x65 (lihat tautan di atas untuk perincian tentang angka).

Warna dalam animasi di bawah ini berarti yang berikut:

• - tidak ada dalam cache,
• - dalam cache,
• - hit cache,
• - baru saja membaca dari RAM,
• - miss cache.

Kasus 64x64:

Perhatikan bagaimana hampir setiap akses ke baris baru menghasilkan cache yang hilang. Dan sekarang tampilannya seperti kasus normal, sebuah matriks 65x65:

Di sini Anda dapat melihat bahwa sebagian besar akses setelah pemanasan awal adalah hit cache. Beginilah cara cache CPU dimaksudkan untuk bekerja secara umum.

_{Kode yang menghasilkan bingkai untuk animasi di atas dapat dilihat di sini .}