Złożoność obliczania dyskretnej transformaty Fouriera?

Jaka jest złożoność (na standardowej całkowitej liczbie pamięci RAM) obliczania standardowej dyskretnej transformaty Fouriera wektora $n$ liczb całkowitych?

Klasyczny algorytm szybkich transformacji Fouriera , niewłaściwie ^[1] przypisany Cooleyowi i Tukeyowi, jest zwykle opisywany jako działający w czasie $O(n \log n)$ . Ale większość operacji arytmetycznych wykonywanych w tym algorytmie rozpoczyna się od złożonych $n$ pierwiastków jedności, które są (dla większości $n$ ) irracjonalne, więc dokładna ocena w stałym czasie nie jest rozsądna. Ten sam problem pojawia się w przypadku naiwnego algorytmu czasu $O(n^2)$ (pomnożenie przez macierz Vandermonde'a złożonych pierwiastków jedności).

Nie jest nawet jasne, jak dokładnie reprezentować wynik DFT (w jakiejkolwiek przydatnej formie). Innymi słowy, nie jest jasne, czy obliczanie DFT jest rzeczywiście możliwe!

Więc załóżmy, musimy tylko $b$ bitów precyzji w każdej wartości wyjściowej. Jaka jest złożoność obliczania dyskretnej transformaty Fouriera w funkcji $n$ i $b$ ? (Dla konkretności możesz założyć, że $n$ jest potęgą $2$ ).

Czy może każdy przypadek „FFT” w literaturze oznacza „szybką transformację teoretyczną ”? ^[2]

Zobacz moje powiązane pytania dotyczące złożoności eliminacji Gaussa i najkrótszych ścieżek euklidesowych .

_{^[1] Naprawdę należy go nazwać (jakimś przedrostkiem) algorytmem Gaussa-Runge-Königa-Yatesa-Stumpfa-Danielsona-Lánczosa-Cooleya-Tukeya.}

_{^[2] A jeśli tak, to dlaczego większość podręczników opisuje tylko algorytm liczb zespolonych?}

Ta odpowiedź jest wariantem analizy pierwszego algorytmu („Methode A”) Schönhage'a i Strassena do mnożenia długich liczb całkowitych.

Załóżmy, że chcemy obliczyć FFT o długości . Skaluj dane wejściowe tak, aby wszystkie wartości były mniejsze niż 1. Załóżmy najpierw, że obliczamy z arytmetyką stałego punktu m -bit ( m bitów po punkcie binarnym). Niech δ = 2 1 / 2 - m jest ( "kompleks"), jednostka najmniej położeniu. Niech ω = exp ( 2 π i / K ) .$K = 2^k$$m$$m$$\delta = 2^{1/2 -m}$$\omega = \exp(2\pi i/K)$

1) Można obliczyć przybliżenia takie, że | ω ′ j - ω j | ≤ ( 2 k - 1 ) δ dla wszystkich 0 ≤ j ≤ K - 1 . Można tego dokonać w czasie O ( K M ( m ) ), gdzie M ( m ) to czas potrzebny do pomnożenia liczby m- bitowej. (patrz Knuth Vol. 2, wydanie trzecie, strona 309).$\omega_j'$$|\omega_j' - \omega^j| \le (2k-1)\delta$$0 \le j \le K-1$$O(K M(m))$$M(m)$$m$

Jeśli standardowa liczba całkowita RAM oznacza koszt logarytmiczny, to . Jeśli standardowa liczba całkowita RAM oznacza słowo RAM, to M ( m ) = O ( m ) . (Schönhage i Strassen pokazują w „Metodzie A”, jak zmniejszyć w czasie liniowym mnożenie liczb m- bitowych do m mnożenia liczb bitów O ( log m ) . To ostatnie można zrobić po kosztach jednostkowych.)$M(m) = O(m \log m)$$M(m) = O(m)$$m$$m$$O(\log m)$

