Ускорение матрично-матричного умножения с использованием векторных инструкций SSE - PullRequest
0 голосов
/ 30 октября 2018

У меня возникли проблемы с векторизацией кода C с использованием векторных инструкций SSE. Код, который я должен победить, -

#define N 1000
void matrix_mul(int mat1[N][N], int mat2[N][N], int result[N][N])
{
   int i, j, k;
   for (i = 0; i < N; ++i)
   {
      for (j = 0; j < N; ++j)
      {
         for (k = 0; k < N; ++k)
         {
              result[i][k] += mat1[i][j] * mat2[j][k];
         }
      }
   }
}

Вот что я получил до сих пор:

void  matrix_mul_sse(int mat1[N][N], int mat2[N][N], int result[N][N])
{
   int i, j, k; int* l;
   __m128i v1, v2, v3;
   v3 = _mm_setzero_si128();
   for (i = 0; i < N; ++i)
   {
       for (j = 0; j < N; j += 4)
       {

           for (k = 0; k < N; k += 4)
           {

               v1 = _mm_set1_epi32(mat1[i][j]);
               v2 = _mm_loadu_si128((__m128i*)&mat2[j][k]);
               v3 = _mm_add_epi32(v3, _mm_mul_epi32(v1, v2));
               _mm_storeu_si128((__m128i*)&result[i][k], v3);
               v3 = _mm_setzero_si128();
           }
       }
   }
}

После исполнения я получил неверный результат. Я знаю, что причина в загрузке из памяти в v2. Я перебираю mat1 в главном порядке строк, поэтому мне нужно загрузить mat2 [0] [0], mat2 [1] [0], mat2 [2] [0], mat2 [3] [0] .... но что фактически загружается mat2 [0] [0], mat2 [0] [1], mat2 [0] [2], mat2 [0] [3] ... потому что mat2 хранится в памяти в главном порядке строк. Я пытался решить эту проблему, но без каких-либо улучшений. Может кто-нибудь помочь мне, пожалуйста.

Ответы [ 2 ]

0 голосов
/ 30 октября 2018

Я изменил ваш код, чтобы сделать адресацию явной [в этом случае это помогает].

#define N 100

Это заглушка для операции векторного умножения и накопления; у вас должна быть возможность заменить NV любым броском, который есть у вашего векторного модуля, и поместить здесь соответствующие коды операций.

#define NV 8
int Vmacc(int *A, int *B) {
   int i = 0;
   int x = 0;
   for (i = 0; i < NV; i++) {
        x += *A++ * *B++;
    }
    return x;
}

Это умножение имеет два заметных отклонения от нормы: 1. Он кэширует столбчатый вектор в смежный. 2. Он пытается вытолкнуть кусочки множителя в вектор-подобный функционал. Даже без использования векторного модуля это занимает половину времени простой версии только из-за лучшего использования кэша / предварительной выборки.

void mm2(int *A, int *B, int n, int *C) {
    int c, r;
    int stride = 0;
    int cache[N];
    for (c = 0; c < n; c++) {
        /* cache cumn i: */
        for (r = 0; r < n; r++) {
            cache[r] = B[c + r*n];
        }
        for (r = 0; r < n; r++) {
            int k = 0;
            int x = 0;
            int *Av = A + r*n;
            for (k = 0; k+NV-1 < n; k += NV) {
                x += Vmacc(Av+k, cache+k);
            }
            while (k < n) {
                x += Av[k] * cache[k];
                k++;
            }
            C[r*n + c] = x;
        }
    }
}
0 голосов
/ 30 октября 2018

Ниже исправлена ​​ваша реализация:

void  matrix_mul_sse(int mat1[N][N], int mat2[N][N], int result[N][N])
{
   int i, j, k;
   __m128i v1, v2, v3, v4; 
   for (i = 0; i < N; ++i)
   {
       for (j = 0; j < N; ++j) // 'j' must be incremented by 1
       {
           // read mat1 here because it does not use 'k' index
           v1 = _mm_set1_epi32(mat1[i][j]); 
           for (k = 0; k < N; k += 4)
           {   
               v2 = _mm_loadu_si128((const __m128i*)&mat2[j][k]);

               // read what's in the result array first as we will need to add it later to our calculations
               v3 = _mm_loadu_si128((const __m128i*)&result[i][k]);

               // use _mm_mullo_epi32 here instead _mm_mul_epi32 and add it to the previous result
               v4 = _mm_add_epi32(v3, _mm_mullo_epi32(v1, v2));

               // store the result
               _mm_storeu_si128((__m128i*)&result[i][k], v4);
           }
       }
   }
}

Короче говоря _mm_mullo_epi32 (требуется SSE4.1) дает 4 x int32 результата, в отличие от _mm_mul_epi32, который дает 2 x int64 результата. Если вы не можете использовать SSE4.1, взгляните на ответ здесь для альтернативного решения SSE2.

Полное описание от Intel Intrinsic Guide :

_mm_mullo_epi32: умножить упакованные 32-разрядные целые числа на a и b, получив промежуточные 64-разрядные целые числа, и сохранить младшие 32 бита промежуточных целых чисел в dst.

_mm_mul_epi32: Умножьте младшие 32-битные целые числа от каждого упакованного 64-битного элемента в a и b и сохраните 64-битные результаты со знаком в формате dst.

...