Матрично-векторное умножение в AVX не пропорционально быстрее, чем в SSE
Я писал матрично-векторное умножение в SSE и AVX, используя следующее:
for(size_t i=0;i<M;i++) {
size_t index = i*N;
__m128 a, x, r1;
__m128 sum = _mm_setzero_ps();
for(size_t j=0;j<N;j+=4,index+=4) {
a = _mm_load_ps(&A[index]);
x = _mm_load_ps(&X[j]);
r1 = _mm_mul_ps(a,x);
sum = _mm_add_ps(r1,sum);
}
sum = _mm_hadd_ps(sum,sum);
sum = _mm_hadd_ps(sum,sum);
_mm_store_ss(&C[i],sum);
}
Я использовал аналогичный метод для AVX, однако в конце, так как AVX не имеет эквивалентной инструкции для _mm_store_ss(), Я использовал:
_mm_store_ss(&C[i],_mm256_castps256_ps128(sum));
Код SSE дает мне ускорение 3,7 по сравнению с последовательным кодом. Тем не менее, код AVX дает мне ускорение всего 4,3 по сравнению с последовательным кодом.
Я знаю, что использование SSE с AVX может вызвать проблемы, но я скомпилировал его с флагом -mavx ‘, используя g ++, который должен удалить коды операций SSE.
Я мог бы также использовать: _mm256_storeu_ps(&C[i],sum) сделать то же самое, но ускорение то же самое.
Любое понимание того, что еще я мог бы сделать, чтобы улучшить производительность? Может ли это быть связано с: performance_memory_bound, хотя я не понял ответ на эту тему ясно.
Кроме того, я не могу использовать инструкцию _mm_fmadd_ps (), даже включив заголовочный файл «immintrin.h». У меня есть FMA и AVX включены.
Решение
Я предлагаю вам пересмотреть свой алгоритм. Смотрите обсуждение Эффективное матричное векторное умножение 4×4 с SSE: горизонтальное сложение и точечный продукт — какой смысл?
Вы делаете один длинный продукт и используете _mm_hadd_ps за итерацию. Вместо этого вы должны делать четыре точечных продукта одновременно с SSE (восемь с AVX) и использовать только вертикальные операторы.
Вам нужно сложение, умножение и трансляция. Все это может быть сделано в SSE с _mm_add_ps, _mm_mul_ps, а также _mm_shuffle_ps (для трансляции).
Если у вас уже есть транспонирование матрицы, это действительно просто.
Но независимо от того, есть ли у вас транспонирование или нет, вам нужно сделать свой код более дружественным к кешу. Чтобы это исправить, я предлагаю циклическое разбиение матрицы. Смотрите это обсуждение Какой самый быстрый способ транспонировать матрицу в C ++? чтобы получить представление о том, как сделать разбиение на петли.
Я бы попытался получить правильное разбиение на петли, прежде чем даже пытаться использовать SSE / AVX. Наибольший прирост в моем умножении матриц был не от SIMD или потоков, а от циклического разбиения. Я думаю, что если вы правильно используете кэш, ваш AVX-код будет работать более линейно по сравнению с SSE.
Другие решения
Как кто-то уже предложил, добавьте
-funroll-петли
Любопытно, что это не установлено по умолчанию.
Используйте __restrict для определения любых указателей с плавающей точкой.
Используйте const для константных ссылок на массивы.
Я не знаю, достаточно ли умен компилятор, чтобы признать, что 3 промежуточных значения внутри цикла не нужно поддерживать в рабочем состоянии от итерации к итерации.
Я бы просто удалил эти 3 переменные или, по крайней мере, сделал их локальными внутри цикла (a, x, r1). Индекс может быть объявлен, где j объявлен, чтобы сделать его более локальным.
Убедитесь, что M и N объявлены как const, и если их значения являются константами времени компиляции, дайте компилятору увидеть их.
Рассмотрим этот код. Я не знаком с версией INTEL, но это быстрее, чем XMMatrixMultiply в DirectX. Дело не в том, сколько математики сделано для каждой инструкции, а в том, чтобы уменьшить количество команд (если вы используете быстрые инструкции, что и делает эта реализация).
// Perform a 4x4 matrix multiply by a 4x4 matrix
// Be sure to run in 64 bit mode and set right flags
// Properties, C/C++, Enable Enhanced Instruction, /arch:AVX
// Having MATRIX on a 32 byte bundry does help performance
struct MATRIX {
union {
float f[4][4];
__m128 m[4];
__m256 n[2];
};
}; MATRIX myMultiply(MATRIX M1, MATRIX M2) {
MATRIX mResult;
__m256 a0, a1, b0, b1;
__m256 c0, c1, c2, c3, c4, c5, c6, c7;
__m256 t0, t1, u0, u1;
t0 = M1.n[0]; // t0 = a00, a01, a02, a03, a10, a11, a12, a13
t1 = M1.n[1]; // t1 = a20, a21, a22, a23, a30, a31, a32, a33
u0 = M2.n[0]; // u0 = b00, b01, b02, b03, b10, b11, b12, b13
u1 = M2.n[1]; // u1 = b20, b21, b22, b23, b30, b31, b32, b33
a0 = _mm256_shuffle_ps(t0, t0, _MM_SHUFFLE(0, 0, 0, 0)); // a0 = a00, a00, a00, a00, a10, a10, a10, a10
a1 = _mm256_shuffle_ps(t1, t1, _MM_SHUFFLE(0, 0, 0, 0)); // a1 = a20, a20, a20, a20, a30, a30, a30, a30
b0 = _mm256_permute2f128_ps(u0, u0, 0x00); // b0 = b00, b01, b02, b03, b00, b01, b02, b03
c0 = _mm256_mul_ps(a0, b0); // c0 = a00*b00 a00*b01 a00*b02 a00*b03 a10*b00 a10*b01 a10*b02 a10*b03
c1 = _mm256_mul_ps(a1, b0); // c1 = a20*b00 a20*b01 a20*b02 a20*b03 a30*b00 a30*b01 a30*b02 a30*b03
a0 = _mm256_shuffle_ps(t0, t0, _MM_SHUFFLE(1, 1, 1, 1)); // a0 = a01, a01, a01, a01, a11, a11, a11, a11
a1 = _mm256_shuffle_ps(t1, t1, _MM_SHUFFLE(1, 1, 1, 1)); // a1 = a21, a21, a21, a21, a31, a31, a31, a31
b0 = _mm256_permute2f128_ps(u0, u0, 0x11); // b0 = b10, b11, b12, b13, b10, b11, b12, b13
c2 = _mm256_mul_ps(a0, b0); // c2 = a01*b10 a01*b11 a01*b12 a01*b13 a11*b10 a11*b11 a11*b12 a11*b13
c3 = _mm256_mul_ps(a1, b0); // c3 = a21*b10 a21*b11 a21*b12 a21*b13 a31*b10 a31*b11 a31*b12 a31*b13
a0 = _mm256_shuffle_ps(t0, t0, _MM_SHUFFLE(2, 2, 2, 2)); // a0 = a02, a02, a02, a02, a12, a12, a12, a12
a1 = _mm256_shuffle_ps(t1, t1, _MM_SHUFFLE(2, 2, 2, 2)); // a1 = a22, a22, a22, a22, a32, a32, a32, a32
b1 = _mm256_permute2f128_ps(u1, u1, 0x00); // b0 = b20, b21, b22, b23, b20, b21, b22, b23
c4 = _mm256_mul_ps(a0, b1); // c4 = a02*b20 a02*b21 a02*b22 a02*b23 a12*b20 a12*b21 a12*b22 a12*b23
c5 = _mm256_mul_ps(a1, b1); // c5 = a22*b20 a22*b21 a22*b22 a22*b23 a32*b20 a32*b21 a32*b22 a32*b23
a0 = _mm256_shuffle_ps(t0, t0, _MM_SHUFFLE(3, 3, 3, 3)); // a0 = a03, a03, a03, a03, a13, a13, a13, a13
a1 = _mm256_shuffle_ps(t1, t1, _MM_SHUFFLE(3, 3, 3, 3)); // a1 = a23, a23, a23, a23, a33, a33, a33, a33
b1 = _mm256_permute2f128_ps(u1, u1, 0x11); // b0 = b30, b31, b32, b33, b30, b31, b32, b33
c6 = _mm256_mul_ps(a0, b1); // c6 = a03*b30 a03*b31 a03*b32 a03*b33 a13*b30 a13*b31 a13*b32 a13*b33
c7 = _mm256_mul_ps(a1, b1); // c7 = a23*b30 a23*b31 a23*b32 a23*b33 a33*b30 a33*b31 a33*b32 a33*b33
c0 = _mm256_add_ps(c0, c2); // c0 = c0 + c2 (two terms, first two rows)
c4 = _mm256_add_ps(c4, c6); // c4 = c4 + c6 (the other two terms, first two rows)
c1 = _mm256_add_ps(c1, c3); // c1 = c1 + c3 (two terms, second two rows)
c5 = _mm256_add_ps(c5, c7); // c5 = c5 + c7 (the other two terms, second two rose)
// Finally complete addition of all four terms and return the results
mResult.n[0] = _mm256_add_ps(c0, c4); // n0 = a00*b00+a01*b10+a02*b20+a03*b30 a00*b01+a01*b11+a02*b21+a03*b31 a00*b02+a01*b12+a02*b22+a03*b32 a00*b03+a01*b13+a02*b23+a03*b33
// a10*b00+a11*b10+a12*b20+a13*b30 a10*b01+a11*b11+a12*b21+a13*b31 a10*b02+a11*b12+a12*b22+a13*b32 a10*b03+a11*b13+a12*b23+a13*b33
mResult.n[1] = _mm256_add_ps(c1, c5); // n1 = a20*b00+a21*b10+a22*b20+a23*b30 a20*b01+a21*b11+a22*b21+a23*b31 a20*b02+a21*b12+a22*b22+a23*b32 a20*b03+a21*b13+a22*b23+a23*b33
// a30*b00+a31*b10+a32*b20+a33*b30 a30*b01+a31*b11+a32*b21+a33*b31 a30*b02+a31*b12+a32*b22+a33*b32 a30*b03+a31*b13+a32*b23+a33*b33
return mResult;
}
Еще раз, если вы хотите построить свой собственный алгоритм умножения матриц, пожалуйста остановись. Я помню, что на форуме Intel AVX один из их инженеров признался, что им потребовалось очень много времени для написания сборок AVX, чтобы достичь теоретической пропускной способности AVX для умножения двух матриц (особенно маленькие матрицы), поскольку инструкции загрузки и сохранения AVX в настоящее время довольно медленные, не говоря уже о трудностях преодоления накладных расходов потоков для параллельной версии.
Пожалуйста, установите Библиотека Intel Math Kernel, потратить полчаса на чтение руководства и строки кода 1, чтобы позвонить в библиотеку, СДЕЛАННЫЙ!
detector