Итак, я пытаюсь, как и многие другие до меня, просеять все небесные выравнивания фазовых сдвигов и коэффициентов нормализации, необходимые для дифференциации с помощью преобразования Фурье.
Я пытаюсь использовать как можно меньше кода, полностью используя функциональность numpy для поддержания чистоты кода при использовании его операций с массивами, которые, как мне сказали, быстрее, чем явные циклы python.
Из строки 106 таблицы Функциональные отношения, одномерные , в разделе о Важных преобразованиях Фурье статьи Википедии о преобразованиях Фурье, я понял, что преобразование Фурье обладает этим свойством, когда дифференцирование в реальном пространстве эквивалентно умножению в частотном пространстве, так что d^nf/dx^n = ifft[F*(i*w)^n]
, где F
- это преобразование Фурье f
.
В поисках более надежной ссылки для определения того, что будет w
в d^nf/dx^n = ifft[F*(i*w)^n]
, я нашел документ PDF по математике дифференцирования с использованием ДПФ. В статье автор обсуждает, что существует бесконечно много так называемых тригонометрических интерполяций для данной функции, и показывает способ получения «менее волнистого» решения, которое является уникальным.
Итак, я использовал его определения и написал небольшой код на Python для проведения различий. Вот код, который я придумал (обратите внимание на последние две строки, так как они изменятся)
N = 51
a = 0
b = 2*np.pi
X = np.linspace(a,b, N)
L = abs(b-a)
TwoPiOverL = 2.0*np.pi/L
freqs = np.fft.fftfreq(N,1./N)
# THESE COME FROM [1] #
D1 = TwoPiOverL*freqs*1j
D2 = -(TwoPiOverL*freqs)**2
#=========================#
S = np.sin(X)
fftS=np.fft.fft(S)
C = np.cos(X)
# COMPUTING THE DERIVATIVES #
Y1 = np.fft.ifft(D1*fftS).real
Y2 = np.fft.ifft(D2*fftS).real
На рисунках сгенерированных графиков запись ~d/dx sin
означает производную функции синуса, вычисленную согласно умножению с преобразованием Фурье, вместо аналитического выражения.
Тест 1
Первый тест, используя только определения, найденные в PDF:
Итак, первый тест был беспорядок. Производные хорошо интерполируются вокруг середины области, но конечные точки сумасшедшие. В самом PDF-документе в начале раздела 6 автор пишет:
Все вышеизложенное относится к спектральному дифференцированию периодических функций y (x). Если функция непериодическая,
артефакты появятся в производных из-за неявных разрывов на конечных точках.
Звучит правдоподобно, и это означает, что я где-то напортачу, потому что, если вся статья написана с обсуждением периодических сигналов, и если мой сигнал попадает в непериодический класс, то есть место для улучшения, если я можно использовать периодический сигнал.
Но я использую DFT, что означает, что я применяю числовую процедуру к массиву. Дискретные данные, которые не несут никакой информации о периодичности. И тогда я подумал:
Ну, как я могу сказать своему дифференциатору Фурье, что сигнал периодический? Что это влечет за собой в коде? Где, в коде, я должен измениться, чтобы передать эту информацию в мою программу?
То, что казалось ответом, пришло из стека переполнения. В этом посте о численном дифференцировании с использованием БПФ (да, избыточно, но я туда доберусь), @dkato дает инструкции, которые я интерпретировал как изменение двух последних строк в моем коде с
Y1 = np.fft.ifft(D1*fftS).real
Y2 = np.fft.ifft(D2*fftS).real
до
Y1 = np.fft.ifft(np.exp(D1)*fftS).real
Y2 = np.fft.ifft(np.exp(D2)*fftS).real
Итак, я сделал это, а затем проверил, что я получу.
Тест 2
Второй тест с использованием поправок @dkato:
Как видно из рисунка выше, при использовании поправок, предоставленных @dkato, покачивание конечной точки исчезает, но оценка производной все еще довольно далека. Я использую 51 точку данных, поэтому я предполагаю, что это не связано с неточностями в расширении. Кажется, что способ вычисления производных неверен.
Источником наибольшей потери точности между первой аналитической производной и ее оценкой преобразования Фурье, по-видимому, является сдвиг вправо, тогда как вторая производная не только смещена, но и плохо масштабируется. Довольно забавно, однако, когда мы используем
Y2 = np.fft.ifft(-np.exp(D2)*fftS).real
для оценки второй производной мы видим, что вся ошибка из-за смещения исчезает, оставляя только масштабирование .Я не знаю, является ли это свойством общих вторичных производных, работающих на БПФ, или это происходит только потому, что я использую функцию синуса в качестве своего y (x), в частности.
The Point, наконец
Я рад, что странные, уродливые покачивания на конечных точках исчезли, а это означает, что я должен был каким-то образом сказать своему Фурье-Дифференциатору, что работаю с периодическими сигналами.Но эти странные сдвиги и вычисления моих дифференциальных оценок все еще надоедливы, и я не могу на всю жизнь положить палец на то, что не так.Как вы могли догадаться, отчасти потому, что у меня мало опыта игры со всеми движущимися частями в DFT.
В этом качестве у меня есть два вопроса к сообществу, если вы будете так любезны, чтобы помочь мнес этой проблемой.
1) Если это правда, что я как-то сказал своему дифференциатору Фурье, что я использую периодическую функцию, как именно я это сделал?В чем причина кода, который работает, и кода, который не работает?
и
2) Что является источником этих странных несоответствий масштаба / сдвига междуаналитические и оценочные производные?Как я могу их исправить?