moltiplicazione matrice sparse AVX2

Question

Sto cercando di sfruttare le nuove istruzioni AVX2 GATHER per accelerare una moltiplicazione di matrice sparsa - vettore. La matrice è in formato CSR (o Yale) con un puntatore a riga che punta a un array di indici di colonne che a sua volta contiene le colonne. Il codice C per un tale mul mat-vec fa apparire come questo:

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

Ora il mio obiettivo è quello di accelerare questo con intrinseche AVX2. Il seguente codice funziona con Intel o GCC più recenti, basato su https://blog.fox-toolkit.org/?p=174. Ho rimosso il resto qui perché le mie righe sono tutte allineate su 4 doppi (colonne% 4 == 0) comunque (fortunato a me). Ho anche il codice che si occupa del resto se qualcuno è interessato, ma il punto è che il codice è in realtà leggermente più lento. Ho controllato lo smontaggio e per la versione precedente vengono generate solo le istruzioni FP e per il mio codice AVX2 tutte le operazioni AVX2 appaiono come previsto. Anche con piccole matrici che si adattano alla cache, la versione AVX2 non va bene. Sono perplesso qui ...

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

Qualsiasi suggerimento è benvenuto.

Grazie mille, Chris

Answer 1

Raccogliere su Haswell è lento. Ho implementato una ricerca LUT a 8 bit indice di valori a 16 bit (per GF16 moltiplicare per par2) in alcuni modi diversi, per scoprire cosa è più veloce. Su Haswell, la versione VPGATHERDD impiegava 1,7 volte la versione movd/pinsrw. (Solo una coppia VPUNPCK/istruzioni per il turno erano necessarie al di là del raddrizzamento.) code here, if anyone wants to run the benchmark.

Come è comune quando un'istruzione viene introdotta per la prima volta, non dedica un'enorme quantità di silicio per renderlo super-veloce. È lì solo per ottenere supporto HW là fuori, quindi il codice può essere scritto per usarlo. Per prestazioni ottimali su tutte le CPU, è necessario fare ciò che x264 ha fatto per pshufb: avere un flag SLOW_SHUFFLE per CPU come Core2 e farlo rientrare nel proprio set-pointer-setting di ricerca delle routine, piuttosto che solo ciò che insns una CPU sostiene.

Per i progetti meno entusiasti di avere le versioni di ASM ottimizzate per ogni CPU su cui possono essere eseguite, introdurre una versione di accelerazione di un'istruzione farà sì che le persone lo usino prima, quindi quando il prossimo progetto arriva e il suo codice più veloce velocizza. Rilasciare un design come Haswell dove si raduna in realtà è un rallentamento è un po 'rischioso. Forse volevano vedere come la gente avrebbe usato? Aumenta la densità del codice, il che aiuta quando il gather non è in un circuito chiuso.

Si suppone che Broadwell abbia un'implementazione più rapida, ma non ne ho accesso. Il manuale Intel che elenca latenza/throughput per le istruzioni dice che la raccolta di Broadwell è di circa 1,6 volte più veloce, quindi sarebbe ancora leggermente più lento di un loop realizzato a mano che sposta/decomprime gli indici nei reg di GP e li usa per i vettori in PINSRW.

Se gather potrebbe sfruttare i casi in cui più elementi avevano lo stesso indice, o anche un indice che puntava allo stesso blocco di recupero 32B, potrebbero esserci grandi accelerazioni a seconda dei dati di input.

Speriamo che Skylake migliorerà ulteriormente. Ho pensato di leggere qualcosa dicendo che sarebbe successo, ma al momento del controllo non ho trovato nulla.

RE: matrici sparse: non esiste un formato che duplica i dati, quindi è possibile eseguire letture contigue per righe o colonne? Non è qualcosa per cui ho dovuto scrivere codice, ma penso di averne parlato in alcune risposte.