Использование инструкций FMA (слитое умножение) для сложного умножения
Я хотел бы использовать имеющиеся объединенные инструкции процессора для сложения и умножения, чтобы помочь в сложном умножении массива приличного размера. По сути, основная математика выглядит так:
void ComplexMultiplyAddToArray(float* pDstR, float* pDstI, const float* pSrc1R, const float* pSrc1I, const float* pSrc2R, const float* pSrc2I, int len)
{
for (int i = 0; i < len; ++i)
{
const float fSrc1R = pSrc1R[i];
const float fSrc1I = pSrc1I[i];
const float fSrc2R = pSrc2R[i];
const float fSrc2I = pSrc2I[i];
// Perform complex multiplication on the input and accumulate with the output
pDstR[i] += fSrc1R*fSrc2R - fSrc1I*fSrc2I;
pDstI[i] += fSrc1R*fSrc2I + fSrc2R*fSrc1I;
}
}
Как вы, вероятно, видите, данные структурированы так, что у нас есть отдельные массивы действительных и мнимых чисел. Теперь предположим, что у меня есть следующие функции, доступные как встроенные функции для отдельных инструкций, которые выполняютб + с и аб-в соответственно:
float fmadd(float a, float b, float c);
float fmsub(float a, float b, float c);
Наивно, я вижу, что я могу заменить 2 умножения, одно сложение и одно вычитание одним fmadd и одним fmsub, вот так:
// Perform complex multiplication on the input and accumulate with the output
pDstR[i] += fmsub(fSrc1R, fSrc2R, fSrc1I*fSrc2I);
pDstI[i] += fmadd(fSrc1R, fSrc2I, fSrc2R*fSrc1I);
Это приводит к очень скромным улучшениям производительности, а также, как я полагаю, точности, но я думаю, что мне действительно не хватает чего-то, где математическая математика может быть изменена так, что я могу заменить еще пару комбинаций «мульт / добавление» или «мульт / суб». В каждой строке есть дополнительное прибавление и дополнительное умножение, которое, как мне кажется, я могу преобразовать в одну FMA, но, к сожалению, я не могу понять, как это сделать, не меняя порядок операций и получая неправильный результат. Есть математические эксперты с идеями?
Ради вопроса, целевая платформа, вероятно, не так важна, так как я знаю, что такого рода инструкции существуют на разных платформах.
Решение
Это хорошее начало. Вы можете уменьшить еще одно дополнение:
// Perform complex multiplication on the input and accumulate with the output
pDstR[i] += fmsub(fSrc1R, fSrc2R, fSrc1I*fSrc2I);
pDstI[i] += fmadd(fSrc1R, fSrc2I, fSrc2R*fSrc1I);
Здесь вы можете использовать другой fmadd в расчете на мнимую часть:
pDstI[i] = fmadd(fSrc1R, fSrc2I, fmadd(fSrc2R, fSrc1I, pDstI[i]));
Точно так же вы можете сделать то же самое с реальной частью, но вам нужно отрицать аргумент. Если это делает вещи быстрее или медленнее, во многом зависит от микросинхронизации архитектуры, над которой вы работаете:
pDstR[i] = fmsub(fSrc1R, fSrc2R, fmadd(fSrc1I, fSrc2I, -pDstR[i]));
Кстати, вы может получить дополнительные улучшения производительности, если вы объявите свои целевые массивы как не сглаживающие с помощью restrict ключевое слово. В данный момент компилятор должен предположить, что pDstR и pDstI могут перекрываться или указывать на один и тот же кусок памяти. Это предотвратит загрузку компилятором pDstI [i] перед записью в pDstR [i].
После этого может помочь некоторая осторожная развертка цикла, если компилятор этого еще не сделал. Проверьте вывод ассемблера вашего компилятора!
Другие решения
Я нашел следующее (с небольшой помощью), кажется, приводит к правильному ответу:
pDstR[i] = fmsub(fSrc1R, fSrc2R, fmsub(fSrc1I, fSrc2I, pDstR[i]));
pDstI[i] = fmadd(fSrc1R, fSrc2I, fmadd(fSrc2R, fSrc1I, pDstI[i]));
Но, как ни странно, на AVX производительность не улучшается так сильно, как если бы в математической математике оставалась часть с реальным результатом, используя половину FMA, а при воображаемом результате используйте полную FMA:
pDstR[i] += fmsub(fSrc1R, fSrc2R, fSrc1I*fSrc2I);
pDstI[i] = fmadd(fSrc1R, fSrc2I, fmadd(fSrc2R, fSrc1I, pDstI[i]));
Спасибо всем за помощь.