Этот раздел посвящен периодическим функциям. Они рассматриваются для случая периода равного $ 2\pi_{} $, к которому можно свести общий случай соответствующей заменой переменной («растяжением» шкалы).

Тригонометрическим полиномом от переменной $ x_{} $ называется выражение вида $$ \begin{matrix} g_n(x)=a_0 & + & a_1 \cos x + a_2 \cos 2x+\dots + a_n \cos nx + \\ \ &+&b_1 \sin x + b_2 \sin 2x+\dots + b_n \sin nx = \end{matrix} $$ $$ =a_0 + \sum_{k=1}^n (a_k \cos\, kx+ b_k \sin \, kx) \ . $$ при фиксированных постоянных $ a_0,a_1,a_2,\dots,a_n,b_1,\dots,b_n $ — вещественных или комплексных. Если $ a_{n} $ или $ b_n $ отличны от нуля, то говорят, что полином $ g_n(x) $ имеет порядок $ n_{} $. В этом случае числа $ a_0,\{a_j,b_j\}_{j=1}^n $ называются коэффициентами тригонометрического полинома, а числа $ a_{n} $ и $ b_n $ — старшими коэффициентами. Если все числа $ \{b_j\}_{j=1}^n $ равны нулю, то о полиноме $ g_n(x) $ говорят как о полиноме по косинусам, если все числа $ \{a_j\}_{j=0}^n $ равны нулю — то полином называют полиномом по синусам.

Чтобы отличать тригонометрический полином от «обычного» $ a_0+a_1x+\dots + a_nx^n $ мы будем называть последний степенным или алгебраическим.

Пример. Функции $ \cos^n x $ и $ \sin^n x $ являются тригонометрическими полиномами $ n_{} $-го порядка; первый из них является полиномом по косинусам:

$$ \cos^3 x= \frac{3}{4}\cos x + \frac{1}{4}\cos 3\,x \ ; $$ второй же — в зависимости от четности $ n_{} $ — либо полиномом по косинусам, либо — по синусам: $$ \sin^4 x=\frac{3}{8}-\frac{1}{2}\cos 2\,x+\frac{1}{8}\cos 4\, x \ , \ \sin^5 x=\frac{5}{8}\sin x-\frac{5}{16} \sin 3\,x+\frac{1}{16}\sin 5\,x \ . $$ Для общего случая доказано ☞ ЗДЕСЬ.

Теорема. Произвольный тригонометрический полином $ n_{} $-го порядка может быть представлен в виде

$$ g_n(x) \equiv e^{-\mathbf i n x} G_{2n} (e^{\mathbf i x}) \ , $$ где $$ G_{2n}(z)=\sum_{k=0}^{2n} u_kz^k $$ — некоторый алгебраический полином степени $ 2n_{} $.

Доказательство. Перейдем к комплексной переменной: в предположении, что $ x_{} $ переменная вещественная , обозначим $$ z= \cos x + \mathbf i \sin x \ , $$ Тогда, на основании свойств операций над комплексными числами (в том числе формулы Муавра ), имеем $$ z^{-1} = \cos x - \mathbf i \sin x \ , z^j= \cos j\,x + \mathbf i \sin j\,x \ , z^{-j}= \cos j\,x - \mathbf i \sin j\,x \quad npu \ \in\{1,\dots,n\} \ , $$ и, следовательно, $$ \cos j\,x = \frac{1}{2}\left(z^j+z^{-j}\right),\quad \sin j\, x=\frac{1}{2\mathbf i} \left(z^j-z^{-j}\right) \ . $$ Подставим эти выражения в $ g_n(x) $ и домножим получившееся на $ z^n $: $$G_{2n}(z)=\underbrace{\frac{a_n+\mathbf i b_n}{2}}_{u_0}+\underbrace{\frac{a_{n-1}+\mathbf i b_{n-1}}{2}}_{u_1}z+\dots+\underbrace{\frac{a_1+\mathbf i b_1}{2}}_{u_{n-1}}z^{n-1}+ \underbrace{a_0}_{u_n}z^n+ $$ $$ +\underbrace{\frac{a_1-\mathbf i b_1}{2}}_{u_{n+1}}z^{n+1}+\dots+\underbrace{\frac{a_{n-1}-\mathbf i b_{n-1}}{2}}_{u_{2n-1}}z^{2n-1} + \underbrace{\frac{a_n-\mathbf i b_n}{2}}_{u_{2n}}z^{2n} \ . $$ Если переменная $ x $ предполагается комплексной, то придется использовать экспоненту от комплексного числа: $ z=e^{\mathbf i x} $ и представления через нее $ \sin x $ и $ \cos x $ (см. ☞ЗДЕСЬ ). ♦

Если все коэффициенты $ a_0,\{a_j,b_j\}_{j=1}^n $ тригонометрического полинома $ g_n(x) $ вещественны, то полином $ G_{2n}(z) $ из теоремы удовлетворяет тождеству:

$$ G_{2n}(z) \equiv z^{2n} \overline{G_{2n}}(z^{-1}) \ . $$ Иными словами, для коэффициентов такого полинома $$ G_{2n}(z)=u_0+u_1z+u_2z^2+\dots+u_{2n-1}z^{2n-1} + u_{2n}z^{2n} $$ справедлива «симметрия с сопряжением» набора его коэффициентов $ (u_0,u_1,\dots,u_{2n-1},u_{2n}) $ относительно середины этого набора: $$ u_0=\overline{u_{2n}},u_1=\overline{u_{2n-1}},\dots, u_j=\overline{u_{2n-j}},\dots $$

Корнем или нулем тригонометрического полинома $ g_n(x) $ называется комплексное число $ \lambda $ такое, что $ g_n(\lambda)=0 $. Договариваются не различать между собой корни различающиеся на целое кратное $ 2\pi_{} $ (т. е. сравнимые по модулю $ 2\pi_{} $). Определение кратного корня тригонометрического полинома вводится аналогично случаю полинома алгебраического.

Следующая теорема является аналогом основной теоремы высшей алгебры для степенных полиномов.

Теорема. Тригонометрический полином $ n_{} $-го порядка, у которого старшие коэффициенты $ a_{n} $ и $ b_n $ удовлетворяют условиям $ a_n+ \mathbf i b_n \ne 0, a_n- \mathbf i b_n \ne 0 $, имеет в точности $ 2n $ комплексных корней с учетом их кратностей.

Доказательство ☞ ЗДЕСЬ. ♦

