Захотелось научиться писать БПФ без арифметики с плавающей точкой - начал изучать по этой статье. Вроде все более-менее понял, написал и даже кое-что работает, но остается непонятными несколько вещей:
Во-первых - пусть A = B*C(в смысле полиномов) и мы посчитали FFT(B), FFT(C). Все значения в этих массивах меньше модуля. Потом мы должны скалярно перемножить (FFT(B), FFT(C)) и применить обратное FFT. Но когда мы перемножим - не все элементы этого массива уже будут меньше модуля - что с этим делать? просто брать их по тому же модулю и считать обратное? Если я так делаю, у меня начинает все неправильно работать(из-за таких переполнений), числа длиной около тысячи знаков умножаются правильно только при базе 10 или 100, что очень печально.
Во-вторых - пусть после FFT у нас получилось: FFT(A) = (1, 3), FFT(B) = (4, 0). После того как мы скалярно перемножим, у нас получится (4, 0), и после обратного FFT мы получим тот же B. Такое бывает(может я не понял чего-то в теории), и если да - что делать?
И последнее - модуль в статье на e-maxx подходит для случая, когда в векторе не более 2^20 элементов. Вряд ли мне понадобится больше - но вопрос - а где-то есть готовые последовательности типа "модуль - образующий элемент" ?. На oeis - только для небольших чисел, можно конечно прочитать и познать еще один алгоритм - но очень уж лень :)
Кстати - мне кажется, или на e-maxx единственная статья про целочисленное БПФ в интернете? Если есть еще что почитать(на русском или английском), очень буду рад ссылкам.
upd: еще на всякий случай код
По первому вопросу -- надо брать модуль, который вместит произведение любой пары коэфициентов из перемножаемых многочленов.
По второму -- такое может быть только если и оригинальные многочлены обладали таким свойством (A*B=B).
По третьему -- на e-maxx это написано. Тебе надо найти простое число вида c*2^k+1 -- у него будет первообразный корень из единицы степени 2^k.
p=167772161 g=3
p=469762049 g=3 (это даже 226k+1)
p=1107296257 g=10
p=1711276033 g=29
p=1811939329 g=13 (это даже 226k+1)
p=2013265921 g=31 (это даже 227k+1)
p=2113929217 g=5
Можно делать еще такой трюк, чтобы избежать переполнения. Посчитать всё по двум модулям. А потом по КТО восстановить. Тогда можно брать хоть 106 за основание системы счисления. Из-за двух модулей все будет в два раза дольше, но за счёт большого основания получается все таки выигрыш во времени.
Да, с информацией по целочисленному ФФТ напряжёнка, но когда-то я находил это в одной книге. Сейчас попытался снова откопать - не смог :(
Но, видимо, ответы на все вопросы и так уже получены
http://arxiv.org/abs/math/0307321