Come chiedere a GCC di srotolare completamente questo ciclo (ad esempio, sbucciare questo ciclo)?

Question

C'è un modo per istruire GCC (versione I utilizzata 4.8.4) per srotolare il ciclo while nella funzione di fondo completamente, ad esempio sbucciare questo ciclo? Il numero di iterazioni del ciclo è noto al momento della compilazione: 58.

Lasciatemi prima spiegare cosa ho provato.

Controllando ouput GAS:

gcc -fpic -O2 -S GEPDOT.c 

12 registri XMM0 - XMM11 vengono utilizzati. Se mi passa la bandiera -funroll-loops a gcc:

gcc -fpic -O2 -funroll-loops -S GEPDOT.c 

il ciclo è srotolato solo due volte. Ho controllato le opzioni di ottimizzazione GCC. GCC dice che -funroll-loop si accende -frame-registers, quindi, quando GCC srotola un ciclo, la sua scelta prioritaria per l'allocazione del registro è quella di utilizzare i registri "rimanenti". Ma ci sono solo 4 rimasti su XMM12 - XMM15, quindi GCC può solo srotolare 2 volte al suo meglio. Se fossero disponibili 48 anziché 16 registri XMM, GCC srotolerà il ciclo while 4 volte senza problemi.

Eppure ho fatto un altro esperimento. Per prima cosa ho srotolato manualmente il ciclo while due volte, ottenendo una funzione GEPDOT_2. Allora non c'è alcuna differenza affatto tra

gcc -fpic -O2 -S GEPDOT_2.c 

e

gcc -fpic -O2 -funroll-loops -S GEPDOT_2.c 

Da GEPDOT_2 già utilizzato tutti i registri, senza srotolamento è eseguita.

GCC si registra la ridenominazione per evitare potenziale false dipendenze introdotte. Ma so per certo che non ci sarà un tale potenziale nel mio GEPDOT; anche se esiste, non è importante. Ho provato a srotolare il loop da solo, e srotolare 4 volte è più veloce dello srotolamento 2 volte, più veloce di quanto non si possa srotolare. Naturalmente posso srotolare manualmente più volte, ma è noioso. GCC può farlo per me? Grazie.

// C file "GEPDOT.c" 
#include <emmintrin.h> 

void GEPDOT (double *A, double *B, double *C) { 
    __m128d A1_vec = _mm_load_pd(A); A += 2; 
    __m128d B_vec = _mm_load1_pd(B); B++; 
    __m128d C1_vec = A1_vec * B_vec; 
    __m128d A2_vec = _mm_load_pd(A); A += 2; 
    __m128d C2_vec = A2_vec * B_vec; 
    B_vec = _mm_load1_pd(B); B++; 
    __m128d C3_vec = A1_vec * B_vec; 
    __m128d C4_vec = A2_vec * B_vec; 
    B_vec = _mm_load1_pd(B); B++; 
    __m128d C5_vec = A1_vec * B_vec; 
    __m128d C6_vec = A2_vec * B_vec; 
    B_vec = _mm_load1_pd(B); B++; 
    __m128d C7_vec = A1_vec * B_vec; 
    A1_vec = _mm_load_pd(A); A += 2; 
    __m128d C8_vec = A2_vec * B_vec; 
    B_vec = _mm_load1_pd(B); B++; 
    int k = 58; 
    /* can compiler unroll the loop completely (i.e., peel this loop)? */ 
    while (k--) { 
    C1_vec += A1_vec * B_vec; 
    A2_vec = _mm_load_pd(A); A += 2; 
    C2_vec += A2_vec * B_vec; 
    B_vec = _mm_load1_pd(B); B++; 
    C3_vec += A1_vec * B_vec; 
    C4_vec += A2_vec * B_vec; 
    B_vec = _mm_load1_pd(B); B++; 
    C5_vec += A1_vec * B_vec; 
    C6_vec += A2_vec * B_vec; 
    B_vec = _mm_load1_pd(B); B++; 
    C7_vec += A1_vec * B_vec; 
    A1_vec = _mm_load_pd(A); A += 2; 
    C8_vec += A2_vec * B_vec; 
    B_vec = _mm_load1_pd(B); B++; 
    } 
    C1_vec += A1_vec * B_vec; 
    A2_vec = _mm_load_pd(A); 
    C2_vec += A2_vec * B_vec; 
    B_vec = _mm_load1_pd(B); B++; 
    C3_vec += A1_vec * B_vec; 
    C4_vec += A2_vec * B_vec; 
    B_vec = _mm_load1_pd(B); B++; 
    C5_vec += A1_vec * B_vec; 
    C6_vec += A2_vec * B_vec; 
    B_vec = _mm_load1_pd(B); 
    C7_vec += A1_vec * B_vec; 
    C8_vec += A2_vec * B_vec; 
    /* [write-back] */ 
    A1_vec = _mm_load_pd(C); C1_vec = A1_vec - C1_vec; 
    A2_vec = _mm_load_pd(C + 2); C2_vec = A2_vec - C2_vec; 
    A1_vec = _mm_load_pd(C + 4); C3_vec = A1_vec - C3_vec; 
    A2_vec = _mm_load_pd(C + 6); C4_vec = A2_vec - C4_vec; 
    A1_vec = _mm_load_pd(C + 8); C5_vec = A1_vec - C5_vec; 
    A2_vec = _mm_load_pd(C + 10); C6_vec = A2_vec - C6_vec; 
    A1_vec = _mm_load_pd(C + 12); C7_vec = A1_vec - C7_vec; 
    A2_vec = _mm_load_pd(C + 14); C8_vec = A2_vec - C8_vec; 
    _mm_store_pd(C,C1_vec); _mm_store_pd(C + 2,C2_vec); 
    _mm_store_pd(C + 4,C3_vec); _mm_store_pd(C + 6,C4_vec); 
    _mm_store_pd(C + 8,C5_vec); _mm_store_pd(C + 10,C6_vec); 
    _mm_store_pd(C + 12,C7_vec); _mm_store_pd(C + 14,C8_vec); 
    } 

