Русская Википедия:Алгоритм Кули — Тьюки
Алгоритм Ку́ли — Тью́ки — наиболее часто используемый алгоритм вычисления быстрого преобразования Фурье. Алгоритм позволяет выразить дискретное преобразование Фурье длины, равной произвольному составному числу <math>N</math>, через определённое количество преобразований меньшей длины с помощью рекурсии, понижая таким образом сложность вычислений до <math>O(N \log N)</math> для гладких <math>N</math>. Назван в честь Шаблон:Нп3 и Дж. Тьюки.
Основной алгоритм
Для произвольного натурального числа <math>N</math> дискретным преобразованием Фурье комплексного вектора <math>x(j)</math>, где <math>j = 0, \ldots, N-1</math>, называется комплексный вектор <math>X(k)</math>, где <math>k = 0, \ldots, N-1</math>, компоненты которого задаются формулой
- <math> X(k) = \sum\limits_{j=0}^{N-1} x(j) \omega^{kj}, \; k = 0, \ldots, N-1, </math>
где <math display="inline">\omega = e^{-\frac{2\pi i}{N}}</math>.
Для построения алгоритма делается предположение, что <math>N = N_1 N_2</math> для некоторых натуральных <math>N_1</math> и <math>N_2</math> и производится следующая замена индексовШаблон:Sfn:
- <math> \begin{array}{lll}
&j = j_1 + N_1 j_2, \ &j_1 = 0, \ldots, N_1-1, \ &j_2 = 0, \ldots, N_2-1, \\
&k = k_2 + N_2 k_1, \ &k_1 = 0, \ldots, N_1-1, \ &k_2 = 0, \ldots, N_2-1, \end{array} </math>
в результате которой получается
- <math> X(k_2 + N_2 k_1) = \sum\limits_{j_1 = 0}^{N_1-1} \sum\limits_{j_2 = 0}^{N_2-1} x(j_1 + N_1 j_2) \omega^{(j_1 + N_1 j_2)(k_2 + N_2 k_1)}. </math>
Далее векторы входных и выходных данных преобразуются в двумерные массивы <math>x'</math> и <math>X'</math>, задающиеся равенствами
- <math> \begin{array}{lll}
&x'(j_1,j_2) = x(j_1 + N_1 j_2), \ &j_1 = 0, \ldots, N_1-1, \ &j_2 = 0, \ldots, N_2-1, \\
&X'(k_1,k_2) = X(k_2 + N_2 k_1), \ &k_1 = 0, \ldots, N_1-1, \ &k_2 = 0, \ldots, N_2-1. \end{array} </math>
Стоит заметить, что компоненты <math>X</math> упорядочены не так, как компоненты <math>x</math>.
Далее пусть <math>\omega^{N_1} = \gamma</math> и <math>\omega^{N_2} = \beta</math>. Очевидно, что <math>\omega^{N_1 N_2 j_2 k_1} = 1</math>. Тогда в терминах двумерных переменных формула преобразуется к видуШаблон:Sfn
- <math> X'(k_1,k_2) = \sum\limits_{j_1=0}^{N_1-1} \beta^{j_1 k_1} \left( \omega^{j_1 k_2} \sum\limits_{j_2=0}^{N_2-1} x'(j_1,j_2) \gamma^{j_2 k_2} \right). </math>
Таким образом, вычисление преобразования длины <math>N = N_1 N_2</math> сводится к выполнению следующих действий:
- Вычисление <math>N_1</math> преобразований длины <math>N_2</math>.
- Умножение на комплексные «поворачивающие» множители.
- Вычисление <math>N_2</math> преобразований длины <math>N_1</math>.
При этом вместо <math>N^2</math> комплексных сложений и <math>N^2</math> комплексных умножений исходной формулы итоговая схема содержит не более <math>N(N_1 + N_2 - 2)</math> комплексных сложений и <math>N(N_1 + N_2 + 1)</math> комплексных умноженийШаблон:Sfn.
Каждое из преобразований длины <math>N_1</math> или <math>N_2</math> можно вычислять с помощью различных быстрых алгоритмов, в частности, снова по вышеприведённой схеме. В этом случае преобразование длины <math display="inline">N = \prod\limits_{i=1}^n N_i</math> может быть представлено в форме, требующей выполнения <math display="inline">N \sum\limits_{i=1}^n N_i</math> комплексных умноженийШаблон:Sfn.
Алгоритм по основанию два
Во многих приложениях длина преобразования равна степени двойки: <math>N = 2^m</math>. Тогда в вышеприведённой схеме возможны варианты <math>N_1 = 2</math> или <math>N_2 = 2</math>. В этом случае говорят об алгоритме Кули — Тьюки по основанию дваШаблон:Sfn (Шаблон:Lang-en).
Если <math>N_1 = 2</math>, то говорят об алгоритме Кули — Тьюки с прореживанием по времениШаблон:Sfn. В этом случае уравнения преобразуются к виду
- <math> \begin{array}{rcl}
X(k) &= &\sum\limits_{j=0}^{N/2 - 1} x(2j) \omega^{2kj} + \omega^k \sum\limits_{j=0}^{N/2 - 1} x(2j+1) \omega^{2kj} \\ X(k+N/2) &= &\sum\limits_{j=0}^{N/2 - 1} x(2j) \omega^{2kj} - \omega^k \sum\limits_{j=0}^{N/2 - 1} x(2j+1) \omega^{2kj}, \ k = 0,\ldots,N/2-1. \end{array} </math>
Если ввести обозначения <math>E(k) = \sum\limits_{j=0}^{N/2 - 1} x(2j) \omega^{2kj}</math> и <math>O(k) = \sum\limits_{j=0}^{N/2 - 1} x(2j+1) \omega^{2kj}</math>, то уравнения можно переписать как
- <math> \begin{array}{rcl}
X(k) &= &E(k) + \omega^k O(k) \\ X(k+N/2) &= &E(k) - \omega^k O(k), \ k = 0,\ldots,N/2-1. \end{array} </math> Такая операция обычно называется «бабочкой».
Данная процедура может быть применена к входному вектору рекурсивно. На каждом шаге <math>N</math>-точечное преобразование разбивается на два <math>(N/2)</math>-точечных преобразования, которые, в свою очередь, разбиваются таким же образом до тех пор, пока длина преобразования не станет равна единице. Затем происходит обратный ход, на каждом шаге длины результатов преобразований удваиваются с помощью «бабочек». При такой реализации выполняется <math>(N/2) \log_2 N</math> комплексных умножений и <math>N \log_2 N</math> комплексных сложений.
Этот процесс является примером применения методики «разделяй и властвуй». При этом во многих реализациях прямой рекурсии избегают и вместо неё дерево вычислений проходится в порядке поиска в ширину.
Если <math>N_2 = 2</math>, то говорят об алгоритме Кули — Тьюки с прореживанием по частотеШаблон:Sfn. В этом случае уравнения преобразуются к виду
- <math> \begin{array}{rcl}
X(2k) &= &\sum\limits_{j=0}^{N/2 - 1} \bigg( x(j) + x(j+N/2) \bigg) \omega^{2kj}, \\ X(2k+1) &= &\sum\limits_{j=0}^{N/2 - 1} \bigg( x(j) - x(j+N/2) \bigg) \omega^{j(2k+1)}, \ k = 0,\ldots,N/2-1. \end{array} </math>
Алгоритм Рейдера — Бреннера
Пусть
- <math> a(k) = \left\lbrace \begin{array}{rl} 0, &k=0, \\ \dfrac{x(k) - x(k+N/2)}{2i \sin\left( 2\pi k/N \right)}, &k=1,\ldots,N/2-1 \end{array} \right. </math>
и пусть
- <math> A(k) = \sum\limits_{j=0}^{N/2-1} a(j) \omega^{2kj}, \ k = 0,\ldots,N/2-1. </math>
С использованием формул алгоритма с прореживанием по частоте нетрудно убедиться, что выполняется следующее соотношение:
- <math> X(2k+1) = A(k) - A(k+1) + (x(0) - x(N/2)), \ k = 0,\ldots,N/2-1. </math>
Такая модификация алгоритма по основанию два называется алгоритмом Рейдера — Бреннера. Она позволяет уменьшить вычислительную сложность за счёт более простых умножений на чисто мнимые константы <math>1/(2i \sin\left( 2\pi k/N \right))</math>Шаблон:Sfn. Аналогичные формулы можно получить с использованием вещественных констант <math>1/(2 \cos\left( 2\pi k/N \right))</math>Шаблон:Sfn.
История
Алгоритм и его рекурсивная реализация были изобретены около 1805 года К. Гауссом при интерполировании траекторий астероидов Паллада и Юнона[1]. Тогда открытие не получило широкого распространения и было опубликовано лишь после смерти учёного на новой латыни. Результат Гаусса несколько раз переоткрывался в различных формах в течение последующих 150 лет и стал популярным после публикации в 1965 году статьи Шаблон:Нп3 из IBM и Дж. Тьюки из Принстона, в которой алгоритм был в очередной раз переоткрыт, а также описывалась удобная реализация для ЭВМ[2].
Тот факт, что первооткрывателем алгоритма является Гаусс, был обнаружен лишь через несколько лет после публикации Кули и Тьюки. В своей статье они ссылались только на работу И. Дж. Гуда, в которой был описан алгоритм Гуда — Томаса.
Выражение преобразования Фурье длины <math>N</math> через два преобразования длины <math>N/2</math> иногда называют леммой Шаблон:Нп3 — Ланцоша, так как оно было получено этими двумя авторами в 1942 году[3].
Примечания
Литература
См. также
