Pertanyaan Bagaimana BLAS mendapatkan kinerja ekstrem seperti itu?


Karena penasaran saya memutuskan untuk mengukur fungsi perkalian matriks saya sendiri versus implementasi BLAS ... Saya harus mengatakan paling tidak terkejut pada hasilnya:

Implementasi Kustom, 10 uji coba   Perkalian matriks 1000x1000:

Took: 15.76542 seconds.

Implementasi BLAS, 10 percobaan dari   Perkalian matriks 1000x1000:

Took: 1.32432 seconds.

Ini menggunakan angka floating point presisi tunggal.

Implementasi Saya:

template<class ValT>
void mmult(const ValT* A, int ADim1, int ADim2, const ValT* B, int BDim1, int BDim2, ValT* C)
{
    if ( ADim2!=BDim1 )
        throw std::runtime_error("Error sizes off");

    memset((void*)C,0,sizeof(ValT)*ADim1*BDim2);
    int cc2,cc1,cr1;
    for ( cc2=0 ; cc2<BDim2 ; ++cc2 )
        for ( cc1=0 ; cc1<ADim2 ; ++cc1 )
            for ( cr1=0 ; cr1<ADim1 ; ++cr1 )
                C[cc2*ADim2+cr1] += A[cc1*ADim1+cr1]*B[cc2*BDim1+cc1];
}

Saya punya dua pertanyaan:

  1. Mengingat bahwa perkalian matriks-matriks mengatakan: nxm * mxn membutuhkan perkalian n * n * m, sehingga dalam kasus di atas operasi 1000 ^ 3 atau 1e9. Bagaimana mungkin pada prosesor 2.6Ghz saya untuk BLAS melakukan operasi 10 * 1e9 dalam 1,32 detik? Bahkan jika multiplcations adalah operasi tunggal dan tidak ada lagi yang dilakukan, perlu waktu ~ 4 detik.
  2. Mengapa penerapan saya jauh lebih lambat?

75
2017-08-19 23:30


asal


Jawaban:


Untuk banyak alasan.

Pertama, kompiler Fortran sangat dioptimalkan, dan bahasa memungkinkan mereka untuk menjadi seperti itu. C dan C ++ sangat longgar dalam hal penanganan susunan (misalnya kasus pointer mengacu pada area memori yang sama). Ini berarti bahwa compiler tidak dapat mengetahui sebelumnya apa yang harus dilakukan, dan dipaksa untuk membuat kode generik. Di Fortran, kasus Anda lebih ramping, dan kompilator memiliki kontrol yang lebih baik terhadap apa yang terjadi, memungkinkan dia untuk mengoptimalkan lebih banyak (misalnya menggunakan register).

Hal lain adalah bahwa Fortran menyimpan barang-barang secara bertingkat, sementara C menyimpan data secara berurutan. Saya belum memeriksa kode Anda, tetapi berhati-hatilah dengan cara Anda melakukan produk. Dalam C Anda harus memindai baris bijaksana: cara ini Anda memindai array Anda bersama memori yang berdekatan, mengurangi cache misses. Cache miss adalah sumber inefisiensi pertama.

Ketiga, itu tergantung implementasi blas yang Anda gunakan. Beberapa implementasi mungkin ditulis dalam assembler, dan dioptimalkan untuk prosesor spesifik yang Anda gunakan. Versi netlib ditulis dalam fortran 77.

Juga, Anda melakukan banyak operasi, sebagian besar diulang dan berlebihan. Semua perkalian untuk mendapatkan indeks merugikan untuk kinerja. Saya tidak benar-benar tahu bagaimana hal ini dilakukan dalam BLAS, tetapi ada banyak trik untuk mencegah operasi yang mahal.

Misalnya, Anda dapat mengolah ulang kode Anda dengan cara ini

template<class ValT>
void mmult(const ValT* A, int ADim1, int ADim2, const ValT* B, int BDim1, int BDim2, ValT* C)
{
if ( ADim2!=BDim1 ) throw std::runtime_error("Error sizes off");

memset((void*)C,0,sizeof(ValT)*ADim1*BDim2);
int cc2,cc1,cr1, a1,a2,a3;
for ( cc2=0 ; cc2<BDim2 ; ++cc2 ) {
    a1 = cc2*ADim2;
    a3 = cc2*BDim1
    for ( cc1=0 ; cc1<ADim2 ; ++cc1 ) {
          a2=cc1*ADim1;
          ValT b = B[a3+cc1];
          for ( cr1=0 ; cr1<ADim1 ; ++cr1 ) {
                    C[a1+cr1] += A[a2+cr1]*b;
           }
     }
  }
} 