---

aggiornamento 1

Grazie al commento di @ user3386109, vorrei estendere questa domanda un po '. @ user3386109 solleva una domanda molto buona. In realtà ho qualche dubbio sulla capacità del compilatore di allocare il registro in modo ottimale, quando ci sono così tante istruzioni parallele da programmare.

Personalmente ritengo che un modo affidabile sia quello di codificare prima il corpo del loop (che è la chiave per HPC) nell'assemblaggio in linea asm, quindi duplicarlo tutte le volte che voglio. Ho avuto un post impopolare all'inizio di quest'anno: inline assembly. Il codice era leggermente diverso perché il numero di iterazioni di loop, j, è un argomento di funzione quindi sconosciuto al momento della compilazione. In quel caso non riesco a srotolarlo completamente, quindi ho solo duplicato il codice assembly due volte e ho convertito il loop in un'etichetta e saltato.Si è scoperto che le prestazioni risultanti del mio assembly scritto sono circa il 5% più alte dell'assemblaggio generato dal compilatore, il che potrebbe suggerire che il compilatore non riesce ad allocare i registri nel modo previsto e ottimale.

Ero (e sono ancora) un bambino nella codifica di assiemi, quindi mi serve un buon caso di studio per imparare un po 'sull'assemblaggio di x86. Ma a lungo termine non sono incline a codificare GEPDOT con una grande proporzione per il montaggio. Ci sono principalmente tre ragioni:

1. asm assembly inline è stato critisized per non essere portabile. Anche se non capisco perché. Forse perché diverse macchine hanno diversi registri danneggiati?
2. Anche il compilatore sta migliorando. Quindi preferirei sempre l'ottimizzazione algoritmica e l'abitudine al codice C per aiutare il compilatore a generare buoni risultati;
3. L'ultimo motivo è più importante. Il numero di iterazioni potrebbe non essere sempre 58. Sto sviluppando una subroutine di fattorizzazione della matrice ad alte prestazioni. Per un fattore di blocco della cache nb, il numero di iterazioni sarà (nb-2). Non ho intenzione di inserire nb come argomento di funzione, come ho fatto nel post precedente. Questo è un parametro specifico della macchina che verrà definito come una macro. Quindi il numero di iterazioni è noto al momento della compilazione, ma può variare da macchina a macchina. Indovina quanto lavoro tedioso devo svolgere nello srotolamento del ciclo manuale per una varietà di nb. Quindi, se c'è un modo per istruire semplicemente il compilatore a sbucciare un ciclo, è grandioso.

Sarei molto apprezzato se si può anche condividere qualche esperienza nella produzione di librerie portatili ad alte prestazioni.

Answer 1

Prova a modificare i parametri di ottimizzatore:

gcc -O3 -funroll-loops --param max-completely-peeled-insns=1000 --param max-completely-peel-times=100 

Questo dovrebbe fare il trucco.

Answer 2

Questa non è una risposta, ma potrebbe essere di interesse per gli altri che cercano di vettorizzare le moltiplicazioni di matrici con GCC.

seguito, assumo c è una matrice 4 × 4 nella riga-ordine maggiore, un ea 4 riga, n -column matrice per colonna principale (trasposto), b è un 4-colonna, n -ROW matrice per righe maggiore, e l'operazione di calcolare è c =un × b + c, dove × denota moltiplicazione matriciale.

