20 年前我写了一个矩阵计算 C++ 库,我愿意使用英特尔 MKL 库来提高它的性能。对于复值向量/矩阵,我的库使用两个拆分数组:一个用于实部,一个用于虚部。
以下是计时结果:
- N=65536, fftw 时间 = 0.005(s), fftpack 时间 = 0.001(s)
- N=100000, fftw 时间 = 0.005(s), fftpack 时间 = 0.003(s)
- N=131072, fftw 时间 = 0.006(s), fftpack 时间 = 0.004(s)
- N=250000, fftw 时间 = 0.013(s), fftpack 时间 = 0.007(s)
- N=262144,fftw 时间 = 0.012(s),fftpack 时间 = 0.008(s)
- N=524288, fftw 时间 = 0.022(s), fftpack 时间 = 0.018(s)
- N=750000, fftw 时间 = 0.037(s), fftpack 时间 = 0.025(s)
- N=1048576, fftw 时间 = 0.063(s), fftpack 时间 = 0.059(s)
- N=1500000,fftw 时间 = 0.114(s),fftpack 时间 = 0.079(s)
- N=2097152, fftw 时间 = 0.126(s), fftpack 时间 = 0.146(s)
- N=4194304,fftw 时间 = 0.241(s),fftpack 时间 = 0.35(s)
- N=8388608, fftw 时间 = 0.433(s), fftpack 时间 = 0.788(s)
对于长度 < 1500000 的向量,双值 fftpack 比 fftw 快。
这是我使用的代码:
Matrix X=randn(M,1); //input vector
//start timer
Matrix Y = MyFFTW(X);
// measure time
//function to compute the FFT
Matrix MyFFTW(Matrix X)
{
int M= X.rows();
int N= X.cols();
Matrix Y(T_COMPLEX,M,N); // output complex to store FFT results
// Input data could also ba matrix
double* in_data = (double*)fftw_malloc(sizeof(double) * M );
fftw_complex* out_data = (fftw_complex*)fftw_malloc(sizeof(fftw_complex) * (M / 2 + 1));
fftw_plan fftplan = fftw_plan_dft_r2c_1d(M, in_data, out_data, FFTW_ESTIMATE);
//one fftplan is used for all the matrix columns
for (int i = 1; i <= N; i++)
{
//copy column number i to in_dataused by the fftplan, arrays indexing is 1-based like matlab
memcpy(in_data, X.pr(1,i), M* sizeof(double));
fftw_execute(fftplan);
//split out_data to real and imag parts
double* pr = Y.pr(1,i), * pi = Y.pi(1,i);
int k = (M - 1) / 2, j;
for (j = 0; j <= k; j++)
{
*pr++ = out_data[j][0];
*pi++ = out_data[j][1];
}
if (M % 2 == 0)
{
*pr++ = out_data[M/2][0];
*pi++ = out_data[M/2][1];
}
for (j = k; j >= 1; j--)
{
*pr++ = out_data[j][0];
*pi++ = out_data[j][1];
}
}
fftw_destroy_plan(fftplan);
fftw_free(in_data);
fftw_free(out_data);
return Y;
}
结果是在 Intel core i7 @ 3.2 GHz 上使用 Visual Studio 2019 作为编译器和最后一个 Intel MKL 库获得的。编译器标志是:
/fp:fast /DWIN32 /O2 /Ot /Oi /Oy /arch:AVX2 /openmp /MD
链接器库是:
mkl_intel_c.lib mkl_intel_thread.lib mkl_core.lib libiomp5md.lib
有没有更好的方法让小尺寸矢量的 fftw 更快?
更新:
我针对使用 MKL fftw 进行 fft 计算的 Matlab 进行了测试:
- N=65536,matlab fft 时间 = 0.071233(s)
- N=100000,matlab fft 时间 = 0.011437(s)
- N=131072,matlab fft 时间 = 0.0074411(s)
- N=250000,matlab fft 时间 = 0.015349(s)
- N=262144,matlab fft 时间 = 0.0082545(s)
- N=524288,matlab fft 时间 = 0.011395(s)
- N=750000,matlab fft 时间 = 0.022364(s)
- N=1048576,matlab fft 时间 = 0.019683(s)
- N=1500000,matlab fft 时间 = 0.033493(s)
- N=2097152,matlab fft 时间 = 0.035345(s)
- N=4194304,matlab fft 时间 = 0.069539(s)
- N=8388608,matlab fft 时间 = 0.1387(s)
除了第一次使用 N=65536 调用 fft 之外,Matlab(64 位)比我使用 fftpack(对于 N > 500000)和使用 MKL fftw 的函数(win32)更快。
谢谢