Как выполнить матричную операцию 8 x 8 с помощью SSE?

Моя первоначальная попытка выглядела так (предположим, мы хотим умножить)

  __m128 mat[n]; /* rows */
  __m128 vec[n] = {1,1,1,1};
  float outvector[n];
   for (int row=0;row<n;row++) {
       for(int k =3; k < 8; k = k+ 4)
       {
           __m128 mrow = mat[k];
           __m128 v = vec[row];
           __m128 sum = _mm_mul_ps(mrow,v);
           sum= _mm_hadd_ps(sum,sum); /* adds adjacent-two floats */
       }
           _mm_store_ss(&outvector[row],_mm_hadd_ps(sum,sum));
 }

Но это явно не работает. Как мне подойти к этому?

Я должен загрузить 4 за раз ....

Другой вопрос: если мой массив очень большой (скажем, n = 1000), как я могу выровнять его по 16 байтам? Это вообще возможно?


person user1012451    schedule 27.11.2011    source источник
comment
Какой результат вы ожидаете? Я не вижу никакой матрицы, только векторное умножение. Кроме того, откуда берутся 3, 8 и 4?   -  person pezcode    schedule 27.11.2011
comment
@ user963889, размеры не имеют никакого смысла. Что ты пытаешься сделать? Умножить вектор 8x1 или массив векторов на матрицу 8x8?   -  person Brett Hale    schedule 27.11.2011
comment
@BrettHale Предположим, у нас есть 8x8, кратные вектору 8x1. Хочу в итоге получить 8х1. Я застрял. Можете ли вы, ребята, вести меня в правильном направлении? Спасибо.   -  person user1012451    schedule 28.11.2011
comment
@user963889 user963889, ОК - дал ответ на [8x8] x [8x1] ...   -  person Brett Hale    schedule 28.11.2011


Ответы (2)


ОК... Я буду использовать соглашение о мажорной матрице. Для каждой строки [m] требуется (2) __m128 элементов, чтобы получить 8 чисел с плавающей запятой. Вектор 8x1 v является вектором-столбцом. Поскольку вы используете инструкцию haddps, я предполагаю, что SSE3 доступен. Находка r = [m] * v :

void mul (__m128 r[2], const __m128 m[8][2], const __m128 v[2])
{
    __m128 t0, t1, t2, t3, r0, r1, r2, r3;

    t0 = _mm_mul_ps(m[0][0], v[0]);
    t1 = _mm_mul_ps(m[1][0], v[0]);
    t2 = _mm_mul_ps(m[2][0], v[0]);
    t3 = _mm_mul_ps(m[3][0], v[0]);

    t0 = _mm_hadd_ps(t0, t1);
    t2 = _mm_hadd_ps(t2, t3);
    r0 = _mm_hadd_ps(t0, t2);

    t0 = _mm_mul_ps(m[0][1], v[1]);
    t1 = _mm_mul_ps(m[1][1], v[1]);
    t2 = _mm_mul_ps(m[2][1], v[1]);
    t3 = _mm_mul_ps(m[3][1], v[1]);

    t0 = _mm_hadd_ps(t0, t1);
    t2 = _mm_hadd_ps(t2, t3);
    r1 = _mm_hadd_ps(t0, t2);

    t0 = _mm_mul_ps(m[4][0], v[0]);
    t1 = _mm_mul_ps(m[5][0], v[0]);
    t2 = _mm_mul_ps(m[6][0], v[0]);
    t3 = _mm_mul_ps(m[7][0], v[0]);

    t0 = _mm_hadd_ps(t0, t1);
    t2 = _mm_hadd_ps(t2, t3);
    r2 = _mm_hadd_ps(t0, t2);

    t0 = _mm_mul_ps(m[4][1], v[1]);
    t1 = _mm_mul_ps(m[5][1], v[1]);
    t2 = _mm_mul_ps(m[6][1], v[1]);
    t3 = _mm_mul_ps(m[7][1], v[1]);

    t0 = _mm_hadd_ps(t0, t1);
    t2 = _mm_hadd_ps(t2, t3);
    r3 = _mm_hadd_ps(t0, t2);

    r[0] = _mm_add_ps(r0, r1);
    r[1] = _mm_add_ps(r2, r3);
}

Что касается выравнивания, переменная типа __m128 должна автоматически выравниваться в стеке. С динамической памятью это небезопасное предположение. Некоторые реализации malloc/new могут возвращать только память, гарантированно выровненную по 8 байтам.

Заголовок встроенных функций предоставляет _mm_malloc и _mm_free. В этом случае параметр align должен быть (16).

person Brett Hale    schedule 28.11.2011
comment
Спасибо. Мой код теперь похож на ваш после 2 дней работы... но ваш намного понятнее. Я узнал. Спасибо. - person user1012451; 01.12.2011

Intel разработала Small Matrix Library для матриц размером от 1×1 до 6. ×6. Примечание по применению Расширения потоковой передачи SIMD для AP-930 — умножение матриц подробно описывает алгоритм умножения двух матриц 6×6. Это должно быть адаптировано к матрицам других размеров с некоторыми усилиями.

person Daniel Trebbien    schedule 27.11.2011