La funzione ingenuo per raggiungere questo obiettivo è

void slow_4(double  *c, 
      const double *a, 
      const double *b, 
      size_t  n) 
{ 
    size_t row, col, i; 

    for (row = 0; row < 4; row++) 
     for (col = 0; col < 4; col++) 
      for (i = 0; i < n; i++) 
       c[4*row+col] += a[4*i+row] * b[4*i+col]; 
} 

GCC genera codice abbastanza buono per SSE2/SSE3 utilizzando

#if defined(__SSE2__) || defined(__SSE3__) 

typedef double vec2d __attribute__((vector_size (2 * sizeof (double)))); 

void fast_4(vec2d  *c, 
      const vec2d *a, 
      const vec2d *b, 
      size_t  n) 
{ 
    const vec2d *const b_end = b + 2L * n; 

    vec2d s00 = c[0]; 
    vec2d s02 = c[1]; 
    vec2d s10 = c[2]; 
    vec2d s12 = c[3]; 
    vec2d s20 = c[4]; 
    vec2d s22 = c[5]; 
    vec2d s30 = c[6]; 
    vec2d s32 = c[7]; 

    while (b < b_end) { 
     const vec2d b0 = b[0]; 
     const vec2d b2 = b[1]; 
     const vec2d a0 = { a[0][0], a[0][0] }; 
     const vec2d a1 = { a[0][1], a[0][1] }; 
     const vec2d a2 = { a[1][0], a[1][0] }; 
     const vec2d a3 = { a[1][1], a[1][1] }; 
     s00 += a0 * b0; 
     s10 += a1 * b0; 
     s20 += a2 * b0; 
     s30 += a3 * b0; 
     s02 += a0 * b2; 
     s12 += a1 * b2; 
     s22 += a2 * b2; 
     s32 += a3 * b2; 
     b += 2; 
     a += 2; 
    } 

    c[0] = s00; 
    c[1] = s02; 
    c[2] = s10; 
    c[3] = s12; 
    c[4] = s20; 
    c[5] = s22; 
    c[6] = s30; 
    c[7] = s32; 
} 

#endif 

Per AVX, GCC può fare ancora meglio con

#if defined(__AVX__) || defined(__AVX2__) 

typedef double vec4d __attribute__((vector_size (4 * sizeof (double)))); 

void fast_4(vec4d  *c, 
      const vec4d *a, 
      const vec4d *b, 
      size_t  n) 
{ 
    const vec4d *const b_end = b + n; 

    vec4d s0 = c[0]; 
    vec4d s1 = c[1]; 
    vec4d s2 = c[2]; 
    vec4d s3 = c[3]; 

    while (b < b_end) { 
     const vec4d bc = *(b++); 
     const vec4d ac = *(a++); 
     const vec4d a0 = { ac[0], ac[0], ac[0], ac[0] }; 
     const vec4d a1 = { ac[1], ac[1], ac[1], ac[1] }; 
     const vec4d a2 = { ac[2], ac[2], ac[2], ac[2] }; 
     const vec4d a3 = { ac[3], ac[3], ac[3], ac[3] }; 
     s0 += a0 * bc; 
     s1 += a1 * bc; 
     s2 += a2 * bc; 
     s3 += a3 * bc; 
    } 

    c[0] = s0; 
    c[1] = s1; 
    c[2] = s2; 
    c[3] = s3; 
} 

#endif 

La versione SSE3 dell'assieme generato utilizzando gcc-4.8.4 (-O2 -march=x86-64 -mtune=generic -msse3) è essenzialmente

fast_4: 
     salq $5, %rcx 
     movapd (%rdi), %xmm13 
     addq %rdx, %rcx 
     cmpq %rcx, %rdx 
     movapd 16(%rdi), %xmm12 
     movapd 32(%rdi), %xmm11 
     movapd 48(%rdi), %xmm10 
     movapd 64(%rdi), %xmm9 
     movapd 80(%rdi), %xmm8 
     movapd 96(%rdi), %xmm7 
     movapd 112(%rdi), %xmm6 
     jnb  .L2 