Cobalah, saya yakin Anda akan menghemat sesuatu.

Pada Anda # 1 pertanyaan, alasannya adalah bahwa skala perkalian matriks sebagai O (n ^ 3) jika Anda menggunakan algoritma trivial. Ada algoritma itu skala jauh lebih baik.


-24
2017-08-19 23:36



Titik awal yang bagus adalah buku yang bagus The Science of Programming Matrix Computations oleh Robert A. van de Geijn dan Enrique S. Quintana-Ortí. Mereka menyediakan versi unduh gratis.

BLAS dibagi menjadi tiga tingkatan:

  • Level 1 mendefinisikan satu set fungsi aljabar linear yang beroperasi pada vektor saja. Fungsi-fungsi ini mendapat manfaat dari vektorisasi (misalnya dari menggunakan SSE).

  • Fungsi Level 2 adalah operasi matriks-vektor, mis. beberapa produk matriks-vektor. Fungsi-fungsi ini dapat diimplementasikan dalam hal fungsi Level1. Namun, Anda dapat meningkatkan kinerja fungsi ini jika Anda dapat memberikan implementasi khusus yang memanfaatkan beberapa arsitektur multiprosesor dengan memori bersama.

  • Fungsi level 3 adalah operasi seperti produk matriks-matriks. Sekali lagi Anda bisa menerapkannya dalam hal fungsi Level2. Tetapi fungsi Level3 melakukan O (N ^ 3) operasi pada O (N ^ 2) data. Jadi, jika platform Anda memiliki hierarki cache, Anda dapat meningkatkan kinerja jika Anda memberikan implementasi khusus cache dioptimalkan / cache ramah. Ini dijelaskan dengan baik dalam buku ini. Dorongan utama dari fungsi Level3 berasal dari optimasi cache. Peningkatan ini secara signifikan melebihi dorongan kedua dari paralelisme dan pengoptimalan perangkat keras lainnya.

By the way, sebagian besar (atau bahkan semua) implementasi BLAS kinerja tinggi TIDAK diimplementasikan di Fortran. ATLAS diimplementasikan dalam C. GotoBLAS / OpenBLAS diimplementasikan dalam C dan bagian-bagian penting kinerjanya di Assembler. Hanya implementasi referensi dari BLAS yang diterapkan di Fortran. Namun, semua implementasi BLAS ini menyediakan antarmuka Fortran sehingga dapat dihubungkan dengan LAPACK (LAPACK mendapatkan semua kinerjanya dari BLAS).

Kompilator yang dioptimalkan memainkan peran kecil dalam hal ini (dan untuk GotoBLAS / OpenBLAS kompiler tidak masalah sama sekali).

Implementasi IMHO no BLAS menggunakan algoritma seperti algoritma Coppersmith – Winograd atau algoritma Strassen. Saya tidak tahu pasti alasannya, tapi ini tebakan saya:

  • Mungkin itu tidak mungkin untuk menyediakan cache dioptimalkan implementasi dari algoritma ini (yaitu Anda akan kehilangan lebih banyak maka Anda akan menang)
  • Algoritma ini secara numerik tidak stabil. Karena BLAS adalah kernel komputasional LAPACK, ini tidak ada jalannya.

Edit / Perbarui:

Makalah baru dan ground breaking untuk topik ini adalah Makalah BLIS. Mereka ditulis dengan sangat baik. Untuk kuliah saya "Dasar-Dasar Software untuk Komputasi Kinerja Tinggi" Saya mengimplementasikan produk matriks-matriks mengikuti kertas mereka. Sebenarnya saya menerapkan beberapa varian dari produk matriks-matriks. Varian yang paling sederhana sepenuhnya ditulis dalam huruf C dan memiliki kurang dari 450 baris kode. Semua varian lain hanya mengoptimalkan loop

    for (l=0; l<MR*NR; ++l) {
        AB[l] = 0;
    }
    for (l=0; l<kc; ++l) {
        for (j=0; j<NR; ++j) {
            for (i=0; i<MR; ++i) {
                AB[i+j*MR] += A[i]*B[j];
            }
        }
        A += MR;
        B += NR;
    }

Kinerja keseluruhan dari produk matriks-matriks hanya tergantung pada loop ini. Sekitar 99,9% waktu dihabiskan di sini. Dalam varian lain saya menggunakan intrinsik dan kode assembler untuk meningkatkan kinerja. Anda dapat melihat tutorial melalui semua varian di sini:

ulmBLAS: Tutorial tentang GEMM (Matrix-Matrix Product)

