1626435587-55f52a4de97976f3c6215fa7c103f544 (Смирнов 2017 - Основы вычислительной физики ч2), страница 10

8 использована формула = 20 log10 (1 /2 ).14мощно-стейамплитуд48(б)0.30.20.1-2010-2-4010-30-4-2024-60-4Частота ω, отн.ед.Рис. 9.10-1децибелы0.41,42 дБ31,5 дБ0.5Спектр |h(ω)|Спектр |h(ω)|(а)-2024Частота ω, отн.ед.Спектр оконной функции Ханна |ℎ̃()| (72) в линейном (а) и логарифмическом (б) масштабестрируемых спектрах. В случае прямоугольного окна отношение высотглавного и первого побочного лепестка составляет всего 13,5 дБ.Аналогично для оконной функции Ханна можем записать:⃒⃒(︂)︂⃒⃒122 ⃒ sinc ( /2) ⃒ℎ() =1 − cos, |ℎ̃()| = 2 ⃒.(72)222( ) − 4 ⃒График спектральной функции (72) показан на рис. 9 в линейном илогарифмическом масштабе. По сравнению со спектром прямоугольного окна (см. рис. 8) видно несколько отличий: соотношение высотыглавного и побочного пика в спектре возросло с 13,5 до 31,5 дБ, чтопозволяет значительно уменьшить эффект частокола при использовании окна Ханна.

Однако при этом мы наблюдаем также увеличениеширины центрального пика, что приводит к снижению разрешения вспектре, содержащем лишь «целые» частоты (период которых делитнацело ширину сетки).Таким образом, мы видим, что использование оконных функцийявляется компромиссным решением. Использование любой оконнойфункции подразумевает внесение некоторых искажений в регистрируемый спектр сигнала, а выбор той или иной оконной функции обусловлен спецификой конкретной задачи и зависит от допустимости тех илииных искажений.2.10. Быстрое преобразование ФурьеПреобразование Фурье играет чрезвычайно важную роль и повсеместно используется для цифровой обработки сигналов, в том числе49потоков аудио- и видеоданных, ставших неотъемлемой частью современной жизни.

Необходимость больших объёмов вычислений (многиеиз которых требуется выполнять в реальном масштабе времени) диктует спрос на экономичные алгоритмы и их эффективные программныеи аппаратные реализации. Заметим, что возникновение практического интереса к эффективным вычислительным алгоритмам возникло заполвека до появления Youtube. Актуальные задачи тех дней включалианализ сейсмоданных для обнаружения испытаний атомного оружиявероятным противником.Вычисление дискретного преобразование Фурье непосредственно поопределению (56) на с.

33 требует (2 ) арифметических операций. Алгоритмы, позволяющие решить данную задачу за меньшее число шагов, принято называть(БПФ). Начало активного развития данного направления связано с публикациейработы [A4] в 1965 году, хотя, как выяснилось уже после выхода статьи,базовые идеи были высказаны значительно раньше.К настоящему времени разработано значительное количество алгоритмов быстрого преобразования Фурье (и ещё большее количество ихразнообразных описаний). Для наших целей, однако, будет вполне достаточным поверхностного знакомства с основными идеями, позволяющими достичь значительной экономии при выполнении дискретногопреобразования Фурье, а также с одной из наиболее эффективных универсальных программных реализаций, о которой пойдёт речь в п.

2.11.Достаточно компактная формальная запись формул, позволяющихуменьшить число арифметических операций за счёт перегруппировкислагаемых в (56), может быть найдена в пионерской работе [A4]15 . Однако для наших целей представляется более полезным и интереснымначать рассмотрение основной идеи алгоритмов БПФ, проведя аналогию с вычислением полиномов, как это сделано в замечательной монографии [A6, гл.

2].быстрым преобразованием Фурье2.10.1. БПФ и вычисление полиномовЗаметим, что матрица × дискретного преобразования Фурье()Λ ∝ , где ≡ exp(2/), с точностью до коэффициента пропорциональности совпадает с матрицей Вандермонда ∆ = при = .Другими словами, выполнение дискретного преобразования Фурье длявектора (0 , 1 , . . . , −1 ) в соответствии с выражением15Бесплатную копию можно найти в Интернет, используя сервис Google Scholar.50⎛01...⎜⎜⎜⎜⎜⎝ −2−1⎞⎛11...⎟ ⎜⎟ ⎜⎟ ⎜⎟=⎜⎟ ⎜⎠ ⎝ 1111.........−2−1......11−2...−1...(−2)(−2)(−1)(−2)(−2)(−1)(−1)(−1)⎞⎛01...⎟⎜⎟⎜⎟⎜⎟⎜⎟⎜⎠ ⎝ −2−1⎞⎟⎟⎟⎟⎟⎠эквивалентно вычислению значений полинома16 () = 0 + 1 1 + . .

. + −1 −1в точках — комплексных корнях -й степени из единицы:{ } = { }, = 0, 1, . . . , − 1, = exp(2/).(73)Непосредственное вычисление значений ( ) полинома потребуетприблизительно 22 арифметических операций над комплексными числами (полагаем, что таблица значений построена предварительно ииспользуется для многократного вычисления полиномов с разными ).Как добиться экономии при вычислении значений полинома? Для простоты ограничимся рассмотрением случая = 2 , где — целое число.Разобьём полином на сумму двух многочленов разной чётности: () = /2 (2 ) + /2 (2 ) =(︀)︀= 0 + 2 2 + .

