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

Однако это выходит за рамки поставленных нами целей, в связис чем мы направляем интересующихся читателей к монографии [A6].2.10.3. Общий случай составных Выше мы рассматривали алгоритм БПФ для частного случая массивов длины = 2 . Указанное ограничение на размер таблиц данныххотя и позволяет достичь наибольшей экономии числа операций, ноне является обязательным. Чтобы построить более общий алгоритмБПФ заметим, что экономия числа операций достигается за счёт использования симметрии расположения комплексных корней из единицы { } = { }. В рассмотренном выше частном случае = 2 мывоспользовались симметрией относительно замены → −, что эквивалентно домножению на , или повороту комплексной плоскости на .

Аналогичным образом можно использовать симметрию относительно поворота на угол 2/, где — делитель числа точек . Рассмотрим данную мысль на примере = 3. Полагая кратным 3 и действуяпо аналогии с (74), разобьём полином на сумму трёх многочленов,каждый из которых является собственной функцией преобразования → 2/3 с собственным значением 1 либо ±2/3 : () = /3 (3 ) + /3 (3 ) + 2 /3 (3 ) =(︀)︀= 0 + 3 3 + . .

. + −3 −3 +(︀)︀+ 1 + 4 3 + . . . + −2 −3 · +(︀)︀+ 2 + 5 3 + . . . + −1 −3 · 2 .Поскольку mod 3 = 0, для каждой точки в наборе (73) присутствует также симметрично расположенные точки = ±2/3 , = ( ± /3) mod . Следовательно, вычисление значений полинома ( ) можно оптимизировать, если разбить набор точек (73) на /3троек.

В каждой тройке для нахождения необходимо вычислитьзначения полиномов /3 , /3 и /3 :54 () = /3 (3 ) + ( (2343) = /3 (3 ) + 3) = /3 ( ) + · /3 (3 ) +2343 · /3 (3 ) + 3 · /3 ( ) + 2 · /3 (3 ),43232 · /3 (3 ),2(81)3 · /3 ( ).Использование соотношений (81) позволяет получить ответ всегоза 3 × 2( 3 )2 + 12 3 ≈ 23 2 операций (вместо 22 при непосредственномвычислении по формулам (56) на с. 33).Аналогичные рассуждения можно произвести для произвольногосоставного числа точек . Краткое формальное рассмотрение общегослучая = 1 · 2 можно найти в статье [A4].

Перегруппировка слагаемых в сумме позволяет свести исходную задачу о выполнении дискретного преобразования Фурье массива длины к вычислению суммы из1 слагаемых, каждое из которых есть результат БПФ от подмножестваисходного массива длины 2 .2.10.4. Cлучай простого Работа описанного выше алгоритма основана на разделении задачидля -точечного дискретного преобразования Фурье на последовательность подзадач с кратно меньшим числом точек. Каждый шаг деления даёт выигрыш в быстродействии, так что общая экономия будет увеличиваться с ростом числа шагов.

Последнее определяется, очевидно,числом сомножителей, на которые раскладывается длина исходного массива. Можно ли достичь какого-либо ускорения в случае, если — простое число, т. е. непредставимо в виде произведения целочисленных сомножителей? Через три года после выхода работы Кули иТьюки было показано [A7], что в случае простых дискретное преобразование Фурье может быть сведено к вычислению циклической корреляции (либо свёртки) двух последовательностей длины ( − 1). Данная задача может быть решена с помощью дискретных преобразованийФурье17 .

Поскольку ( − 1) не является простым числом (предполагаем, что простое и > 3), задача сводится к рассмотренному вышеалгоритму быстрого преобразования Фурье, что позволяет решить еёза ( log ) арифметических операций, хотя коэффициент пропорциональности может быть значительно больше, чем в случае составных и тем более = 2 .17Фурье-образ от свёртки * равен произведению Фурье-образов и .55Рассмотрим идею данного алгоритма более подробно. Данный параграф предназначен лишь для самых увлечённых и любознательныхчитателей, тогда как при первом знакомстве с предметом рекомендуется переходить к п.

2.11. Для наглядности, помимо общих формул,выпишем также матрицы для частного случая = 5. Прежде всего,обратим внимание, что все элементы первой строки и первого столбцаматрицы дискретного преобразования Фурье Λ = exp(2) равны1 (см. матрицу на с. 51). В этой связи для коэффициента Фурье с индексом 0 имеем:−1∑︁0 = .(82)=0Для остальных коэффициентов можем записать:(︂)︂−1∑︁2 exp , = 1, . . . , − 1.

= 0 +(83)=1В частном случае = 5, обозначая для краткости ≡ exp(2/5),имеем:⎛⎞ ⎛ 1⎞⎛⎞1 − 0 2 3 41⎜ 2 − 0 ⎟ ⎜ 2 4 1 3 ⎟ ⎜ 2 ⎟⎜⎟ ⎜⎟⎜⎟⎝ 3 − 0 ⎠ = ⎝ 3 1 4 2 ⎠ ⎝ 3 ⎠4 − 04 3 2 14 .Выполним перестановку строк и столбцов полученной матрицы, заменив индексы и в (83) так, чтобы получилось взаимно однозначноеотображение = 0, 1, . . ., − 2 на = 1, 2, .

