Algoritmo FFT de Rader

Transformada discreta de Fourier para tamaños primos

El algoritmo de Rader (1968), ^[1] llamado así por Charles M. Rader del Laboratorio Lincoln del MIT , es un algoritmo de transformada rápida de Fourier (FFT) que calcula la transformada discreta de Fourier (DFT) de tamaños primos reexpresando la DFT como una transformación cíclica. convolución (el otro algoritmo para FFT de tamaños primos, el algoritmo de Bluestein , también funciona reescribiendo la DFT como una convolución).

Dado que el algoritmo de Rader solo depende de la periodicidad del núcleo DFT, es directamente aplicable a cualquier otra transformada (de orden primo) con una propiedad similar, como una transformada de teoría de números o la transformada discreta de Hartley .

El algoritmo se puede modificar para obtener un factor de ahorro de dos para el caso de DFT de datos reales, utilizando una reindexación/permutación ligeramente modificada para obtener dos convoluciones cíclicas de datos reales de tamaño medio; ^[2] una adaptación alternativa para DFT de datos reales utiliza la transformada discreta de Hartley . ^[3]

Winograd amplió el algoritmo de Rader para incluir tamaños DFT de potencia primaria , ^{[4] [5]} y hoy en día el algoritmo de Rader se describe a veces como un caso especial del algoritmo FFT de Winograd , también llamado algoritmo de transformada multiplicativa de Fourier (Tolimieri et al., 1997), ^[6] que se aplica a una clase de tamaños aún mayor. Sin embargo, para tamaños compuestos como las potencias primarias, el algoritmo FFT de Cooley-Tukey es mucho más simple y práctico de implementar, por lo que el algoritmo de Rader generalmente solo se usa para casos base de primos grandes de la descomposición recursiva de la DFT de Cooley-Tukey . ^[3] ${\displaystyle p^{m}}$

Algoritmo

Representación visual de una matriz DFT en el algoritmo FFT de Rader. La matriz consta de relojes de colores que representan una matriz DFT de tamaño 11. Al permutar filas y columnas (excepto la primera de cada una) según secuencias generadas por las potencias de la raíz primitiva de 11, la matriz DFT original se convierte en una matriz circulante . Multiplicar una secuencia de datos con una matriz circulante es equivalente a la convolución cíclica con el vector fila de la matriz. Esta relación es un ejemplo del hecho de que el grupo multiplicativo es cíclico: . ${\displaystyle (\mathbb {Z} /p\mathbb {Z} )^{\times }\cong C_{p-1}}$

Comience con la definición de la transformada discreta de Fourier:

${\displaystyle X_{k}=\sum _{n=0}^{N-1}x_{n}e^{-{\frac {2\pi i}{N}}nk}\qquad k=0,\dots ,N-1.}$

Si N es un número primo, entonces el conjunto de índices distintos de cero forma un grupo bajo el módulo de multiplicación N. Una consecuencia de la teoría de números de tales grupos es que existe un generador del grupo (a veces llamado raíz primitiva , que puede encontrarse mediante una búsqueda exhaustiva o algoritmos ligeramente mejores ^[7] ). Este generador es un número entero g tal que para cualquier índice distinto de cero n y para un único (formando una biyección de q a n distinto de cero ). De manera similar, para cualquier índice k distinto de cero y para un único , donde el exponente negativo denota el inverso multiplicativo de . Eso significa que podemos reescribir la DFT usando estos nuevos índices p y q como: ${\displaystyle n\in {}\{1,\dots ,N-1\}}$ ${\displaystyle n=g^{q}{\pmod {N}}}$ ${\displaystyle q\in {}\{0,\dots ,N-2\}}$ ${\displaystyle k=g^{-p}{\pmod {N}}}$ ${\displaystyle p\in {}\{0,\dots ,N-2\}}$ ${\displaystyle g^{p}\mod N}$

${\displaystyle X_{0}=\sum _{n=0}^{N-1}x_{n},}$

${\displaystyle X_{g^{-p}}=x_{0}+\sum _{q=0}^{N-2}x_{g^{q}}e^{-{\frac {2\pi i}{N}}g^{-(p-q)}}\qquad p=0,\dots ,N-2.}$