Если числа $ \{\lambda_1,\dots, \lambda_{2n} \} $ при $ \{ \mathfrak{Re} (\lambda_k)\}_{k=1}^{2n} \subset [0,2\pi[ $ --- корни тригонометрического полинома $ g_n(x) $, у которого $ a_n+ \mathbf i b_n \ne 0, a_n- \mathbf i b_n \ne 0 $, то этот полином можно представить в виде

$$ g_n(x) \equiv A \prod_{k=1}^{2n} \sin \frac{x-\lambda_k}{2} \quad \mbox{ при } \quad A=\pm 2^{2n-1}\sqrt{a_n^2+b_n^2} \ . $$

Имеет место тождество:

$$ \sin \, nx \equiv 2^{n-1} \prod_{j=0}^{n-1} \sin\left(x-\frac{\pi j}{n} \right) \equiv 2^{n-1} \sin \, x \sin \left(x-\frac{\pi}{n}\right) \sin \left(x-\frac{2\pi}{n} \right)\times \dots \times \sin \left(x-\frac{(n-1)\pi}{n} \right) \ . $$

Доказательство обоих следствий ☞ ЗДЕСЬ. ♦

Задача. Построить тригонометрический полином порядка $ n_{} $, принимающий значения по следующей таблице $$ \begin{array}{c|ccccc} x & x_1 & x_2 & \dots & x_{2\,n+1} \\ \hline y & y_1 & y_2 &\dots & y_{2\,n+1} \end{array} $$ т.е. $ g_n(x_j) = y_j $ при $ j\in \{1,\dots,2\,n+1\} $. Здесь узлы интерполяции $ \{x_j\}_{j=1}^{2\,n+1} $ предполагаются вещественными, различными и расположенными в интервале $ [0, 2\pi [ $.

Теорема. Функция

$$ \begin{matrix} g_n(x)&=&y_1\frac{\displaystyle \sin\frac{x-x_2}{2}\times\dots \times\sin \frac{x-x_{2n+1}}{2}}{\displaystyle \sin \frac{x_1-x_2}{2}\times \dots \times \sin \frac{x_1-x_{2n+1}}{2}} + \\ &+&y_2\frac{\displaystyle \sin\frac{x-x_1}{2}\sin\frac{x-x_3}{2}\times \dots \times \sin\frac{x-x_{2n+1}}{2}}{\displaystyle \sin\frac{x_2-x_1}{2}\sin\frac{x_2-x_3}{2}\times \dots \times \sin\frac{x_2-x_{2n+1}}{2}} + \\ & + & \dots + \\ &+& y_{2n+1}\frac{\displaystyle \sin \frac{x-x_1}{2} \times \dots \times \sin \frac{x-x_{2n}}{2}}{ \displaystyle \sin \frac{x_{2n+1}-x_1}{2} \times \dots \times \sin \frac{x_{2n+1}-x_{2n}}{2}} \end{matrix} $$ удовлетворяет условиям $ g_n(x_j)=y_j $ при $ j\in \{1,\dots,2n+1\} $. Эта функция является тригонометрическим полиномом порядка не выше $ n_{} $. При таком ограничении на порядок, полином, удовлетворяющий заданным интерполяционным условиям, определяется единственным образом.

Доказательство факта, что функция $ g_{n}(x) $ удовлетворяет интерполяционной таблице очевидно: имеем полный аналог интерполяционного степенного полинома в форме Лагранжа. Образно говоря, каждое слагаемое $ y_jW_j(x) $ в этой сумме «отвечает» исключительно только за свой узел интерполяции (т.е. обеспечивает в нем нужное значение) — и при этом «не портит» остальные узлы (обращается в них в нуль).

Теперь нужно показать, что функция $ g_{n}(x) $ является именно тригонометрическим полиномом, т.е. допускает представление вида $ \displaystyle a_0 + \sum_{k=1}^n (a_k \cos \, kx + b_k \sin \, kx) $. Это достигается путем преобразования произведений синусов, стоящих в числителях слагаемых. В каждом таком произведении содержится четное число сомножителей, объединяем их попарно и используем формулы приведения: $$ \sin \frac{x-x_1}{2} \sin \frac{x-x_2}{2} \equiv \frac{1}{2} \left( \cos \frac{x_2-x_1}{2} - \cos(x-\frac{x_1+x_2}{2} \right)\equiv $$ $$ \equiv \frac{1}{2} \cos \frac{x_2-x_1}{2} - \frac{1}{2} \cos \frac{x_1+x_2}{2} \cos\, x - \frac{1}{2} \sin \frac{x_1+x_2}{2} \sin \,x \ . $$ В результате, мы избавляемся от половинного аргумента под синусами. Дальнейшие преобразования аналогичны. Окончательная формула ☞ ЗДЕСЬ.

Единственность этого полинома при условии ограничения на порядок следует из теоремы предыдущего пункта: если $ h_n(x) $ — другой полином, удовлетворяющий условию теоремы, то полином $ g_n(x)-h_n(x) $ порядка $ \le 2\,n $ имеет $ \ge (2n+1) $-го корней на интервале $ [0,2\, \pi[ $, что невозможно. ♦

Итак, теоретическая база решения задачи тригонометрической интерполяции построена. Осталась практическая часть: интересует представление полинома из теоремы в каноническом виде $$ g_n(x)\equiv a_0 + \sum_{k=1}^n (a_k \cos\, kx+ b_k \sin \, kx) \ , $$ т.е. его коэффициенты. Стандартный способ их поиска — метод неопределенных коэффициентов — приводит к системе из $ (2n+1) $-го линейного уравнения $$ \left(\begin{array}{llllllll} 1 & \cos x_1 & \sin x_1 & \cos \, 2\, x_1 & \sin \, 2\, x_1 & \dots & \cos \, n\, x_1 & \sin \, n\, x_1 \\ 1 & \cos x_2 & \sin x_2 & \cos \, 2\, x_2 & \sin \, 2\, x_2 & \dots & \cos \, n\, x_2 & \sin \, n\, x_2 \\ \dots & & & & & & & \dots \\ 1 & \cos x_{2n+1} & \sin x_{2n+1} & \cos \, 2\, x_{2n+1} & \sin \, 2\, x_{2n+1} & \dots & \cos \, n\, x_{2n+1} & \sin \, n\, x_{2n+1} \end{array} \right) \left(\begin{array}{l} a_0 \\ a_1 \\ \vdots \\ b_n \end{array} \right)= \left(\begin{array}{l} y_1 \\ y_2 \\ \vdots \\ y_{2n+1} \end{array} \right) \ . $$ относительно $ (2n+1) $-го коэффициента $ a_0,a_1,b_1,a_2,b_2,\dots,a_n,b_n $. Согласно теореме, эта система должна иметь решение относительно столбца коэффициентов, причем решение единственное. Из этого факта, в частности, следует, что определитель матрицы этой системы должнен быть отличен от нуля (см. ☞ формулы Крамера ). Так оно и есть.

Теорема. Определитель матрицы системы равен

$$ 2^{2n^2} \prod_{0\le k < j \le 2n} \sin \frac{x_k-x_j}{2} \ . $$

Теперь обсудим возможные упрощения представления тригонометрического полинома из теоремы, которые получатся за счет удачного выбора узлов интерполяции.

Теорема. Если в качестве узлов интерполяции взять равноотстоящие

$$ x_j=\frac{2(j-1)\pi}{2n+1} \quad npu \quad j\in \{1,2,\dots,2n+1\} \ , $$ то коэффициенты интерполяционного тригонометрического полинома находятся по формулам $$ a_0=\frac{1}{2n+1}\sum_{j=1}^{2n+1}y_j, $$ $$ a_k= \frac{2}{2n+1} \sum_{j=1}^{2n+1} y_j \cos k\,x_j ,\quad b_k= \frac{2}{2n+1} \sum_{j=1}^{2n+1} y_j \sin k\,x_j \quad npu \ k\in \{1,\dots,n \} . $$

Доказательство ☞ ЗДЕСЬ.

Формулы последней теоремы остаются справедливыми и в случае «сдвига» системы узлов интерполяции на величину, кратную $ \pi_{} $, например, при выборе системы узлов $$ x_j=\frac{2\pi j}{2n+1} \quad npu \quad j\in \{-n,-n+1,\dots,-1,0,1,2,\dots,n\} $$ на интервале $ [-\pi,\pi[ $: $$ a_0=\frac{1}{2n+1}\sum_{j=-n}^{n}y_j,\quad a_k= \frac{2}{2n+1} \sum_{j=-n}^n y_j \cos k\,x_j ,\quad b_k= \frac{2}{2n+1} \sum_{j=-n}^n y_j \sin k\,x_j \quad npu \ k\in \{1,\dots,n \} . $$

Задачу тригонометрической интерполяции можно интерпретировать как построение некоторого отображения. «На вход» этого отображения подается набор значений $ \{y_1,\dots,y_{2n+1}\} $ интерполируемой функции, «на выходе» должны получить коэффициенты $ \{a_0,a_1,b_1,\dots,a_n,b_n\} $ интерполирующего тригонометрического полинома. Обсудим свойства этого отображения. Первым делом условимся, что узлы интерполяции $ \{x_1,\dots,x_{2n+1} \} $ фиксированы, а значения в них интерполируемых функций могут меняться от эксперимента к эксперименту. Запишем каждый набор значений в виде вектора, и в том же стиле запишем коэффициенты тригонометрического полинома. Для определенности будем рассматривать эти векторы как столбцы в пространстве $ \mathbb C^n $ (т.е. и значения функции и коэффициенты тригонометрического полинома могут быть комплексными), получаем некоторое соответствие $$ \left( \begin{array}{l} y_1\\ y_2 \\ \vdots \\ y_{2n+1} \end{array} \right) \mapsto \left( \begin{array}{l} a_0\\ a_1 \\ b_1 \\ \vdots \\ a_{n} \\ b_n \end{array} \right) \ . $$

Итак, вектор пространства $ \mathbb C^n $ отображается в вектор того же пространства. Это соответствие является именно функцией, причем функцией взаимно-однозначной: поскольку, по доказанному выше, каждому набору значений соответствует единственный набор коэффициентов соответствующего тригонометрического полинома. Таким образом, мы имеем дело с преобразованием пространства $ \mathbb C^n $.

Теорема. Это преобразование линейно.

Доказательство этого факта разбивается на проверку двух свойств линейности преобразования, которые мы переведем «на язык интерполяции». Умножение всех значений $ y_j $ на одну и ту же константу $ \alpha_{} $ приводит к умножению на ту же константу интерполяционного полинома $$ \begin{array}{c} \begin{array}{c|ccccc} x & x_1 & x_2 & \dots & x_{2\,n+1} \\ \hline y & y_1 & y_2 &\dots & y_{2\,n+1} \end{array} \\ \downarrow \\ g_n(x) \end{array} \qquad \qquad \begin{array}{c} \begin{array}{c|ccccc} x & x_1 & x_2 & \dots & x_{2\,n+1} \\ \hline y & \alpha y_1 & \alpha y_2 &\dots & \alpha y_{2\,n+1} \end{array} \\ \downarrow \\ \alpha f_n(x) \end{array} $$ а сложение двух интерполяционных таблиц приводит к сложению соответствующих интерполяционных полиномов: $$ \begin{array}{c} \begin{array}{c|ccccc} x & x_1 & x_2 & \dots & x_{2\,n+1} \\ \hline y & y_1 & y_2 &\dots & y_{2\,n+1} \end{array} \\ \downarrow \\ g_n(x) \\ \searrow \end{array} \qquad \qquad \begin{array}{c} \begin{array}{c|ccccc} x & x_1 & x_2 & \dots & x_{2\,n+1} \\ \hline y & \tilde y_1 & \tilde y_2 &\dots & \tilde y_{2\,n+1} \end{array} \\ \downarrow \\ \tilde g_n(x) \\ \swarrow \end{array} $$ $$ \begin{array}{c} \begin{array}{c|ccccc} x & x_1 & x_2 & \dots & x_{2\,n+1} \\ \hline y & y_1+ \tilde y_1 & y_2 + \tilde y_2 &\dots & y_{2\,n+1} + \tilde y_{2\,n+1} \end{array} \\ \downarrow \\ g_n(x)+\tilde g_n(x) \end{array} $$ Оба утверждения следуют непосредственно из доказательства теоремы о существовании тригонометрического полинома: $$ g_n(x) = \sum_{j=1}^{2n+1} y_j W_j(x) \, $$ где тригонометрические полиномы $ \{ W_j(x) \}_{j=1}^{2n+1} $ не зависят от $ \{y_1,y_2,\dots,y_{2n+1}\} $. ♦

Теперь наша задача — наиболее удобным способом формализовать это преобразование и выписать его матрицу. В общем случае, решение этой задачи громоздко и мы будем заниматься случаем равноотстоящих узлов $$ x_j=\frac{2(j-1)\pi}{2n+1} \quad npu \quad j\in \{1,2,\dots,2n+1\} \ , $$ для которого мы получили компактные формулы в предыдущем пункте: $$ a_0=\frac{1}{2n+1}\sum_{j=1}^{2n+1}y_j, $$ $$ a_k= \frac{2}{2n+1} \sum_{j=1}^{2n+1} y_j \cos k\,x_j ,\quad b_k= \frac{2}{2n+1} \sum_{j=1}^{2n+1} y_j \sin k\,x_j \ npu \ k\in \{1,\dots,n \} . $$ Обратим внимание на симметрию вхождения косинусов и синусов в две последние формулы. Вспоминая тригонометрическую форму комплексного числа, хочется объединить эти формулы в одну, с использованием мнимой единицы $ \mathbf i $: $$ a_k+ \mathbf i b_k = \frac{2}{2n+1} \sum_{j=1}^{2n+1} y_j \left( \cos k\,x_j + \mathbf i \sin k\,x_j \right) = $$ Теперь еще одно приятное воспоминание о замечательной формуле Муавра позволит нам привести выражение к виду $$ =\frac{2}{2n+1} \sum_{j=1}^{2n+1} y_j \left( \cos \, x_j + \mathbf i \sin \, x_j \right)^k \ . $$ Но это ещё не всё, мы можем продолжить и дальше, если обратим внимание на величины $ x_{j} $: $$ \cos \, x_j + \mathbf i \sin \, x_j = \cos \frac{2(j-1)\pi}{2n+1} + \mathbf i \sin \frac{2(j-1)\pi}{2n+1} $$ и в правой части формулы стоит объект, известный как корень из единицы степени $ 2n+1 $. Если обозначить его традиционным для теории комплексных чисел обозначением $$ \varepsilon_{j-1} = \cos \frac{2\pi (j-1)}{2n+1} + \mathbf i \sin \frac{2\pi (j-1)}{2n+1} \quad npu \quad j\in \{1,2,\dots,2n+1 \} \, $$ то мы получаем окончательную формулу $$ a_k+ \mathbf i b_k = \frac{2}{2n+1} \sum_{j=1}^{2n+1} y_j \varepsilon_{j-1}^k \quad npu \quad k\in \{1,2,\dots,n \} \ . $$

Теперь объединим формулы для $ a_{k} $ и $ b_{k} $ по-другому и преобразуем выражение $ a_k- \mathbf i b_k $: $$ a_k- \mathbf i b_k = \frac{2}{2n+1} \sum_{j=1}^{2n+1} y_j \left( \cos k\,x_j - \mathbf i \sin k\,x_j \right) = \frac{2}{2n+1} \sum_{j=1}^{2n+1} y_j \left( \cos (-k\,x_j) + \mathbf i \sin (-k\,x_j) \right)= $$ $$ = \frac{2}{2n+1} \sum_{j=1}^{2n+1} y_j \left( \cos (-\,x_j) + \mathbf i \sin (-\,x_j) \right)^k=\frac{2}{2n+1} \sum_{j=1}^{2n+1} y_j \varepsilon_{-j+1}^k \quad npu \quad k\in \{1,2,\dots,n \} \ . $$ Здесь комплексное число $$ \varepsilon_{-j+1} =\cos \frac{2(j-1)\pi}{2n+1} - \mathbf i \sin \frac{2(j-1)\pi}{2n+1} $$ также является корнем степени $ 2n+1 $ из единицы.

Числа $ (a_k+ \mathbf i b_k)/2 $ и $ (a_k- \mathbf i b_k)/2 $ уже встречались нам ВЫШЕ и являлись коэффициентами алгебраического полинома $$G_{2n}(z)=\underbrace{\frac{a_n+\mathbf i b_n}{2}}_{u_0}+\underbrace{\frac{a_{n-1}+\mathbf i b_{n-1}}{2}}_{u_1}z+\dots+\underbrace{\frac{a_1+\mathbf i b_1}{2}}_{u_{n-1}}z^{n-1}+ \underbrace{a_0}_{u_n}z^n+ $$ $$ +\underbrace{\frac{a_1-\mathbf i b_1}{2}}_{u_{n+1}}z^{n+1}+\dots+\underbrace{\frac{a_{n-1}-\mathbf i b_{n-1}}{2}}_{u_{2n-1}}z^{2n-1} + \underbrace{\frac{a_n-\mathbf i b_n}{2}}_{u_{2n}}z^{2n} \ , $$ связанного с тригонометрическим полиномом $ g_n(x) $ соотношением $$ g_n(x) \equiv e^{-\mathbf i n x} G_{2n} (e^{\mathbf i x}) \ . $$ Оказывается, что наша задача тригонометрической интерполяции удобнее всего записывается именно в терминах коэффициентов полинома $ G_{2n}(z) $. Собственно, мы почти дошли до цели: мы собираемся рассматривать отображение вектора значений интерполяционной таблицы в вектор коэффициентов полинома $ G_{2n}(z) $. Осталось сделать последнее упрощение: с не очень понятными пока последствиями мы переставим местами эти коэффициенты: $$ \left( \begin{array}{l} y_1\\ y_2 \\ \vdots \\ y_{2n+1} \end{array} \right) \mapsto \left( \begin{array}{l} u_n\\ u_{n+1} \\ \vdots \\ u_{2n} \\ u_0 \\ u_1 \\ \vdots \\ u_{n-1} \end{array} \right) \ . $$ Теперь выпишем формулы, связывающие эти столбцы. Первые $ n_{} $ компонент второго столбца $$ u_n=a_0=(y_1+\dots+y_{2n+1})/(2n+1) $$ $$ 2 u_{n+1}=a_1-\mathbf i b_1=(y_1+y_2\varepsilon_{-1}+y_3\varepsilon_{-2}+\dots+y_{2n+1}\varepsilon_{-2n})/(2n+1) $$ $$ 2 u_{n+2}=a_2-\mathbf i b_2=(y_1+y_2\varepsilon_{-1}^2+y_3\varepsilon_{-2}^2+\dots+y_{2n+1}\varepsilon_{-2n}^2)/(2n+1) $$ $$ \dots $$ $$ 2 u_{2n}=a_n-\mathbf i b_n=(y_1+y_2\varepsilon_{-1}^n+y_3\varepsilon_{-2}^n+\dots+y_{2n+1}\varepsilon_{-2n}^n)/(2n+1) $$ Теперь рассмотрим следующую компоненту столбца коэффициентов $$ 2u_0=a_n+\mathbf i b_n=(y_1+y_2\varepsilon_{1}^n+y_3\varepsilon_{2}^n+\dots+y_{2n+1}\varepsilon_{2n}^n)/(2n+1) \ . $$ Заметим, что корни степени $ 2n+1 $ из единицы удовлетворяют легко проверяемым соотношениям $$ \varepsilon_{1}^n=\varepsilon_{-1}^{n+1}, \varepsilon_{2}^n=\varepsilon_{-2}^{n+1},\dots \varepsilon_{2n}^n=\varepsilon_{-2n}^{n+1} \ , $$ и, таким образом, последняя формула эквивалентна $$ 2u_0=a_n+\mathbf i b_n=(y_1+y_2\varepsilon_{-1}^{n+1}+y_3\varepsilon_{-2}^{n+1}+\dots+y_{2n+1}\varepsilon_{-2n}^{n+1})/(2n+1) \ . $$ Пойдем дальше аналогичными преобразованиями $$ 2u_1=a_{n+1}+\mathbf i b_{n+1}=(y_1+y_2\varepsilon_{-1}^{n+2}+y_3\varepsilon_{-2}^{n+2}+\dots+y_{2n+1}\varepsilon_{-2n}^{n+2})/(2n+1) \ , $$ $$ \dots $$ $$ 2u_{n-1}=a_{2n}+\mathbf i b_{2n}=(y_1+y_2\varepsilon_{-1}^{2n}+y_3\varepsilon_{-2}^{2n}+\dots+y_{2n+1}\varepsilon_{-2n}^{2n})/(2n+1) \ . $$ Теперь переписываем формулы в матричном виде $$ (2n+1)\left( \begin{array}{rl} & u_n \\ 2 & u_{n+1} \\ 2 & u_{n+2} \\ & \vdots \\ 2 & u_{2n} \\ 2 & u_0 \\ 2 & u_1 \\ & \vdots \\ 2 & u_{n-1} \end{array} \right) = \left( \begin{array}{lllll} 1 & 1 & 1 & \dots & 1 \\ 1 & \varepsilon_{-1} & \varepsilon_{-2} & \dots & \varepsilon_{-2n} \\ 1 & \varepsilon_{-1}^2 & \varepsilon_{-2}^2 & \dots & \varepsilon_{-2n}^2 \\ & \dots & & & \dots \\ 1 & \varepsilon_{-1}^n & \varepsilon_{-2}^n & \dots & \varepsilon_{-2n}^n \\ 1 & \varepsilon_{-1}^{n+1} & \varepsilon_{-2}^{n+1} & \dots & \varepsilon_{-2n}^{n+1} \\ & \dots & & & \dots \\ 1 & \varepsilon_{-1}^{2n} & \varepsilon_{-2}^{2n} & \dots & \varepsilon_{-2n}^{2n} \end{array} \right) \left( \begin{array}{l} y_1 \\ y_2 \\ \vdots \\ y_{2n+1} \end{array} \right) $$ и видим, что в правой части возникла матрица Вандермонда так характерная для задачи алгебраической интерполяции.

Пришло время давать самое главное определение.

Дискретное преобразование Фурье

Пусть $ N_{} $ — произвольное натуральное число. Говорят, что cистема (последовательность) чисел $ \{ \mathfrak X_0,\dots, \mathfrak X_{N-1} \} $ получается из системы (последовательности) чисел $ \{ \mathfrak x_0,\dots, \mathfrak x_{N-1} \} $ в результате дискретного преобразования Фурье (ДПФ)¹⁾, если $$ \left( \begin{array}{l} \mathfrak X_0 \\ \mathfrak X_1 \\ \vdots \\ \mathfrak X_{N-1} \end{array} \right) = \left( \begin{array}{lllll} 1 & 1 & 1 & \dots & 1 \\ 1 & \varepsilon_{-1} & \varepsilon_{-2} & \dots & \varepsilon_{-N+1} \\ 1 & \varepsilon_{-1}^2 & \varepsilon_{-2}^2 & \dots & \varepsilon_{-N+1}^2 \\ & \dots & & & \dots \\ 1 & \varepsilon_{-1}^{N-1} & \varepsilon_{-2}^{N-1} & \dots & \varepsilon_{-N+1}^{N-1} \\ \end{array} \right)_{N\times N} \left( \begin{array}{l} \mathfrak x_0 \\ \mathfrak x_1 \\ \vdots \\ \mathfrak x_{N-1} \end{array} \right) \quad npu \quad \varepsilon_{-j}= \cos \frac{2\pi j}{N} - \mathbf i \sin \frac{2\pi j}{N} \ . $$ или, в эквивалентном виде, $$ \mathfrak X_k=\sum_{j=0}^{N-1} \mathfrak x_j \varepsilon_{-j}^k \quad npu \quad k\in\{0,\dots, N-1\} \ . $$ Будем обозначать ДПФ буквой $ \mathcal F $: $$ \mathcal F(\mathfrak x_0,\dots, \mathfrak x_{N-1}) = (\mathfrak X_0,\dots, \mathfrak X_{N-1}) \ ; $$ оно, очевидно, является линейным преобразованием пространства $ \mathbb C^N $. Естественно считать наборы чисел $ (\mathfrak x_0,\dots, \mathfrak x_{N-1}) $ и $ (\mathfrak X_0,\dots, \mathfrak X_{N-1}) $ векторами этого пространства; но в ТЕОРИИ СИГНАЛОВ принято говорить о последовательности $ \{ \mathfrak x_j \}_{j=0}^{N-1} $ как о входной последовательности; а о последовательности $ \{ \mathfrak X_k \}_{k=0}^{N-1} $ — как о выходной последовательности. Еще раз подчеркнем смысл чисел $ \{\varepsilon_{-j}\}_{j=0}^{N-1} $: это — система всех различных корней $ N_{} $-й степени из $ 1_{} $ ; только в отличие от стандартного определения, они пронумерованы в обратном порядке (на единичной окружности комплексной плоскости движемся по часовой стрелке). Поскольку имеет место равенство $$ \varepsilon_{-j}=\varepsilon_{-1}^j = \varepsilon_{1}^{-j} \ , $$ то дискретное преобразование Фурье можно переписать еще в двух эквивалентных видах: $$ \mathfrak X_k=\sum_{j=0}^{N-1} \mathfrak x_j \varepsilon_{-1}^{jk}= \sum_{j=0}^{N-1} \mathfrak x_j \varepsilon_{1}^{-jk} \quad npu \quad k\in\{0,\dots, N-1\} \ . $$ Наконец, использование экспоненциального представления корня из единицы $ \varepsilon_{-j} =e^{-\mathbf i (2\pi j/N)} $ дает еще одну часто используемую формулу: $$ \mathfrak X_k=\sum_{j=0}^{N-1} \mathfrak x_j e^{-\mathbf i(2\pi jk/N)} \quad npu \quad k\in\{0,\dots, N-1\} \ . $$

Найти ДПФ набора $ (-1, 2,1,-3) $.

Решение. Здесь $ N_{}=4 $, а корни четвертой степени из единицы хорошо известны: $ \{\varepsilon_{-j}\}_{j=0}^3=\{1,-\mathbf i,-1, \mathbf i\} $. $$ \left( \begin{array}{rrrr} 1 & 1 & 1 & 1 \\ 1 & -\mathbf i & -1 & \mathbf i \\ 1 & -1 & 1 & -1 \\ 1 & \mathbf i & -1 & -\mathbf i \end{array} \right) \left( \begin{array}{r} -1 \\ 2 \\ 1 \\ -3 \end{array} \right) = \left( \begin{array}{c} -1 \\ -2-5\, \mathbf i \\ 1 \\ -2+5\, \mathbf i \end{array} \right) . $$

Ответ. $ ( -1,\ -2-5\, \mathbf i, \ 1, \ -2+5\, \mathbf i ) $.

Введенное определение позволяет описать проблему тригонометрической интерполяции как вычисление ДПФ системы значений интерполируемой функции: если $$ \mathcal F(y_1,\dots, y_{2n+1})=(A_0,\dots, A_{2n}) \ , $$ то интерполяционный тригонометрический полином строится по формуле $$ g_n(x) \equiv \frac{1}{2n+1} \left[ A_0 +A_1 z+A_2z^2+\dots+A_nz^n+A_{n+1}z^{-n}+A_{n+2}z^{-n+1}+\dots+A_{2n}z^{-1}\right] \quad npu \quad z=\cos \, x + \mathbf i \sin \, x \ . $$

Такая нумерация коэффициентов полинома позволяет развить альтернативный подход к изложению теории ДПФ. Вернемся к исходной точке: будем искать коэффициенты полинома из системы линейных уравнений. Если $ x_j=2\pi j /(2n+1) $, то соответствующее значение $ z_{j} $ будет снова корнем из $ 1_{} $ степени $ 2n+1 $, как и любая положительная степень $ z_{j} $: $$ z_j= \cos \frac{2\pi j}{2n+1} +\mathbf i \sin \frac{2\pi j}{2n+1}=\varepsilon_j, z_j^n=\varepsilon_j^n \ , $$ а отрицательные степени $ z_{j} $ переписываются через положительные степени $ \varepsilon_j $: $$ z_j^{-n}=\varepsilon_j^{-n}=\varepsilon_j^{n+1}, z_j^{-n+1}= \varepsilon_j^{n+2},\dots, z_j^{-1}= \varepsilon_j^{2n} \ . $$ Следовательно, система линейных уравнений имеет вид: $$ \frac{1}{2n+1} \left( \begin{array}{lllll} 1 & 1 & 1 & \dots & 1 \\ 1 & \varepsilon_{1} & \varepsilon_{1}^2 & \dots & \varepsilon_{1}^{2n} \\ 1 & \varepsilon_{2} & \varepsilon_{2}^2 & \dots & \varepsilon_{2}^{2n} \\ \vdots & & & & \vdots \\ 1 & \varepsilon_{2n} & \varepsilon_{2n}^2 & \dots & \varepsilon_{2n}^{2n} \\ \end{array} \right) \left( \begin{array}{l} A_0 \\ A_1 \\ \vdots \\ A_{2n} \end{array} \right)= \left( \begin{array}{l} y_1 \\ y_2 \\ \vdots \\ y_{2n+1} \end{array} \right) \ . $$ При таком подходе матрица Вандермонда возникает в еще более естественной манере: фактически мы сводим задачу тригонометрической интерполяции к задаче интерполяции алгебраической. Теперь надо решить получившуюся систему. Однако ответ нам уже известен! И этот ответ — то есть ДПФ — записан снова с помощью матрицы Вандермонда, которая отличается от только что полученной только сменой знаков у индексов. Если обозначить эту последнюю через $ V(1,\varepsilon_{1},\varepsilon_{2},\dots, \varepsilon_{2n}) $, а первую — через $ V(1,\varepsilon_{-1},\varepsilon_{-2},\dots, \varepsilon_{-2n}) $ то получаем два соотношения $$ \frac{1}{2n+1} V(1,\varepsilon_{1},\varepsilon_{2},\dots, \varepsilon_{2n}) \left( \begin{array}{l} A_0 \\ A_1 \\ \vdots \\ A_{2n} \end{array} \right) =\left( \begin{array}{l} y_1 \\ y_2 \\ \vdots \\ y_{2n+1} \end{array} \right) \quad u \quad V(1,\varepsilon_{-1},\varepsilon_{-2},\dots, \varepsilon_{-2n}) \left( \begin{array}{l} y_1 \\ y_2 \\ \vdots \\ y_{2n+1} \end{array} \right) = \left( \begin{array}{l} A_0 \\ A_1 \\ \vdots \\ A_{2n} \end{array} \right) \ , $$ откуда выводим изящную формулу: $$ V(1,\varepsilon_{1},\varepsilon_{2},\dots, \varepsilon_{2n})^{-1} = \frac{1}{2n+1} V(1,\varepsilon_{-1},\varepsilon_{-2},\dots, \varepsilon_{-2n}) \ . $$ Обращение матрицы интерполяционной системы оказывается крайне простым!

Опять-таки, отвлекаясь на время от задачи интерполяции, введем для такой удивительной матрицы специальное название. Для произвольного натурального $ N_{} $ назовем матрицу $$ F=\left[ \varepsilon_j^{k} \right]_{j,k=0}^{N-1}= \left( \begin{array}{lllll} 1 & 1 & 1 & \dots & 1 \\ 1 & \varepsilon_1 & \varepsilon_1^2 & \dots & \varepsilon_1^{N-1} \\ 1 & \varepsilon_2 & \varepsilon_2^2 & \dots & \varepsilon_2^{N-1} \\ 1 & \varepsilon_3 & \varepsilon_3^2 & \dots & \varepsilon_3^{N-1} \\ \vdots & & & & \vdots \\ 1 & \varepsilon_{N-1} & \varepsilon_{N-1}^{2} & \dots & \varepsilon_{N-1}^{N-1} \end{array} \right)_{N\times N} \quad npu \quad \varepsilon_j = \cos \frac{2 \pi j}{N} + {\mathbf i} \, \sin \frac{2 \pi j}{N} $$ матрицей дискретного преобразования Фурье. Для нее действительно оказывается справедливым правило вычисления обратной матрицы: $$ F^{-1}=\frac{1}{N}\left[ \varepsilon_{-j}^{k} \right]_{j,k=0}^{N-1} \ . $$

Непосредственное доказательство этого свойства, а также вывод других свойств матрицы ДПФ ☞ ЗДЕСЬ.

Дискретное преобразование Фурье, рассмотренное в предыдущем пункте, представляет собой только лишь альтернативную форму записи решения задачи тригонометрической интерполяции с равноотстоящими узлами, приведенного в конце ☞ ПУНКТА — его введение не дает ни какого вычислительного преимущества.

Подсчитаем количество элементарных операций, необходимых для вычисления ДПФ в общем случае: $$ \mathfrak X_k=\sum_{j=0}^{N-1} \mathfrak x_j \varepsilon_{-j}^k \quad npu \quad k\in\{0,\dots, N-1\} ; $$ числа $ \{\mathfrak x_j\}_{j=0}^{N-1} $ предполагаются комплексными. Вычисление $ \mathfrak X_k $ «в лоб» требует $ N_{} $ умножений и $ N-1 $ сложений. Следовательно, вычисление всех чисел $ \{\mathfrak X_k\}_{k=0}^{N-1} $ требует $ N^2 $ умножений и $ N(N-1) $ сложений.

Попробуем теперь оптимизировать вычисления — например, группировкой слагаемых, входящих в формулы.

Пример. Для случая $ N_{}=6 $ сначала перегруппируем суммы, составляющие ДПФ, отделив слагаемые, содержащие $ \mathfrak x_j $ с четными индексами от слагаемых, содержащих $ \mathfrak x_j $ с нечетными индексами; потом воспользуемся свойствами корней $ 6 $-й степени из единицы:

$$ \varepsilon_{-3}=-1, \varepsilon_{-4}=-\varepsilon_{-1}, \varepsilon_{-5}=-\varepsilon_{-2} \ . $$ $$ \begin{array}{lccccc} \mathfrak X_0=& (\mathfrak x_0+\mathfrak x_2+\mathfrak x_4) & +(\mathfrak x_1+\mathfrak x_3+\mathfrak x_5) & & \\ \mathfrak X_1=& (\mathfrak x_0+\varepsilon_{-1}^2 \mathfrak x_2+\varepsilon_{-1}^4\mathfrak x_4) & +\varepsilon_{-1}(\mathfrak x_1+\varepsilon_{-1}^2\mathfrak x_3+\varepsilon_{-1}^4\mathfrak x_5) &= (\mathfrak x_0+\varepsilon_{-2} \mathfrak x_2+\varepsilon_{-2}^2\mathfrak x_4) & + \varepsilon_{-1}(\mathfrak x_1+\varepsilon_{-2} \mathfrak x_3+\varepsilon_{-2}^2\mathfrak x_5) \\ \mathfrak X_2=& (\mathfrak x_0+\varepsilon_{-2}^2 \mathfrak x_2+\varepsilon_{-2}^4\mathfrak x_4) & +\varepsilon_{-2}(\mathfrak x_1+\varepsilon_{-2}^2\mathfrak x_3+\varepsilon_{-2}^4\mathfrak x_5) &= (\mathfrak x_0+\varepsilon_{-4} \mathfrak x_2+\varepsilon_{-4}^2\mathfrak x_4) &+ \varepsilon_{-2}(\mathfrak x_1+\varepsilon_{-4} \mathfrak x_3+\varepsilon_{-4}^2\mathfrak x_5) \\ \mathfrak X_3=& (\mathfrak x_0+\mathfrak x_2+\mathfrak x_4) & - (\mathfrak x_1+\mathfrak x_3+\mathfrak x_5) & & \\ \mathfrak X_4=& (\mathfrak x_0+\varepsilon_{-4}^2 \mathfrak x_2+\varepsilon_{-4}^4\mathfrak x_4) & +\varepsilon_{-4}(\mathfrak x_1+\varepsilon_{-4}^2\mathfrak x_3+\varepsilon_{-4}^4\mathfrak x_5) &= (\mathfrak x_0+\varepsilon_{-2} \mathfrak x_2+\varepsilon_{-2}^2\mathfrak x_4) &- \varepsilon_{-1}(\mathfrak x_1+\varepsilon_{-2} \mathfrak x_3+\varepsilon_{-2}^2\mathfrak x_5) \\ \mathfrak X_5=& (\mathfrak x_0+\varepsilon_{-5}^2 \mathfrak x_2+\varepsilon_{-5}^4\mathfrak x_4) &+ \varepsilon_{-5}(\mathfrak x_1+\varepsilon_{-5}^2\mathfrak x_3+\varepsilon_{-5}^4\mathfrak x_5) &= (\mathfrak x_0+\varepsilon_{-4} \mathfrak x_2+\varepsilon_{-4}^2\mathfrak x_4) &- \varepsilon_{-2}(\mathfrak x_1+\varepsilon_{-4} \mathfrak x_3+\varepsilon_{-4}^2\mathfrak x_5) \end{array} $$ В получившихся формулах наблюдается следующая закономерность: вычисления $ \mathfrak X_0 $ и $ \mathfrak X_3 $ фактически дублируют друг друга (с различием лишь в знаке одного слагаемого); то же самое можно утверждать относительно пар $ \mathfrak X_1 $ и $ \mathfrak X_4 $, а также $ \mathfrak X_2 $ и $ \mathfrak X_5 $. Иными словами, грамотная организация процесса вычисления $ \mathfrak X_0, \mathfrak X_1, \mathfrak X_2 $ позволит фактически «даром» получить значения $ \mathfrak X_3, \mathfrak X_4, \mathfrak X_5 $.

Рассмотрим теперь повнимательней каждую из скобок $ ( \ ) $. Все скобки, содержащие $ \mathfrak x_0,\mathfrak x_2,\mathfrak x_4 $, имеют одинаковую структуру: они представляют суммы этих величин, домноженных на степени $$ \varepsilon_{-2} = \cos \frac{4\pi}{6} - \mathbf i \, \sin \frac{4\pi}{6} = \cos \frac{2\pi}{3} - \mathbf i \, \sin \frac{2\pi}{3} \ . $$ Последняя величина, с одной стороны является корнем $ 6_{} $-й степени из единицы, а, с другой стороны, является корнем $ 3_{} $-й степени из единицы. Если ввести обозначение этих корней: $$\tilde \varepsilon_{0}=1,\ \tilde \varepsilon_{-1}=\cos \frac{2\pi}{3} - \mathbf i \, \sin \frac{2\pi}{3},\ \tilde \varepsilon_{-2}=\cos \frac{4\pi}{3} - \mathbf i \, \sin \frac{4\pi}{3} \ , $$ то скобки,содержащие $ \mathfrak x_0,\mathfrak x_2,\mathfrak x_4 $, переписываются в виде $$ \begin{array}{rl} \mathfrak Y_0 = & \mathfrak x_0 +\mathfrak x_2\, \tilde \varepsilon_{0}+\mathfrak x_4\, \tilde \varepsilon_{0}^2, \\ \mathfrak Y_1 = & \mathfrak x_0 +\mathfrak x_2\, \tilde \varepsilon_{-1}+\mathfrak x_4\, \tilde \varepsilon_{-1}^2, \\ \mathfrak Y_2 = & \mathfrak x_0 +\mathfrak x_2\, \tilde \varepsilon_{-2}+\mathfrak x_4\, \tilde \varepsilon_{-2}^2, \end{array} $$ т.е. они представляют ДПФ набора $ \mathfrak x_0,\mathfrak x_2,\mathfrak x_4 $. Отличие от исходной задачи заключается в уменьшении размерности: исходная размерность была шестой, получившаяся — третьей. Совершенно аналогично расправляемся со скобками, содержащими $ \mathfrak x_1,\mathfrak x_3,\mathfrak x_5 $: $$ \begin{array}{rl} \mathfrak Z_0 = & \mathfrak x_1 +\mathfrak x_3\, \tilde \varepsilon_{0}+\mathfrak x_5\, \tilde \varepsilon_{0}^2, \\ \mathfrak Z_1 = & \mathfrak x_1 +\mathfrak x_3\, \tilde \varepsilon_{-1}+\mathfrak x_5\, \tilde \varepsilon_{-1}^2, \\ \mathfrak Z_2 = & \mathfrak x_1 +\mathfrak x_3\, \tilde \varepsilon_{-2}+\mathfrak x_5\, \tilde \varepsilon_{-2}^2; \end{array} $$ они также представляют ДПФ набора $ \mathfrak x_1,\mathfrak x_3,\mathfrak x_5 $. Окончательно: $$ \begin{array}{ccc} \mathfrak X_0=\mathfrak Y_0+\mathfrak Z_0, & \mathfrak X_1= \mathfrak Y_1+\varepsilon_{-1}\mathfrak Z_1, & \mathfrak X_2= \mathfrak Y_2+\varepsilon_{-2}\mathfrak Z_2, \\ \mathfrak X_3=\mathfrak Y_0-\mathfrak Z_0, & \mathfrak X_4= \mathfrak Y_1-\varepsilon_{-1}\mathfrak Z_1, & \mathfrak X_5= \mathfrak Y_2-\varepsilon_{-2}\mathfrak Z_2. \end{array} $$ Полученные формулы сводят вычисление ДПФ от шестого порядка к третьему. ♦

Идея очевидно распространяется на случай произвольного четного числа $ N_{} $: $$ \mathfrak X_k = \left\{ \begin{array}{llll} \mathfrak Y_k & +\varepsilon_{-k} & \times \mathfrak Z_k & npu \ k\in \{0,1,\dots,N/2-1\}, \\ \mathfrak Y_{k-N/2}& -\varepsilon_{-(k-N/2)}& \times \mathfrak Z_{k-N/2} & npu \ k\in \{N/2,N/2+1,\dots,N-1\}; \end{array} \right. $$ при $$ \varepsilon_{-k}=\cos \frac{2\pi k}{N}- \mathbf i \sin \frac{2\pi k}{N}=e^{-\mathbf i (2\pi k/N)} \ . $$ Здесь набор $ \mathfrak Y_0,\mathfrak Y_1, \dots , \mathfrak Y_{N/2-1} $ представляет собой ДПФ от набора $ \mathfrak x_0,\mathfrak x_2,\dots, \mathfrak x_{N-2} $, а набор $ \mathfrak Z_0,\mathfrak Z_1, \dots , \mathfrak Z_{N/2-1} $ — ДПФ от набора $ \mathfrak x_1,\mathfrak x_3,\dots, \mathfrak x_{N-1} $; число $ \varepsilon_{-k} $ называется дополнительным множителем²⁾. Если число $ N/2 $ само является четным, то процедура повторяется. Если число $ N_{} $ представляет собой степень двойки: $ N=2^{\eta} $, то рекурсивный алгоритм сведет вычисление ДПФ порядка $ N_{} $ к вычислению ДПФ порядка $ 2_{} $. Этот алгоритм вычисления ДПФ называется бинарным прореживанием по времени³⁾.

Пример. Для случая $ N_{}=8 $ схема вычислений будет следующей:

♦

Подсчитаем теперь стоимость бинарного прореживания в элементарных операциях. Если считать, что ДПФ порядка $ N/2 $ уже вычислены, то для вычисления ДПФ порядка $ N_{} $ требуется $ 2\, N $ операций ( $ N_{} $ умножений величин $ \mathfrak Z_{\ell} $ на дополнительные сомножители и $ N_{} $ сложений). Аналогично, для вычисления каждого ДПФ порядка $ N/2 $ потребуется $ 2(N/2) $ операций; и нам требуется вычисление двух таких ДПФ (см. предыдущий пример). Таким образом, в случае $ N=2^{\eta} $ общее количество операций будет равно $$ (2\,N)+2\left(2\frac{N}{2}\right)+4\left(2\frac{N}{4}\right)+\dots \ , $$ и сумма завершается через $ \eta $ шагов когда мы добираемся до ДПФ вторых порядков, каждое из которых стоит две операции⁴⁾. В результате общая стоимость вычислений $ 2\,N \, \eta = 2\, N \log_2 N $, что существенно меньше приведенной в начале пункта стоимости вычислений «в лоб». Так, при $ N=2^{10} $ имеем $$ N^2=1\, 048\, 576, \quad N \log_2 N = 10\, 240 \ . $$

Другой подход к уменьшению стоимости вычислений ДПФ проиллюстрируем снова на примере.

Пример. Для случая $ N_{}=6 $ сначала сгруппируем отдельно входящие в состав ДПФ равенства для $ \mathfrak X_k $ с четными индексами и равенства для $ \mathfrak X_k $ с нечетными индексами; потом объединим скобками слагаемые, содержащие пары $ \mathfrak x_0 $ и $ \mathfrak x_3 $, $ \mathfrak x_1 $ и $ \mathfrak x_4 $, $ \mathfrak x_2 $ и $ \mathfrak x_5 $:

$$ \left\{\begin{array}{lrrrr} \mathfrak X_0& &=(\mathfrak x_0+\mathfrak x_3) &+(\mathfrak x_1+\mathfrak x_4) &+(\mathfrak x_2 +\mathfrak x_5) \\ \mathfrak X_2&=\mathfrak x_0+\varepsilon_{-2} \mathfrak x_1+\varepsilon_{-2}^2 \mathfrak x_2+\varepsilon_{-2}^3 \mathfrak x_3+\varepsilon_{-2}^4 \mathfrak x_4+\varepsilon_{-2}^5 \mathfrak x_5 &= (\mathfrak x_0+\mathfrak x_3)&+\varepsilon_{-2}(\mathfrak x_1+\mathfrak x_4)&+ \varepsilon_{-2}^2(\mathfrak x_2+\mathfrak x_5)\\ \mathfrak X_4&=\mathfrak x_0+\varepsilon_{-4} \mathfrak x_1+\varepsilon_{-4}^2 \mathfrak x_2+\varepsilon_{-4}^3 \mathfrak x_3+\varepsilon_{-4}^4 \mathfrak x_4+\varepsilon_{-4}^5 \mathfrak x_5 & = (\mathfrak x_0+\mathfrak x_3)&+\varepsilon_{-4}(\mathfrak x_1+\mathfrak x_4)&+ \varepsilon_{-4}^2(\mathfrak x_2+\mathfrak x_5) \end{array} \right. $$ Видим, что получившиеся формулы имеют снова структуру ДПФ — только, в отличие от исходного — не $ 6 $-го, а $ 3_{} $-го порядка: величина $$ \varepsilon_{-2} = \cos \frac{4\pi}{6} - \mathbf i \, \sin \frac{4\pi}{6} = \cos \frac{2\pi}{3} - \mathbf i \, \sin \frac{2\pi}{3} $$ является корнем $ 3_{} $-й степени из единицы, а входной последовательностью является $ \mathfrak x_0+\mathfrak x_3, \mathfrak x_1+\mathfrak x_4, \mathfrak x_2+\mathfrak x_5 $.

Для $ \mathfrak X_k $ c нечетными индексами, имеем: $$ \left\{ \begin{array}{lrrrr} \mathfrak X_1 &=\mathfrak x_0+\varepsilon_{-1}\mathfrak x_1+\varepsilon_{-1}^2\mathfrak x_2 +\varepsilon_{-1}^3\mathfrak x_3+\varepsilon_{-1}^4\mathfrak x_4+\varepsilon_{-1}^5\mathfrak x_5 &=(\mathfrak x_0-\mathfrak x_3)&+\varepsilon_{-1}(\mathfrak x_1-\mathfrak x_4)&+\varepsilon_{-1}^2 (\mathfrak x_2-\mathfrak x_5) \\ \mathfrak X_3&=\mathfrak x_0+\varepsilon_{-3}\mathfrak x_1+\varepsilon_{-3}^2\mathfrak x_2 +\varepsilon_{-3}^3\mathfrak x_3+\varepsilon_{-3}^4\mathfrak x_4+\varepsilon_{-3}^5\mathfrak x_5 & =(\mathfrak x_0-\mathfrak x_3)&+\varepsilon_{-3}(\mathfrak x_1-\mathfrak x_4)&+\varepsilon_{-3}^2 (\mathfrak x_2-\mathfrak x_5) \\ \mathfrak X_5&=\mathfrak x_0+\varepsilon_{-5}\mathfrak x_1+\varepsilon_{-5}^2\mathfrak x_2 +\varepsilon_{-5}^3\mathfrak x_3+\varepsilon_{-5}^4\mathfrak x_4+\varepsilon_{-5}^5\mathfrak x_5 & =(\mathfrak x_0-\mathfrak x_3)&+\varepsilon_{-5}(\mathfrak x_1-\mathfrak x_4)&+ \varepsilon_{-5}^2(\mathfrak x_2-\mathfrak x_5) \end{array} \right. $$ Здесь труднее увидеть ДПФ, но всё-таки можно: если за входную последовательность взять $ \mathfrak x_0-\mathfrak x_3,\varepsilon_{-1}(\mathfrak x_1-\mathfrak x_4),\varepsilon_{-1}^2(\mathfrak x_2-\mathfrak x_5) $, то снова получим ДПФ третьего порядка! Таким образом, перегруппировка формул позволяет свести вычисление ДПФ от шестого порядка к третьему. ♦

Идея очевидно распространяется на случай произвольного четного числа $ N_{} $: на основе входной последовательности $ \mathfrak x_0,\mathfrak x_1,\dots, \mathfrak x_{N-1} $ формируются две новые последовательности из $ N/2 $ элементов: $$ \begin{array}{lr} \mathfrak y_j= & \mathfrak x_j+\mathfrak x_{N/2+j} \\ \mathfrak z_j= & \varepsilon_{-j}(\mathfrak x_j-\mathfrak x_{N/2+j}) \end{array} \quad npu \quad j\in \{1,2,\dots, N/2\} ;\quad \varepsilon_{-j}=\cos \frac{2\pi j}{N}- \mathbf i \sin \frac{2\pi j}{N}=e^{-\mathbf i (2\pi j/N)} . $$ К каждой из них применяется ДПФ порядка $ N/2 $, из выходных последовательностей формируется $ \mathfrak X_0,\mathfrak X_1,\dots, \mathfrak X_{N-1} $.

В отличие от предыдущего метода, основанного на измельчении входной последовательности $ \{\mathfrak x_j\}_{j=0}^{N-1} $, в настоящем методе предлагается измельчать выходную последовательность $ \{\mathfrak X_k\}_{k=0}^{N-1} $; он называется бинарным прореживанием по частоте⁵⁾.

Идею метода может прояснить и аналогия с задачей алгебраической (степенной) интерполяции ☞ ЗДЕСЬ.

Оба рассмотренных метода относятся к классу методов вычисления ДПФ, известных под названием быстрого преобразования Фурье⁶⁾.

Задача. Построить тригонометрический полином порядка $ n_{} $, имеющий в узлах интерполяции $ x_{1},\dots,x_n $ заданные значения и заданные значения своей производной: $$ g_n(x_j)=y_j, g_n^{\prime}(x_j)=y_j^{\prime} \quad npu \quad j\in\{1,\dots,n\} \ . $$

Вспомним, что тригонометрический полином порядка $ n_{} $ в общем случае имеет $ 2n+1 $ коэффициент. Поставленная задача формирует условия на эти коэффициенты в виде линейных уравнений. Поскольку число $ 2n_{} $ этих уравнений меньше числа определяемых коэффициентов, то задача, скорее всего (см. ☞ ЗДЕСЬ ), будет иметь бесконечное множество решений. Так оно и оказывается.

Оппенгейм А., Шафер Р. Цифровая обработка сигналов. М.Техносфера. 2009.

Полиа Г., Сеге Г. Задачи и теоремы из анализа. М.Наука. 1978. Часть 2, отдел 6, сс. 85-90.

Турецкий А.Х. Теория интерполирования в задачах. Минск. Вышэйшая школа. 1968.