.L3: 
     movddup (%rsi), %xmm5 
     addq $32, %rdx 
     movapd -32(%rdx), %xmm1 
     addq $32, %rsi 
     movddup -24(%rsi), %xmm4 
     movapd %xmm5, %xmm14 
     movddup -16(%rsi), %xmm3 
     movddup -8(%rsi), %xmm2 
     mulpd %xmm1, %xmm14 
     movapd -16(%rdx), %xmm0 
     cmpq %rdx, %rcx 
     mulpd %xmm0, %xmm5 
     addpd %xmm14, %xmm13 
     movapd %xmm4, %xmm14 
     mulpd %xmm0, %xmm4 
     addpd %xmm5, %xmm12 
     mulpd %xmm1, %xmm14 
     addpd %xmm4, %xmm10 
     addpd %xmm14, %xmm11 
     movapd %xmm3, %xmm14 
     mulpd %xmm0, %xmm3 
     mulpd %xmm1, %xmm14 
     mulpd %xmm2, %xmm0 
     addpd %xmm3, %xmm8 
     mulpd %xmm2, %xmm1 
     addpd %xmm14, %xmm9 
     addpd %xmm0, %xmm6 
     addpd %xmm1, %xmm7 
     ja  .L3 
.L2: 
     movapd %xmm13, (%rdi) 
     movapd %xmm12, 16(%rdi) 
     movapd %xmm11, 32(%rdi) 
     movapd %xmm10, 48(%rdi) 
     movapd %xmm9, 64(%rdi) 
     movapd %xmm8, 80(%rdi) 
     movapd %xmm7, 96(%rdi) 
     movapd %xmm6, 112(%rdi) 
     ret 

La versione AVX del complesso generato (-O2 -march=x86-64 -mtune=generic -mavx) è essenzialmente

fast_4: 
     salq  $5, %rcx 
     vmovapd (%rdi), %ymm5 
     addq  %rdx, %rcx 
     vmovapd 32(%rdi), %ymm4 
     cmpq  %rcx, %rdx 
     vmovapd 64(%rdi), %ymm3 
     vmovapd 96(%rdi), %ymm2 
     jnb  .L2 
.L3: 
     addq  $32, %rsi 
     vmovapd -32(%rsi), %ymm1 
     addq  $32, %rdx 
     vmovapd -32(%rdx), %ymm0 
     cmpq  %rdx, %rcx 
     vpermilpd $0, %ymm1, %ymm6 
     vperm2f128 $0, %ymm6, %ymm6, %ymm6 
     vmulpd  %ymm0, %ymm6, %ymm6 
     vaddpd  %ymm6, %ymm5, %ymm5 
     vpermilpd $15, %ymm1, %ymm6 
     vperm2f128 $0, %ymm6, %ymm6, %ymm6 
     vmulpd  %ymm0, %ymm6, %ymm6 
     vaddpd  %ymm6, %ymm4, %ymm4 
     vpermilpd $0, %ymm1, %ymm6 
     vpermilpd $15, %ymm1, %ymm1 
     vperm2f128 $17, %ymm6, %ymm6, %ymm6 
     vperm2f128 $17, %ymm1, %ymm1, %ymm1 
     vmulpd  %ymm0, %ymm6, %ymm6 
     vmulpd  %ymm0, %ymm1, %ymm0 
     vaddpd  %ymm6, %ymm3, %ymm3 
     vaddpd  %ymm0, %ymm2, %ymm2 
     ja   .L3 
.L2: 
     vmovapd %ymm5, (%rdi) 
     vmovapd %ymm4, 32(%rdi) 
     vmovapd %ymm3, 64(%rdi) 
     vmovapd %ymm2, 96(%rdi) 
     vzeroupper 
     ret 

La schedulazione registro non è ottimale, immagino, ma non sembra atroce neanche. Sono personalmente soddisfatto di quanto sopra, senza tentare di ottimizzarlo a mano a questo punto.

Su un processore Core i5-4200U (compatibile con AVX2), le versioni veloci delle funzioni sopra riportate calcolano il prodotto di due matrici 4 × 256 in 1843 cicli CPU (mediana) per SSE3 e 1248 cicli per AVX2. Questo si riduce a 1,8 e 1,22 cicli per immissione matrice. La versione lenta non programmata richiede circa 11 cicli per immissione matrice, per confronto.

(I conteggi di ciclo sono valori medi, vale a dire la metà dei test erano più veloci. Ho corso solo alcune analisi comparativa di massima con ~ ripetizioni 100k o giù di lì, in modo da prendere questi numeri con un grano di sale.)

Su questo CPU, gli effetti di cache sono tali che AVX2 a 4 × 512 dimensioni della matrice è ancora a 1,2 cicli per voce, ma a 4 × 1024, scende a 1,4, a 4 × 4096 a 1,5, a 4 × 8192 a 1,8 e a 4 × 65536 a 2,2 cicli per entrata. La versione SSE3 rimane a 1,8 cicli per voce fino a 4 × 3072, a quel punto inizia a rallentare; a 4 × 65536 anche questo è di circa 2,2 cicli per entrata. Credo che questa CPU (portatile!) Abbia una larghezza di banda della cache limitata a questo punto.