AVX2 Sparse Matrix Multiplikation

Question

Ich versuche, die neuen AVX2 GATHER Anweisungen zu verwenden, um eine spärliche Matrix - Vektormultiplikation zu beschleunigen. Die Matrix befindet sich im CSR-Format (oder Yale-Format) mit einem Zeilenzeiger, der auf ein Spaltenindexarray zeigt, das wiederum die Spalten enthält. Der C-Code für eine solche Matte-vec mul ist wie folgt aussehen:

for (int row = 0; row < n_rows - 1; row++) { 
    double rowsum = 0; 
    for (int col = row_ptr[row]; col < row_ptr[row + 1]; col++) { 
     rowsum += values[col] * x[col_indices[col]]; 
    } 
    result[row] = rowsum; 
} 

Jetzt ist mein Ziel in diesem mit AVX2 intrinsics zu beschleunigen. Der folgende Code funktioniert mit dem neuesten Intel oder GCC, basierend auf https://blog.fox-toolkit.org/?p=174. Ich habe den Rest hier entfernt, weil meine Zeilen alle auf 4 Doubles (Spalten% 4 == 0) ausgerichtet sind (lucky me). Ich habe den Code auch mit Rest, wenn jemand interessiert ist, aber der Punkt ist, der Code ist tatsächlich etwas langsamer. Ich habe die Disassemblierung überprüft und für die obige Version werden nur FP-Anweisungen generiert und für meinen AVX2-Code erscheinen alle AVX2-Ops wie erwartet. Selbst mit kleinen Matrizen, die in den Cache passen, ist die AVX2-Version nicht gut. Ich bin verwirrt hier ...

double* value_base = &values[0]; 
double* x_base = &x[0]; 
int* index_base = &col_indices[0]; 


for (int row = 0; row < n_rows - 1; row++) { 
    int col_length = row_ptr[row + 1] - row_ptr[row]; 

    __m256d rowsum = _mm256_set1_pd(0.); 
    for (int col4 = 0; col4 < col_length; col4 += 4) { 
     // Load indices for x vector(const __m128i*) 
     __m128i idxreg  = _mm_load_si128((const __m128i*)index_base); 
     // Load 4 doubles from x indexed by idxreg (AVX2) 
     __m256d x_  = _mm256_i32gather_pd(x_base, idxreg, 8); 
     // Load 4 doubles linear from memory (value array) 
     __m256d v_  = _mm256_load_pd(value_base); 
     // FMA: rowsum += x_ * v_ 
     rowsum = _mm256_fmadd_pd(x_, v_, rowsum); 

     index_base += 4; 
     value_base += 4; 
    } 
    __m256d s = _mm256_hadd_pd(rowsum, rowsum); 
    result[row] = ((double*)&s)[0] + ((double*)&s)[2]; 
    // Alternative (not faster): 
    // Now we split the upper and lower AVX register, and do a number of horizontal adds 
    //__m256d hsum = _mm256_add_pd(rowsum, _mm256_permute2f128_pd(rowsum, rowsum, 0x1)); 
    //_mm_store_sd(&result[row], _mm_hadd_pd(_mm256_castpd256_pd128(hsum), _mm256_castpd256_pd128(hsum))); 
} 

Alle Vorschläge sind willkommen.

Vielen Dank, Chris

Answer 1

auf Haswell Sammeln Sie ist langsam. Ich implementierte eine 8-Bit-Index-LUT-Suche von 16-Bit-Werten (für GF16 multipliziere für Par2) auf verschiedene Arten, um herauszufinden, was am schnellsten ist. Auf Haswell nahm die VPGATHERDD Version 1,7x so lange wie die movd/pinsrw Version. (Nur ein paar VPUNPCK/Schicht Anweisungen wurden über die Raffung hinaus benötigt.) code here, if anyone wants to run the benchmark.

Wie üblich, wenn eine Anweisung zuerst eingeführt wird, widmen sie nicht eine riesige Menge an Silizium, um es superschnell zu machen. Es ist nur da, um HW-Support da draußen zu bekommen, also kann Code geschrieben werden, um ihn zu benutzen. Für optimale Leistung auf allen CPUs, dann müssen Sie tun, was x264 für pshufb getan hat: ein SLOW_SHUFFLE Flag für CPUs wie Core2, und Faktor, dass in Ihre beste Routine-Suche Funktion-Zeiger-Einstellung, anstatt nur, was eine CPU unterstützt.

Für Projekte, die weniger fanatisch sind, asm-Versionen für jede CPU abgestimmt zu haben, auf der sie laufen können, führt die Einführung einer nicht beschleunigten Version einer Anweisung dazu, dass die Leute früher damit arbeiten, wenn das nächste Design kommt und schneller, mehr Code wird schneller. Ein Design wie Haswell zu veröffentlichen, bei dem es sich um eine Verlangsamung handelt, ist ein wenig heikel. Vielleicht wollten sie sehen, wie die Leute es benutzen würden? Es erhöht die Codedichte, was hilft, wenn die Sammlung nicht in einer engen Schleife ist.

Broadwell soll eine schnellere Gather-Implementierung haben, aber ich habe keinen Zugriff darauf. Das Intel-Handbuch, das die Latenz/den Durchsatz für Anweisungen auflistet, sagt, dass Broadwells Sammlung etwa 1,6 mal schneller ist, also immer noch etwas langsamer als eine handgefertigte Schleife, die die Indizes in GP-Registern verschiebt/entpackt und sie für PINSRW in Vektoren verwendet.

Wenn gather Fälle ausnutzen könnte, in denen mehrere Elemente denselben Index oder sogar einen Index haben, der auf denselben 32B-Abrufblock verweist, kann es abhängig von den Eingabedaten zu großen Beschleunigungen kommen.

Hoffentlich wird sich Skylake noch weiter verbessern. Ich dachte, ich hätte etwas gelesen, dass es sagen würde, aber als ich nachgesehen habe, habe ich nichts gefunden.

RE: dünn besetzte Matrizen: Gibt es kein Format, das die Daten dupliziert, so dass zusammenhängende Lesevorgänge für Zeilen oder Spalten möglich sind? Es ist nicht etwas, für das ich Code schreiben musste, aber ich denke, ich habe es in einigen Antworten erwähnt.