1

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)更快。

谢谢

4

1 回答 1

0

关于 fftw、AFAIK,MKL 没有具体的性能提示可以帮助加快小型案例的性能。实际上,从 mkl 使用 fftw 的开销可以忽略不计。
关于您的工作台:我看到您测量分配/解除分配部分,创建 fftw 计划,以及内存复制操作。但是,该基准测试中唯一的一个例程(fftw_execute)由 mkl 优化。这可能是这个管道的问题。您可以添加 MKL_VERBOSE 模式来检查 fftw_execute 的执行时间...

于 2020-07-21T06:44:50.887 回答