Bersama dengan makalah BLIS itu menjadi cukup mudah untuk memahami bagaimana perpustakaan seperti Intel MKL dapat memperoleh kinerja seperti itu. Dan mengapa tidak masalah apakah Anda menggunakan penyimpanan utama baris atau kolom!

Tolok ukur terakhir ada di sini (kami menyebut proyek kami ulmBLAS):

Tolok ukur untuk ulmBLAS, BLIS, MKL, openBLAS dan Eigen

Pengeditan / Pembaruan Lainnya:

Saya juga menulis beberapa tutorial tentang bagaimana BLAS digunakan untuk masalah aljabar linear numerik seperti memecahkan sistem persamaan linear:

Faktor LU Kinerja Tinggi

(Faktorisasi LU ini misalnya digunakan oleh Matlab untuk memecahkan sistem persamaan linear.)

Saya berharap menemukan waktu untuk memperluas tutorial untuk menggambarkan dan menunjukkan bagaimana mewujudkan implementasi paralel yang sangat scalable dari faktorisasi LU seperti dalam PLASMA.

Ok, ini dia: Pengkodean Cache Optimized Parallel LU Factorization

P.S .: Saya juga melakukan beberapa eksperimen untuk meningkatkan kinerja uBLAS. Sebenarnya sangat mudah untuk meningkatkan (ya, mainkan kata-kata :)) kinerja uBLAS:

Eksperimen tentang uBLAS.

Di sini proyek serupa dengan API:

Eksperimen tentang BLAZE.


98
2017-07-10 20:23



Jadi pertama-tama BLAS hanyalah antarmuka dari sekitar 50 fungsi. Ada banyak implementasi antarmuka yang bersaing.

Pertama saya akan menyebutkan hal-hal yang sebagian besar tidak berhubungan:

  • Fortran vs C, tidak ada bedanya
  • Algoritma matriks canggih seperti Strassen, implementasi tidak menggunakannya karena mereka tidak membantu dalam praktek

Sebagian besar implementasi memecah setiap operasi menjadi matriks kecil atau operasi vektor dengan cara yang lebih atau kurang jelas. Sebagai contoh, perkalian matriks 1000x1000 yang besar dapat dipecah menjadi urutan perkalian matriks 50x50.

Operasi ukuran kecil berukuran tetap ini (disebut kernel) dikodekan ulang dalam kode assembly khusus CPU menggunakan beberapa fitur CPU dari targetnya:

  • Instruksi gaya SIMD
  • Paralelisme Tingkat Instruksi
  • Kesadaran cache

Lebih lanjut, kernel ini dapat dijalankan secara paralel dengan memperhatikan satu sama lain menggunakan beberapa thread (core CPU), dalam pola desain peta-mengurangi khas.

Lihatlah ATLAS yang merupakan implementasi open source BLAS yang paling sering digunakan. Ini memiliki banyak kernel yang bersaing, dan selama proses pembuatan perpustakaan ATLAS, ia menjalankan kompetisi di antara mereka (beberapa bahkan diparameterisasi, sehingga kernel yang sama dapat memiliki pengaturan yang berbeda). Ini trys konfigurasi yang berbeda dan kemudian memilih yang terbaik untuk sistem target tertentu.

(Tip: Itulah mengapa jika Anda menggunakan ATLAS, Anda lebih baik membangun dan menyetel pustaka dengan tangan untuk mesin khusus Anda, kemudian menggunakan prebuilt.)


19
2018-06-26 14:10



Pertama, ada algoritma yang lebih efisien untuk perkalian matriks daripada yang Anda gunakan.

Kedua, CPU Anda dapat melakukan lebih dari satu instruksi dalam satu waktu.

CPU Anda mengeksekusi 3-4 instruksi per siklus, dan jika unit SIMD digunakan, setiap instruksi memproses 4 pelampung atau 2 kali lipat. (tentu saja angka ini tidak akurat juga, karena CPU biasanya hanya dapat memproses satu instruksi SIMD per siklus)

Ketiga, kode Anda jauh dari optimal:

  • Anda menggunakan pointer mentah, yang berarti bahwa kompiler harus menganggap mereka mungkin alias. Ada kata kunci atau flag khusus kompilator yang dapat Anda tentukan untuk memberi tahu compiler bahwa mereka tidak alias. Sebagai alternatif, Anda harus menggunakan tipe lain selain raw pointer, yang menangani masalah.
  • Anda meronta-ronta cache dengan melakukan traversal naif dari setiap baris / kolom dari input matriks. Anda dapat menggunakan pemblokiran untuk melakukan sebanyak mungkin pekerjaan pada blok matriks yang lebih kecil, yang cocok dalam cache CPU, sebelum pindah ke blok berikutnya.
  • Untuk tugas-tugas numerik murni, Fortran cukup banyak tak terkalahkan, dan C ++ membutuhkan banyak bujukan untuk mendapatkan kecepatan yang sama. Ini bisa dilakukan, dan ada beberapa perpustakaan yang menunjukkannya (biasanya menggunakan templat ekspresi), tapi itu tidak sepele, dan tidak hanya terjadi.