(Recuerde que x _n y X _k son implícitamente periódicos en N , y también que ( identidad de Euler ). Por lo tanto, todos los índices y exponentes se toman en módulo N como lo requiere la aritmética de grupos). ${\displaystyle e^{2\pi i}=1}$

La suma final, arriba, es precisamente una convolución cíclica de las dos secuencias a _q y b _q (de longitud N –1, porque ) definida por: ${\displaystyle q\in {}\{0,\dots ,N-2\}}$

${\displaystyle a_{q}=x_{g^{q}}}$

${\displaystyle b_{q}=e^{-{\frac {2\pi i}{N}}g^{-q}}.}$

Evaluando la convolución

Dado que N –1 es compuesto, esta convolución se puede realizar directamente mediante el teorema de convolución y algoritmos FFT más convencionales. Sin embargo, esto puede no ser eficiente si N –1 tiene factores primos grandes, lo que requiere el uso recursivo del algoritmo de Rader. En cambio, se puede calcular una convolución cíclica de longitud ( N –1) exactamente rellenándola con ceros hasta una longitud de al menos 2( N –1)–1, digamos a una potencia de dos , que luego se puede evaluar en O ( N log N ) tiempo sin la aplicación recursiva del algoritmo de Rader.

Este algoritmo, entonces, requiere O( N ) adiciones más O( N log N ) tiempo para la convolución. En la práctica, las adiciones O( N ) a menudo se pueden realizar absorbiendo las adiciones en la convolución: si la convolución se realiza mediante un par de FFT, entonces la suma de x _n viene dada por la salida DC (0ª) de la FFT de a _q más x ₀ , y x ₀ se pueden sumar a todas las salidas sumándolo al término DC de la convolución antes de la FFT inversa. Aún así, este algoritmo requiere intrínsecamente más operaciones que las FFT de tamaños compuestos cercanos y, en la práctica, normalmente tarda entre 3 y 10 veces más.

Si el algoritmo de Rader se realiza utilizando FFT de tamaño N –1 para calcular la convolución, en lugar de relleno con ceros como se mencionó anteriormente, la eficiencia depende en gran medida de N y del número de veces que el algoritmo de Rader debe aplicarse de forma recursiva. El peor caso sería si N –1 fuera 2 N ₂ donde N ₂ es primo, con N ₂ –1 = 2 N ₃ donde N ₃ es primo, y así sucesivamente. En tales casos, suponiendo que la cadena de números primos se extendiera hasta algún valor acotado, la aplicación recursiva del algoritmo de Rader en realidad requeriría tiempo O( N ² ) ^{[ dudoso – discutir ]} . Tales N _j se denominan primos de Sophie Germain , y dicha secuencia de ellos se denomina cadena de Cunningham de primer tipo. Sin embargo, se observa que las longitudes de las cadenas de Cunningham crecen más lentamente que log ₂ ( N ), por lo que el algoritmo de Rader aplicado de esta manera probablemente no sea O ( N ² ), aunque posiblemente sea peor que O ( N log N ) para los peores casos. Afortunadamente, se puede lograr una garantía de complejidad O ( N log N ) mediante relleno cero.

Referencias

1. CM Rader, "Transformaciones discretas de Fourier cuando el número de muestras de datos es primo", Proc. IEEE 56, 1107–1108 (1968).
2. S. Chu y C. Burrus, "Un algoritmo FTT [ sic ] de factor primo que utiliza aritmética distribuida", IEEE Transactions on Acoustics, Speech, and Signal Processing 30 (2), 217–227 (1982).
3. Matteo Frigo y Steven G. Johnson , "El diseño y la implementación de FFTW3", Actas del IEEE 93 (2), 216–231 (2005).
4. S. Winograd, "Sobre el cálculo de la transformada discreta de Fourier", Proc. Academia Nacional de Ciencias de EE. UU. , 73 (4), 1005–1006 (1976).
5. S. Winograd, "Sobre la computación de la transformada discreta de Fourier", Matemáticas de la Computación , 32 (141), 175-199 (1978).
6. R. Tolimieri, M. An y C.Lu, Algoritmos para convolución y transformada de Fourier discreta , Springer-Verlag, 2ª ed., 1997.
7. Donald E. Knuth, El arte de la programación informática, vol. 2: Algoritmos Seminuméricos , 3ª edición, sección 4.5.4, p. 391 (Addison-Wesley, 1998).