. . + −2 −2 +(︀)︀+ 1 + 3 2 + . . . + −1 −2 · .(74)Несложно увидеть (рис. 10), что при чётных для каждой точки в наборе (73) присутствует также симметрично расположенная точка = − , = ( +/2) mod . Следовательно, вычисление значенийполинома ( ) можно оптимизировать, разбив набор точек (73) на/2 пар: ±0 , ±1 , . . . , ± 2 −1 . В каждой паре для нахождения (±)необходимо вычислить значения полиномов /2 и /2 : ()= /2 (2 ) + · /2 (2 ), (−)= /2 (2 ) − · /2 (2 ).(75)Использование соотношений (75) позволяет получить ответ всего за2 × 2( 2 )2 + 32 ≈ 2 операций, что почти вдвое меньше первоначальнойоценки 22 .16 В данной задаче в качестве нижнего индекса в обозначении полиномов удобнее использовать количество полиномиальных коэффициентов (размерностьпространства), а не степень полиномов.51(а)Im(б)xjImReRe-x j-x jРис.

10.xjСимметричное расположение точек = при = 2 на примере(а) = 16, (б) = 8Однако самое замечательное свойство использованного подхода состоит в том, что его можно применять многократно, что позволяетне просто уменьшить коэффициент в (2 ), но даже изменить асимптотику числа операций. Действительно, формулы (75) сводят задачу о вычислении полинома в точках = к двум задачам овычислении многочленов /2 и /2 в точках 2 = 2 . Поскольку2 = exp(2/(/2)) = /2 , легко заметить, что использование формул (75) приводит нас к условиям первоначальной задачи с заменой на ′ = /2 = 2−1 , ср. рис. 10 (а) и (б). Это позволяет использоватьданный подход многократно, на каждом шаге удваивая число задач иуменьшая вдвое размерность пространства.Посчитаем полное число арифметических операций, которые необходимо совершить для вычисления в = 2 точках { } (73) прирекурсивном использовании формул (75). Поскольку на каждом шагерекурсии число полиномов удваивается, а число коэффицентов уменьшается вдвое, на -м шаге ( = 0, 1, .

. . , −1) будем иметь 2 полиномовс = /2 коэффициентами в каждом. На последнем шаге алгоритмабудем иметь = − 1, −1 = 2, 2 = 2/2 = −1, а для вычисления каждого из 2−1 полиномов 2 в точках {20 , 21 } потребуется двеарифметических операции:2 (+1)= 0 + 1 ,2 (−1)= 0 − 1 .(76)На остальных шагах ( = 0, . . . , − 2) в соответствии с выражением(75) потребуется 3 операции над полиномами и (умножение · ивычисление суммы и разности ±·) для каждого из 2 полиномов вкаждой из /2 = /2+1 пар точек. Следовательно, всего на -м шагеалгоритма необходимо сделать 3 × 2 × (/2+1 ) = 23 арифметическихопераций. Учитывая, что алгоритм включает ( − 1) таких шагов плюс52самый последний шаг при = − 1 (76), получаем для полного числаопераций ( × ) = ( log ).2.10.2.

Запись БПФ через матрицыСформулируем описанный выше алгоритм также на языке матриц,более привычном в контексте данной задачи. Запишем дискретное преобразование Фурье -мерного вектора a = (0 , 1 , . . . , −1 ) в матричной форме, по-прежнему полагая = 2 :b = Λ() a,где()Λ = , = exp(2/).(77)Для использования симметрии матрицы Λ() дискретного преобразования Фурье перейдём теперь к другому базису в пространстве исходныхданных, изменив порядок элементов.

Расположим вначале /2 чётныхкомпонент вектора a, а затем оставшиеся /2 нечётных, что являетсяаналогом записи полинома в виде суммы чётной и нечётной частей(74). Изменение порядка элементов базиса приведёт к соответствующейперестановке столбцов матрицы Λ() . Записав матрицу Λ() в новомбазисе в блочном виде (78), получим, что блоки (/2) и (/2) будутвключать чётные столбцы Λ() , а блоки (/2) и (/2) — нечётные:⎛⎞ ⎛⎞⎛⎞00⎜⎟ ⎜⎟⎜⎟(/2) (/2)⎜ 1 ⎟ ⎜⎟ ⎜ 2 ⎟⎜ .

⎟ ⎜⎟⎜ . ⎟⎜ . ⎟=⎜⎟⎜ . ⎟(78)⎜ . ⎟ ⎜⎟⎜ . ⎟⎜⎟ ⎜⎟⎜⎟(/2)(/2)⎝ −2 ⎠ ⎝ ⎠ ⎝ −3 ⎠−1−1Для матричных элементов блока (/2) имеем:(/2)()= Λ,2 = 2 = /2,где, = 0, . . . ,− 1.2(79)/2Из определения (77) следует, что = −1, а = 1, что позволяетвыразить три других блока через (/2) :(/2)(/2)(/2)()= Λ+/2,2==()Λ,2+1()Λ+/2,2+12(+/2)= ==·(2+1)==(+/2)(2+1)=(/2),(/2) ,(/2)− .(80)Формулы (78–80) сводят задачу о дискретном преобразовании Фурье -мерного вектора a к двум задачам вдвое меньшей размерности53для чётных и нечётных компонент вектора a, что даёт матричное описание рекурсивного алгоритма, рассмотренного выше в терминах полиномов.Несложно построить альтернативную версию рассмотренного алгоритма без использования рекурсии и показать, что быстрое преобразование Фурье позволяет сократить число операций без использования дополнительной памяти для хранения промежуточных вычислений.