О векторных расширениях gcc/clang (2)

[alextutubalin]

Кросспост из блога автора. Комментировать лучше там, но можно и тут

В комментариях к одному из предыдущих постов про оптимизацию матричного преобразования цвета нам предлагают немножко подумать над алгоритмом.
К сожалению, предложенное там решение (офигенно быстрое!) считает неправильно, но направление движение указано верно и мы приходим к такому варианту:

1. транспонируем матрицу, на которую умножаем, дополним нулями правую колонку, чтобы вышло 4x4
2. Каждое из (четырех) входных значений - размножим на вектор.
3. Нужный нам результат - это SIMD-сумма SIMD-произведений вышеупомянутых векторов на строки вышеупомянутой транспонированной матрицы.

Короче, проще кодом:

1. 
   typedef float __v4sf __attribute__ ((__vector_size__ (16)));
2. 
   __v4sf xmat2[4] =
3. 
   {{0.17f, 0.55f, 1.01f, 0.0f},
4. 
   {0.22f, 0.66f, 1.02f, 0.0f},
5. 
   {0.33f, 0.77f, 1.03f, 0.0f},
6. 
   {0.44f, 0.88f, 1.04f, 0.0f}};
7. 
    
8. 
   void dotp_vecT (float *d, int sz)
9. 
   {
10. 
      int i;
11. 
      __v4sf *data = (__v4sf *)d;
12. 
      __v4sf x0,x1,x2,x3,m0,m1,m2,m3;
13. 
      for(i=0;i
14. 
      {
15. 
            x0[0] = x0[1] = x0[2] = x0[3] = data[i][0];
16. 
            x1[0] = x1[1] = x1[2] = x1[3] = data[i][1];
17. 
            x2[0] = x2[1] = x2[2] = x2[3] = data[i][2];
18. 
            x3[0] = x3[1] = x3[2] = x3[3] = data[i][3];
19. 
            m0 = x0 * xmat2[0];
20. 
            m1 = x1 * xmat2[1];
21. 
            m2 = x2 * xmat2[2];
22. 
            m3 = x3 * xmat2[3];
23. 
            data[i] = m0+m1+m2+m3;
24. 
      }
25. 
    }

Результаты

Этот код компилировался gcc 4.6.2 и clang 3.0 (из svn) и запускался на Core2 Q9300 2.5Ghz. Результаты, мягко скажем, разные:

• gcc: 29 Mpix/sec. Это в 3.3 раза хуже чем целочисленный вариант на той же машине и в 5.7 раз хуже чем наилучший (на SSE4.1 dot product) из исследованных на данный момент.
• clang: 164 Mpix/sec, то есть так же, как написанный вручную SSE4.1-вариант.

В отличие от SSE4.1 варианта, этот код работает и для SSE2, причем скорость (на том же core2) не падает
Обсуждение

Разница в производительности объясняется разницей в коде (кто бы сомневался).
Clang нарисовал код, близкий к идеальному (кусочек делает первое умножение):

1. 
           movdqa  (%rdi), %xmm4
2. 
           pshufd  $85, %xmm4, %xmm5       # xmm5 = xmm4[1,1,1,1]
3. 
           mulps   %xmm0, %xmm5            # в xmm0 - вторая строка матрицы
4. 
           pshufd  $0, %xmm4, %xmm6        # xmm6 = xmm4[0,0,0,0]
5. 
           mulps   %xmm1, %xmm6            # в xmm1 - первая строка
6. 
           addps   %xmm5, %xmm6
7. 
   ...

Одна загрузка сразу 4 значений, дальше через shuffle размножаем, умножаем и сложаем. Я бы руками так же написал, ну может в другом порядке, но и только.
В gcc все безобразно. Загрузка данных происходит через стек, что все и объясняет:

1. 
           movl    (%rdi), %eax
2. 
           addl    $1, %edx
3. 
           movl    %eax, -60(%rsp)
4. 
           movl    %eax, -64(%rsp)
5. 
           movl    %eax, -68(%rsp)
6. 
           movl    %eax, -72(%rsp)
7. 
           movl    4(%rdi), %eax
8. 
           movaps  -72(%rsp), %xmm1
9. 
           mulps   %xmm4, %xmm1

И так двенадцать четыре раза (в xmm4 - строка матрицы). Более того, на AVX-машине, где есть прекрасная инструкция vbroadcastss, gcc ей не пользуется. Зато регистры экономятся! C или C++?

Еще один прикол gcc в том, что как C-текст вышеприведенный сниппет компилируется, а как C++ - нет: error: invalid types '__v4sf {aka __vector(4) float}[int]' for array subscript. Удобно, да. Мораль

Мораль, к сожалению, неутешительная. Из векторных расширений (если они используются хоть капельку нетривиально, понятно что SIMD-сложения/умножения работают) gcc может сделать такое, что лучше не надо. Вполне возможно, что из более сложного кода и clang может сделать эдакое, но пока не видел. Но если мы хот-спот векторно "оптимизируем", отчего это место начинает работать в разы хуже, то жить так тоже нельзя.
А значит - писать таки на SIMD-ассемблере. Под 3-4 разные архитектуры (и это я еще AMD не щупал). Если, конечно, интересует результат.