. . , − 1 и = 0, 1, . . . , − 2на = 1, 2, . . . , − 1: = mod , = 0, . . . , − 2, = mod , = 0, . . . , − 2.(84)(Числа , обеспечивающие указанное взаимно однозначное отображение, называют.) Например, при = 5 можновыбрать = 2, что даст следующее отображение индексов:примитивными корнями, , 01122433Выполнив замену индексов (84) в выражении (83), получим:(︂)︂−2∑︁2 + mod − 0 = mod exp·mod ,=0 = 0, . .

. , − 2.56(85)Выражение (85) задаёт произведение некоторой матрицы и вектораa, компоненты которого переставлены в соответствии с (84). Обратимвнимание на структуру матрицы: заменяя в (85) → + 1, несложно увидеть, что каждая следующая строка содержит те же элементы,что и предыдущая, со сдвигом столбцов на 1 ( → − 1). Сдвиг является циклическим ввиду того, что + входит в (85) по модулю .В частном случае = 5 имеем:⎛⎞ ⎛ 1⎞⎛⎞1 − 0 2 4 31⎜ 2 − 0 ⎟ ⎜ 2 4 3 1 ⎟ ⎜ 2 ⎟⎜⎟ ⎜⎟⎜⎟⎝ 4 − 0 ⎠ = ⎝ 4 3 1 2 ⎠ ⎝ 4 ⎠3 − 03 1 2 43 .корреляциисвёрткиПолученное выражение носит названиепары векторов( 1 , 2 , 4 , 3 ) и (1 , 2 , 4 , 3 ), либовектора ( 3 , 4 , 2 , 1 ) с(1 , 2 , 4 , 3 ). Как уже отмечалось в начале данного параграфа, даннаяоперация может быть сведена к обратному дискретному преобразованию Фурье от произведения Фурье-образов сворачиваемых векторов.Всего требуется три преобразования Фурье, однако одно из них — отвектора ( 3 , 4 , 2 , 1 ) — может быть выполнено предварительно, сохранено в памяти и использовано многократно для различных a.

Поскольку свёртка вычисляется в пространстве C−1 , а число ( − 1) является составным для простых ≥ 5, можно использовать описанныев пп. 2.10.1-2.10.3 алгоритмы быстрого преобразования Фурье для составных .Таким образом, мы свели задачу о вычислении дискретного преобразования Фурье для простых к использованию алгоритмов быстрогопреобразования Фурье плюс () дополнительных операций, связанных с компонентами 0 и 0 , см. (82) и (85).

Совместное использованиеалгоритмов БПФ для простых и составных формально позволяетвыполнить дискретное преобразование Фурье для произвольного за( log ) арифметических операций. Для желающих изучить вопросболее глубоко, рекомендуем обратиться к оригинальным работам и монографии [A8].2.11. Библиотека FFTWВ заключение познакомимся с одной из наиболее эффективныхкросс-платформенных библиотек для выполнения дискретного преобразования Фурье — FFTW18 [A9].18 От англ.

Fastest Fourier Transform in the West, самое быстрое преобразованиеФурье на Западе.57Библиотека FFTW является бесплатной (в том числе для коммерческого использования), распространяется по лицензии GNU GPL иимеет целый ряд достоинств, в числе которых:∙ высокая скорость работы на различных аппаратных платформах;∙ возможность дополнительного повышения быстродействия засчёт использования различных наборов векторных (SIMD) инструкций и параллельных вычислений в системах с симметричной мультипроцессорностью (SMP) и массивно-параллельной архитектурой (MPP);∙ выполнение одно- и многомерных преобразований Фурье;∙ возможность работы с массивами произвольной длины (хотя использование таблиц длины = 2 обеспечивает лучшее быстродействие);∙ возможность использования как комплексных, так и вещественных входных / выходных данных;∙ интерфейс на языках Си и Фортран;∙ наличие подробной, хорошо структурированной документации спримерами использования.Указанные достоинства библиотеки FFTW обусловливают её широкое применение.

Использование библиотеки FFTW вместо самописного программного кода в некоторых случаях может дать ускорениевплоть до порядка величины даже без использования параллельныхвычислений в SMP/MPP системах. В данном параграфе мы краткопознакомимся с основными принципами использования FFTW.2.11.1. УстановкаДля установки библиотеки необходимо скачать дистрибутив FFTWсо страницы Download официального сайта [A9].Пользователям операционной системы Microsoft Windows следуетперейти с основной страницы загрузок по ссылке Go here for Windowsи скачать архив с уже готовыми к использованию dll-файлами19 библиотеки FFTW.

Обратим внимание на наличие в архиве сразу трёхdll-файлов: libfftw3-3.dll, libfftw3f-3.dll и libfftw3l-3.dll, содержащих реализацию для типа данных double (файл без суффикса19 Сокращение от, динамически подключаемая библиотекав операционных системах семейства MS Windows.Dynamic Link Library58в имени), float (файл с суффиксом f) и long double (файл с суффиксом l). Требуемый программе dll-файл необходимо разместить врабочем каталоге рядом с создаваемым exe-файлом либо в одном изстандартных каталогов20 , используемых для размещения динамическиподключаемых библиотек (C:\Windows, C:\Windows\System32).Пользователям Microsoft Visual Studio для сборки программы, использующей библиотеку FFTW, необходимо дополнительно создатьlib-файлы с помощью команд21lib /def:libfftw3-3.deflib /def:libfftw3f-3.deflib /def:libfftw3l-3.defДля запуска программы lib.exe нужно открыть консоль разработчика Microsoft Visual Studio (Developer Command Prompt for VS).