离散傅里叶变换，加速计算

Question

输入： hermitian矩阵\rho_{i,j}和i,j=0,1,..d-1

输出： neg=\sum |W(mu,m)|-W(mu,m)，所有mu,m=0,..d-1和W(mu,m)=\sum exp(-4i\pi mu n /d) \rho_{(m-n)%d,(m+n)%d}之和，其中n=0,..d-1

问题： 1)对于大型d (d>5 000)，直接方法(参见代码片段1)相当慢。

2)使用'np.fft.fft()‘要快得多，但在定义中使用的指数是2，而不是4 https://docs.scipy.org/doc/numpy/reference/routines.fft.html#module-numpy.fft

是否有可能使用代码片段2来改进代码段1以提高计算速度？可以用二维fft吗？

片段1:

W=np.zeros([d,d])
neg=0
        for mu in range(d):
            for m in range(d):
                x=0            
                for n in range(d): 
                    x+=np.exp(-4*np.pi*1.0j*mu*n/N)*rho[(m-n)%d,(m+n)%d]
                W[mu,m]=x.real
                neg+=np.abs(W[mu,m])-W[mu,m]

代码片段2:

# create matrix \rho
psi=np.random.rand(500)
psi=psi/np.linalg.norm(psi) # normalize it
d=len(psi) 
rho=np.outer(psi,np.conj(psi)) 

#
start_time=time.time()
m=1 # for example, use particular m 

a=np.array([rho[(m-nn)%d,(m+nn)%d] for nn in range(d)])
ft=np.fft.fft(a)

end_time=time.time()
print(end_time-start_time)

Answer 1

通过使用numpy的数组算法删除嵌套循环。

import numpy as np

def my_sins(x):
    return np.sin(2*x) + np.cos(4*x) + np.sin(3*x)

def dft(x, n=None):
    if n is None:
        n = len(x)   
    k = len(x)    
    cn = np.sum(x*np.exp(-2*np.pi*1j*np.outer(np.arange(k),np.arange(n))/n),1)
    return cn

对于一些样本数据

x = np.linspace(0,2*np.pi,1000)
y = my_sins(x)

%timeit dft(y)

在我的系统中，这会产生：

145 ms ± 953 µs per loop (mean ± std. dev. of 7 runs, 10 loops each)