SSE _mm_load_pd działa przy _mm_store_pd segfaults

Question

Próbuję nauczyć się sznurków SSE w C. Mam fragment kodu, w którym ładuję dwukomponentowy wektor podwójnych danych, dodam do niego coś, a następnie próbuję go zapisać z powrotem do pamięci. Wszystko działa: mogę załadować moje dane do rejestrów SEE, mogę operować na moich danych w tych rejestrach SSE, ale w chwili, gdy spróbuję zapisać te przetworzone dane z powrotem do oryginalnej tablicy (gdzie czytam moje dane z pierwsze miejsce!) Otrzymuję błąd segmentacji.

Czy ktoś może mi doradzić w tej sprawie - to doprowadza mnie do szaleństwa.

double res[2] __attribute__((aligned(16))); 

for(int k=0; k<n; k++){ 
int i=0; 
for(; i+1<n; i+=2) 
    { 
    __m128d cik = _mm_load_pd(&C[i+k*n]); 
    int j = 0; 
    for(; j+1<n; j+=2) 
     { 
     __m128d aTij = _mm_load_pd(&A_T[j+i*n]); 
     __m128d bjk = _mm_load_pd(&B[j+k*n]); 
     __m128d dotpr = _mm_dp_pd(aTij, bjk,2); 
     cik = _mm_add_pd(cik, dotpr); 
     } 
    _mm_store_pd(res, cik); 
    //C[i+k*n] = res[0]; 
    } 
} 

Jak powiedziałem wyżej, wszystko działa w tym kodzie wyjątkiem gdzie przechowywać moje wyniki z powrotem do tej jednowymiarowej tablicy „C” gdzie mogę czytać dane ze w pierwszej kolejności. Oznacza to, że kiedy usunąć znaki komentarza przed

//C[i+k*n] = res[0]; 

otrzymuję winy segmentacji.

Jak to jest możliwe, że mogę odczytać z C z wyrównaną wersją pamięci _mm_load_pd (więc C musi być wyrównane w pamięci!) Podczas pisania z powrotem do niego nie działa? "C" musi być wyrównane, a jak widać "res" musi być wyrównany.

Zastrzeżenie: Mój oryginalny kod odczytać

_mm_store_pd(&C[i+k*n], cik); 

który wyprodukował także usterki segmentacji i zacząłem wprowadzania „RES” z wyraźną wyrównania w moim próba rozwiązania problemu.

Uzupełnienie

A, B, C są zadeklarowane w następujący sposób:

buf = (double*) malloc (3 * nmax * nmax * sizeof(double)); 
double* A = buf + 0; 
double* B = A + nmax*nmax; 
double* C = B + nmax*nmax; 

Próba rozwiązania z posix_memalign

Starając się rozwiązać problem błędu segmentacji pisząc do oryginalna tablica jednowymiarowa, teraz używam buforów dla odpowiednich macierzy. Jednak wciąż się odkłada podczas próby zapisu z powrotem na C_buff!

double res[2] __attribute__((aligned(16))); 

double * A_T; 
posix_memalign((void**)&A_T, 16, n*n*sizeof(double)); 

double * B_buff; 
posix_memalign((void**)&B_buff, 16, n*n*sizeof(double)); 

double * C_buff; 
posix_memalign((void**)&C_buff, 16, n*n*sizeof(double)); 

for(int y=0; y<n; y++) 
    for(int x=0; x<n; x++) 
    A_T[x+y*n] = A[y+x*n]; 

for(int x=0; x<n; x++) 
    for(int y=0; y<n; y++) 
    B_buff[y+x*n] = B[y+x*n]; 

for(int x=0; x<n; x++) 
    for(int y=0; y<n; y++) 
    C_buff[y+x*n] = C[y+x*n]; 

for(int k=0; k<n; k++){ 
    int i=0; 
    for(; i+1<n; i+=2) 
    { 
     __m128d cik = _mm_load_pd(&C_buff[i+k*n]); 
     int j = 0; 
     for(; j+1<n; j+=2) 
     { 
      __m128d aTij = _mm_load_pd(&A_T[j+i*n]); 
      __m128d bjk = _mm_load_pd(&B_buff[j+k*n]); 
      __m128d dotpr = _mm_dp_pd(aTij, bjk,2); 
      cik = _mm_add_pd(cik, dotpr); 
     } 
     _mm_store_pd(&C_buff[i+k*n], cik); 

    //_mm_store_pd(res, cik); 
     //C_buff[i+k*n] = res[0]; 
    //C_buff[i+1+k*n] = res[1]; 
    } 
} 