12
2017-08-20 12:12



Saya tidak tahu secara spesifik tentang implementasi BLAS tetapi ada alogorithms yang lebih efisien untuk Multiplikasi Matriks yang memiliki kompleksitas lebih dari O (n3). Yang juga tahu itu Algoritma Strassen 


9
2017-08-20 00:15



Kebanyakan argumen untuk pertanyaan kedua - assembler, membelah menjadi blok, dll. (Tetapi bukan algoritma kurang dari N ^ 3, mereka benar-benar berkembang) - memainkan peran. Tapi kecepatan rendah dari algoritma Anda pada dasarnya disebabkan oleh ukuran matriks dan pengaturan yang tidak menguntungkan dari tiga loop bersarang. Matriks Anda begitu besar sehingga tidak cocok sekaligus dalam memori cache. Anda dapat mengatur ulang loop sedemikian rupa sehingga sebanyak mungkin akan dilakukan pada baris dalam cache, cara ini secara dramatis mengurangi penyegaran cache (pemisahan BTW ke dalam blok kecil memiliki efek analog, paling baik jika loop di atas blok disusun sama). Implementasi model untuk matriks persegi berikut. Di komputer saya konsumsi waktunya sekitar 1:10 dibandingkan dengan penerapan standar (seperti milik Anda). Dengan kata lain: jangan pernah memprogram perkalian matriks di sepanjang skema "baris kali kolom" yang kita pelajari di sekolah. Setelah mengatur ulang loop, perbaikan lebih banyak diperoleh dengan membuka gulungan, kode assembler, dll.

    void vector(int m, double ** a, double ** b, double ** c) {
      int i, j, k;
      for (i=0; i<m; i++) {
        double * ci = c[i];
        for (k=0; k<m; k++) ci[k] = 0.;
        for (j=0; j<m; j++) {
          double aij = a[i][j];
          double * bj = b[j];
          for (k=0; k<m; k++)  ci[k] += aij*bj[k];
        }
      }
    }

Satu lagi komentar: Implementasi ini bahkan lebih baik di komputer saya daripada mengganti semua oleh BLAS rutin cblas_dgemm (coba di komputer Anda!). Tetapi jauh lebih cepat (1: 4) memanggil dgemm_ perpustakaan Fortran secara langsung. Saya pikir rutinitas ini sebenarnya bukan Fortran tetapi kode assembler (saya tidak tahu apa yang ada di perpustakaan, saya tidak punya sumbernya). Benar-benar tidak jelas bagi saya adalah mengapa cblas_dgemm tidak secepat karena sepengetahuan saya itu hanya pembungkus untuk dgemm_.


4
2017-11-30 20:11



Ini adalah kecepatan yang realistis. Untuk contoh apa yang bisa dilakukan dengan SIMD assembler melalui kode C ++, lihat beberapa contoh Fungsi matriks iPhone - ini lebih dari 8x lebih cepat daripada versi C, dan bahkan tidak "dioptimalkan" perakitan - belum ada pipa-lapisan dan ada operasi stack yang tidak perlu.

Juga kode Anda tidak "batasi benar"- bagaimana kompilator tahu bahwa ketika memodifikasi C, itu tidak memodifikasi A dan B?


3
2017-08-20 00:10



Sehubungan dengan kode asli dalam MM multiply, referensi memori untuk sebagian besar operasi adalah penyebab utama kinerja yang buruk. Memori berjalan pada 100-1000 kali lebih lambat dari cache.

Sebagian besar mempercepat berasal dari menggunakan teknik pengoptimalan loop untuk fungsi loop tiga ini dalam MM kalikan. Dua teknik pengoptimalan loop utama digunakan; membuka gulungan dan memblokir. Sehubungan dengan membuka gulungan, kami membuka gulungan dua loop paling luar dan memblokirnya untuk digunakan kembali dalam cache. Outer loop unrolling membantu mengoptimalkan akses data secara temporal dengan mengurangi jumlah referensi memori ke data yang sama pada waktu yang berbeda selama keseluruhan operasi. Memblokir indeks loop pada nomor tertentu, membantu mempertahankan data dalam cache. Anda dapat memilih untuk mengoptimalkan cache L2 atau cache L3.

https://en.wikipedia.org/wiki/Loop_nest_optimization


1
2018-05-02 12:07