在 Python 中使用 DFT 求一阶导数

Finding first derivative using DFT in Python

我想使用离散傅立叶变换在区间 [0, 2/pi] 上找到 exp(sin(x)) 的一阶导数。基本思想是首先在给定区间上计算 exp(sin(x)) 的 DFT,给出 v_k,然后计算 ikv_k 的逆 DFT,给出所需的答案。实际上,由于傅里叶变换在编程语言中的实现,您可能需要在某处重新排序输出 and/or 乘以不同的因子。

我首先在 Mathematica 中进行了此操作,其中有一个选项 FourierParameters,可让您指定转换的约定。首先,我获得了一个高斯分布的傅立叶级数,以便了解我必须乘以哪些归一化因子,然后继续求导数。不幸的是,此后将我的 Mathematica 代码翻译成 Python(因此我又第一次做了高斯的傅里叶级数——这是成功的),我没有得到相同的结果。这是我的代码:

N=1000
xmin=0
xmax=2.0*np.pi
step = (xmax-xmin)/(N)
xdata = np.linspace(xmin, xmax-step, N)
v = np.exp(np.sin(xdata))
derv = np.cos(xdata)*v
vhat = np.fft.fft(v)
kvals1 = np.arange(0, N/2.0, 1)
kvals2 = np.arange(-N/2.0, 0, 1)
what1 = np.zeros(kvals1.size+1)
what2 = np.empty(kvals2.size)
it = np.nditer(kvals1, flags=['f_index'])
while not it.finished:
    np.put(what1, it.index, 1j*(2.0*np.pi)/((xmax-xmin))*it[0]*vhat[[int(it[0])]])
    it.iternext()
it = np.nditer(kvals2, flags=['f_index'])
while not it.finished:
    np.put(what2, it.index, 1j*(2.0*np.pi)/((xmax-xmin))*it[0]*vhat[[int(it[0])]])
    it.iternext()
xdatafull = np.concatenate((xdata, [2.0*np.pi]))
what = np.concatenate((what1, what2))
w = np.real(np.fft.ifft(what))

fig = plt.figure()
ax = plt.gca()
ax.spines['right'].set_color('none')
ax.spines['top'].set_color('none')
ax.xaxis.set_ticks_position('bottom')
ax.spines['bottom'].set_position(('data',0))
ax.yaxis.set_ticks_position('left')
ax.spines['left'].set_position(('data',0))

plt.plot(xdata, derv, color='blue')
plt.plot(xdatafull, w, color='red')
plt.show()

我可以 post Mathematica 代码,如果有人要我的话。

原来问题是 np.zeros 给了你一个实数数组而不是复数数组,因此之后的赋值不会改变任何东西,因为它们是虚构的。

因此解决方案非常简单

import numpy as np

N=100
xmin=0
xmax=2.0*np.pi
step = (xmax-xmin)/(N)
xdata = np.linspace(step, xmax, N)
v = np.exp(np.sin(xdata))
derv = np.cos(xdata)*v
vhat = np.fft.fft(v)

what = 1j*np.zeros(N)
what[0:N/2.0] = 1j*np.arange(0, N/2.0, 1)
what[N/2+1:] = 1j*np.arange(-N/2.0 + 1, 0, 1)
what = what*vhat
w = np.real(np.fft.ifft(what))

# Then plotting

np.zeros1j*np.zeros

取代