Matrix Multiplication on CUDA Leander Sturm Daniel Gran Hardware-Software-Co-Design

1. Matrix Multiplication on CUDA Leander Sturm Daniel Gran Hardware-Software-Co-Design Universität Erlangen-Nürnberg leander.sturm@mb.stud.uni-erlangen.de daniel.gran@gmx.de

2. Übersicht • Allgemein • Hintergründe und Eigenschaften • allgemeine Optimierung • CUDA • simpler Ansatz und optimierter Ansatz • praktische Ergebnisse • CuBLAS • Eigenschaften & praktische Ergebnisse • CPU • Optimierungsmöglichkeiten & praktische Ergebnisse • Fazit • pro und contra CPU und GPU • direkter Leistungsvergleich

3. Die Matrix-Multiplikation • Matrizen • Schlüsselkonzept der Linearen Algebra • lineare Gleichungssysteme und Abbildungen • Matrix-Multiplikation • Transformationen, geometrische Rotationen • verketten von Abbildungen

4. Eigenschaften • berechnungsintensiv • drei verschachtelte Schleifen • O(n3) ( n2.7 für Strassen-Algorithmus ) n…..Größe der Matrix • für eine Ergebnis-Zeile Zugriff auf eine ganze Matrix nötig • performante Implementierung nicht trivial • naive Implementierung: unzusammenhängende Speicherzugriffe • Strassen-Algorithmus numerisch instabil • blocking in Größe der lokalen Puffer/Caches nötig • Vektorisierung & Parallelisierung von Teilen der Schleifen nötig

5. allgemeine Optimierung • “blocking” • Optimierung des Speicherzugriffs • Größenordnung • CPU: Matrixgröße ca. 4-60 (für L1 bzw. L2 Cache) • GPU mit CUDA: Matrixgröße 16 (16kB shared memory) • massive Erhöhung der cache-hit-Rate -> hoher Speedup!

6. Implementierung in CUDA (I) N • Simpler Ansatz • Kein Blocking – lediglich ein Block • Matrix M und N liegen im “Global Memory” • Jeder Thread berechnet ein Element von P WIDTH M P • Offensichtlich suboptimal • Keine Nutzung des schnellen “Shared Memory” • Lediglich ein Shader-Cluster aktiv ty WIDTH tx WIDTH WIDTH

7. Implementierung in CUDA (II) • Blocking • Jeder Block berechnet eine Submatric Csub • Jeder Thread berechnet ein Element von Csub • Arbeit aufteilbar auf viele Shader-Cluster • Shared Memory • kleine Daten-Portionen von “Global Memory” in “Shared Memory” kopieren • Jeder Block läd Submatrix Asub und Bsub • Jeder Thread läd ein Element von Asub und Bsub

8. Source-Code (I) //Kopieren der Matrizen von Host auf Device (Global Memory) cudaMemcpy(aDevice.elements, aHost.elements, aSize, cudaMemcpyHostToDevice); cudaMemcpy(bDevice.elements, bHost.elements, bSize, cudaMemcpyHostToDevice); //Dimensionierung von Block und Grid dim3 blockDim(BLOCK_SIZE, BLOCK_SIZE); dim3 gridDim(cHost.height/BLOCK_SIZE,cHost.width/BLOCK_SIZE); //Kernel-Aufruf matrixMulKernel<<<gridDim, blockDim>>>(aDevice, bDevice, cDevice); //Zurückkopieren der Ergebnisse cudaMemcpy(cHost.elements, cDevice.elements, cSize, cudaMemcpyDeviceToHost);

9. Source-Code (II) __shared__ float aSubMatrix[BLOCK_SIZE][BLOCK_SIZE]; __shared__ float bSubMatrix[BLOCK_SIZE][BLOCK_SIZE]; float cValue = 0; //Loop über Blocks der Matrizen A und B for (int i=0; i<(aDevice.width/BLOCK_SIZE); i++){ //Daten in den Shared-Memory laden aSubMatrix[ty][tx] = aDevice.elements[IDX( (ty+(by*BS)), (tx+(i*BS)), aD.width)]; bSubMatrix[ty][tx] = bDevice.elements[IDX( (ty+(i*BS)), (tx+(bx*BS)), bD.width)]; __syncthreads(); //Multiplikation der Elemente for (int k=0; k<BLOCK_SIZE; k++){ cValue += aSubMatrix[ty][k]*bSubMatrix[k][tx]; } __syncthreads(); } cDevice.elements[IDX( (ty+(by*BS)), (tx+(bx*BS)), cD.width)] = cValue;

10. praktische Ergebnisse • starke Abhängigkeit von Größe • Große Matrizen = viele Threads • GPU benötigt möglichst viele Threads für gute Performance • Daten-Transfer ist Flaschenhals bei wenigen Berechnungen

11. praktische Ergebnisse (II)

12. CuBLAS • BLAS • “Basic Linear Algebra Subprograms” • hochoptimierte Algorithmen • verfügbar von allen Herstellern von GPU’s und CPU’s • SGEMM • S… Single precision • GEMM… GEneric Matrix Multiplication • CuBLAS • mittels CUDA optimierte BLAS-Implementierung • im CUDA-SDK enthalten

13. praktische Ergebnisse - CuBLAS

14. Matrix-Multiplikation auf der CPU I • moderne CPU’s • 4 Kerne • pro Kern zwei 128bit Vektor-Rechenwerke -> 32 Gleitkomma-Operationen pro Takt • Speicheranbindung • viel Speicher, aber hohe Latenz • große Caches und Hardware-Prefetch -> besondere Sorgfalt bei Zugriffen nötig -> Zugriff möglichst sequentiell

15. Matrix-Multiplikation auf der CPU II • Parallelisierung für Multi-Core • separate Prozessorkerne • Threading • Prozeßkommunikation • getrennte Speicherbereiche • Inter-Prozeß-Kommunikation • Vektorisierung • 128bit Register: 4 sp Werte • Programmierung mit SSE/SSE2-Befehlen • optimierende Compiler • Assembler

16. Matrix-Multiplikation auf der CPU III • hierarchisches Speichermodell • 2 bis 3 separate Caches • in Größe, Bandbreite und Latenz gestaffelte • Optimierung des Speicherzugriffes • gestaffelte Schleifen: blocking • Schleifengrößen jeweils in Größenordnung der Caches

17. praktische Ergebnisse • Intel Core 2 Duo 3,2 GHz • automatische Optimierung • Intel C++ 10.1 • Parallelisierung • Vektorisierung

18. praktische Ergebnisse II • Übergangspunkt CPU-GPU • oberhalb Matrix-Größe ~256 gewinnt GPU • ~ 30’000’000 Operationen • 65536 Threads auf der GPU • Matrix-Größe <=256 gewinnt CPU • kein Transfer-Overhead • zuwenige Threads für GPU

19. Fazit • jede Implementierung hat Vor- und Nachteile: • Problemgröße • Optimierungs-Overhead • Transfer-Overhead • Optimierungsaufwand auf CPU vergleichbar mit CUDA • Vektorisierung • SSE/SSE2 mit Assembler • Parallelisierung • separate Speicherbereiche • Prozeßkommunikation • kann aber vom Compiler erledigt werden!

20. The End • Wir sagen Danke für’s Zuhören! Fragen?