Использование инструкций FMA для алгоритма FFT - PullRequest
5 голосов
/ 26 марта 2020

У меня есть немного кода на C ++, который со временем стал несколько полезной библиотекой FFT, и он был создан для приличной скорости работы с использованием инструкций SSE и AVX. Конечно, все это основано только на алгоритме radix-2, но оно все еще работает. Моя последняя проблема - заставить вычисления бабочки работать с инструкциями FMA. Бабочка basi c radix-2 состоит из 4 умножений и 6 сложений или вычитаний. Простой подход подразумевал бы замену 2 сложений и вычитаний и 2 умножения на 2 инструкции FMA, что привело бы к математически идентичной бабочке, но, очевидно, есть более эффективные способы сделать это:

https://books.google.com/books?id=2HG0DwAAQBAJ&pg=PA56&lpg=PA56&dq=radix+2+fft+fma&source=bl&ots=R5XDWyYBVv&sig=ACfU3U0S2n1hcgiP63LTKMxI5Oc85eEZaQ&hl=en&sa=X&ved=2ahUKEwiz_I3PsrToAhVoHzQIHYmVDGIQ6AEwDXoECAoQAQ#v = onepage & q = radix% 202% 20fft% 20fma & f = false

ci1 = ci1 / cr1
u0 = zinr(0)
v0 = zini(0)
r = zinr(1)
s = sini(1)
u1 = r - s * ci1
v1 = r * ci1 + s
zoutr(0) = u0 + u1 * cr1
zouti(0) = v0 + v1 * cr1
zoutr(1) = u0 - u1 * cr1
zouti(1) = v0 - v1 * cr1

Автор заменяет все 10 сложений, подпрограмм и мультов на 6 FMA при условии, что мнимая часть тиддл-фактора делится на реальная часть. Часть текста гласит «Обратите внимание, что cr1! = 0». По сути, это моя проблема в двух словах. Математика, кажется, работает так же, как рекламируется для всех факторов твидла, за исключением случаев, когда реальный твидл равен нулю, и в этом случае мы делим на ноль. Там, где эффективность абсолютно необходима, ветвление кода, когда cr1 == 0 для другой бабочки, не является хорошим вариантом, особенно когда мы используем SIMD для одновременной обработки нескольких скручиваний и бабочек, где, возможно, только один элемент cr1 == 0. То, что говорит мне моя интуиция, должно быть так: когда cr1 == 0, cr1 и ci1 должны быть полностью другими значениями, и код FMA все равно приведет к правильному ответу, но я не могу понять это , Если бы я мог понять это, было бы относительно просто изменить предварительно вычисленные коэффициенты твида для бабочек FMA, и мы также, конечно, могли бы избежать операции деления в начале бабочки.

1 Ответ

1 голос
/ 28 марта 2020

Книга, кажется, предполагает, что cr1 != 0 всегда верно. Но, к сожалению, это не всегда так (когда угол поворота равен PI / 2).

Я не думаю, что вы можете решить эту проблему, регулируя коэффициенты твида. Единственный вариант, который я вижу, это использовать какое-то очень маленькое число вместо нуля. Это может сработать, но это некрасиво и в некоторых случаях может привести к неточностям.

Возможные решения:

  • Разделите l oop на две части и обработайте этот центральный регистр ( где происходит деление на ноль) специально
  • Вместо деления на cr1, делим на ci1 и соответствующим образом модифицируем форум. В этом случае все еще есть деление на ноль, но это произойдет на первой итерации l oop. Таким образом, вместо центра вы должны обрабатывать первую итерацию специально (поэтому требуется только одна l oop).
  • Используйте другую формулировку FMA:

Обратите внимание, что :

zoutr(1) = u0 - u1 
         = u0 - u1 - (u0 + u1) + (u0 + u1) 
         = u0 - u1 - zoutr(0) + u0 + u1 
         = 2*u0 - zoutr(0)

Итак, эту операцию можно выполнить в 1 FMA.

И если подставить u1 в выражение zoutr(0):

zoutr(0) = u0 + u1
         = u0 + r*cr1 - s*ci1

Это может быть сделано с 2 FMA.

Расчет zouti может быть выполнен так же, как zoutr. Таким образом, вам нужно использовать 6 операций FMA, то есть столько же операций, сколько в книге.

(Обратите внимание, это не означает, что этот вариант будет работать быстрее автоматически, так как он имеет другой цепочка зависимостей данных)

...