2) Klasyczna metoda FFT Cooleya-Tukeya oblicza operacje w postaci . Stosujemy arytmetykę stałoprzecinkową m- bit, te opcje stają się a ′ = t r u n c a t e ( b ′ + ω ′ j c ′ ) . Jeśli wiemy$a = b + \omega^j c$$m$$a' = truncate(b' + \omega_j' c')$ i C ' aż do błędu ε , dostajemy ' do błędu 2 ε + $b'$$c'$$\epsilon$$a'$ .$2\epsilon + 2k\delta$

3) Za pomocą indukcji łatwo zauważyć, że otrzymujemy końcowy wynik z błędem $(2^k - 1) \cdot 2k\delta$ . Aby uzyskać precyzję na końcu, m ≥ k + log k + b + O ( 1 ) . $b$$m \ge k + \log k + b + O(1)$

4) Zatem końcowy czas pracy wynosi .$O(K k M(k+b))$

Powinno to również działać z liczbami zmiennoprzecinkowymi: 1) nadal można to zrobić za pomocą arytmetyki stałych punktów, 2) jest również prawdziwe w przypadku liczb zmiennoprzecinkowych.

Wydaje mi się, że w arytmetyki punktu stałego można to zrobić nawet szybciej. Najpierw redukujemy obliczenia FFT do mnożenia wielomianów za pomocą sztuczki Bluesteina. Długość współczynników potrzebnych do uzyskania pożądanej precyzji powinna wynosić . Następnie zmniejszamy mnożenie wielomianów do mnożenia długich liczb całkowitych. (Dołącz współczynniki do długiej liczby i rozdziel je blokami zerowymi o długości O ( k + b ) .) Długość liczb całkowitych wynosi O ( K ( k + b ) ) .$O(k + b)$$O(k+b)$$O(K(k+b))$

To nie jest pełna odpowiedź, ale mogę wskazać kilka istotnych artykułów, a także częściowo wyjaśnić, dlaczego nie jest tak łatwo wydobyć odpowiedź z konkretnego pytania z literatury.

Zacznę od pytania, dlaczego chcesz poznać odpowiedź na to pytanie? Zazwyczaj ludzie, którym zależy na tego rodzaju problemach, to ci, którzy faktycznie muszą wdrożyć wysokowydajny FFT do praktycznego zastosowania. Tacy ludzie mniej troszczą się o asymptotyczną złożoność w wyidealizowanym modelu obliczeniowym niż o maksymalizację wydajności przy określonych ograniczeniach sprzętowych i programowych. Na przykład twórcy najszybszej transformacji Fouriera na Zachodzie piszą w swoim artykule:

Najlepszy wybór zależy od szczegółów sprzętowych, takich jak liczba rejestrów, opóźnienie i przepustowość instrukcji, rozmiar i powiązanie pamięci podręcznych, struktura potoku procesora itp.

Są to problemy, z którymi teoretycy zazwyczaj nie chcą się zmrużyć, ale mają ogromne znaczenie w rzeczywistych implementacjach. Jeśli teoretyk oświadczy: „Odkryłem absolutnie najlepszą asymptotyczną złożoność bitów w modelu RAM”, praktykujący może powiedzieć: „To miłe”, ale może uznać taki teoretyczny wynik za bezużyteczny dla swoich celów.

Powiedziawszy to, uważam, że najlepiej jest spojrzeć na literaturę analityczną. Na przykład, Tasche i Zeuner dokładnie przyjrzeli się stabilności liczbowej algorytmu FFT. To może wciąż nie być dokładnie to, czego chcesz, ponieważ ogólny konsensus wśród praktyków wydaje się, że aby osiągnąć określoną liczbę liczbową precyzję, najlepszym praktycznym podejściem jest wstępne obliczenie pewnych liczb zwanych „współczynnikami zmienności” z wysoką dokładnością. Jeśli wykonujesz tylko jedną FFT, nie będzie to najszybsze podejście, ponieważ nie będziesz w stanie amortyzować kosztu jednorazowego wstępnego obliczenia w stosunku do dużej liczby obliczeń FFT. Mimo to ich analiza najgorszego przypadku błędu zaokrąglenia powinna być nadal odpowiednia dla twojego pytania.