Answer 1

Nawet z __attribute__((aligned(32))), byłem coraz ten sam błąd (miał 50% szans na mis-alingment). Następnie użyłem następującą funkcję dostać% 100 szanse na wyrównanie (a zasilania powinny być dwa):

void * malloc_float_align(size_t n, unsigned int a/*alignment*/, float *& output) 
{ 
    void * adres=NULL; 
    void * adres2=NULL; 
    adres=malloc(n*sizeof(float)+a); 
    size_t adr=(size_t)adres; 
    size_t adr2=adr+a-(adr&(a-1u)); // a valid address for a alignment 
    adres2=(void *) adr2; 
    output=(float *)adres2; 
    return adres;    //pointer to be used in free() 
} 

następnie wykorzystanie w głównej:

int main() 
{ 


    float * res=NULL; 
    void * origin=malloc_float_align(1024,32u,res); 
    //use res for sse/avx 
    free(origin); // actual allocation is more than 1024 elements 
    return 0; 
} 

Oczywiście to jest w C++, dzięki czemu można sprawiają, że działa tylko zmieniając niektóre style funkcji.

Answer 2

Prosty trik byłoby wykonać dochodzić i zobacz czy wyzwala:

ASSERT(((size_t)(&C_buff[i+k*n]) & 0xF) == 0); 

assert będzie ogień, gdy adres nie jest SSE wyrównany. 64-bitowe kompilacje powinny domyślnie zapewniać wyrównanie 16B. Jeśli planujesz mieć kod 32-bitowy, użyj jednej z powyższych funkcji align_malloc. Musisz użyć odpowiedniego align_free, bo się zawiesisz.

Answer 3

Po usunięciu _mm_store_pd(&C_buff[i+k*n], cik); cała pętla jest zoptymalizowana i usunięta. Kompilator wywnioskuje, że cała pętla for nie prowadzi do żadnej znaczącej pracy i jej nie usuwa. Właśnie dlatego nie masz już błędu segmentacji.
Jestem pewny, że błąd segmentacji wynika z rozmiaru tablicy. Rozważmy prosty program oparty na przykład:

#include <stdio.h> 
#include "emmintrin.h" 

int main(){ 
int n = 15; 
int y,x,k,i,j; 

double * A; 
posix_memalign((void**)&A, 16, n*n*sizeof(double)); 

double * B; 
posix_memalign((void**)&B, 16, n*n*sizeof(double)); 

double * C; 
posix_memalign((void**)&C, 16, n*n*sizeof(double)); 

for(y=0; y<n; y++) 
    for(x=0; x<n; x++) 
    A[x+y*n] = 0.1; 

for(x=0; x<n; x++) 
    for(y=0; y<n; y++) 
    B[y+x*n] = 0.1; 

for(x=0; x<n; x++) 
    for(y=0; y<n; y++) 
    C[y+x*n] = 0.1; 

for(k=0; k<n; k++){ 
    i=0; 
    for(; i+1<n; i+=2) 
    { 
     __m128d cik = _mm_load_pd(&C[i+k*n]); 
     j = 0; 
     for(; j+1<n; j+=2) 
     { 
      __m128d aTij = _mm_load_pd(&A[j+i*n]); 
      __m128d bjk = _mm_load_pd(&B[j+k*n]); 
      __m128d dotpr = _mm_add_pd(aTij, bjk); 
      cik = _mm_add_pd(cik, dotpr); 
     } 
     _mm_store_pd(&C[i+k*n], cik); 
    } 
} 
printf("C[15]: %f\n", C[15]); 
printf("C[14]: %f\n", C[14]); 

To daje winy segmentacji bo n jest nieparzysta. Teraz zmień n = 15 na n = 16 i wszystko będzie działać zgodnie z oczekiwaniami. W związku z tym dopełnienie tablic do parzystej liczby (a nawet lepiej, do rozmiaru linii pamięci podręcznej -> 64 bajty == 8 elementów DP lub 16 elementów SP) zapobiegnie takim problemom i doprowadzi do lepszej wydajności.