c++ - fftw 拆分示例崩溃

Question

我正在尝试在拆分数组上使用 fftw (3.3.4)。我从 C/MEX 指南 (PDF: http://www.researchgate.net/publictopics.PublicPostFileLoader.html?id=551a2e59cf57d7620c8b463b&key=3a91499f-ec8f-41f4-8e87-fff8bf00e54a ) 中举了一个例子。

程序运行时崩溃(段错误)。运行 valgrind 显示内存违规。

我在这里做错了什么？

#include <vector>
#include <fftw3.h>

void DivideArray(double *Data, int NumEl, double Divisor)
{
  int n;
  for(n = 0; n < NumEl; n++)
    Data[n] /= Divisor;
}

void FFTNDSplit(int NumDims, const int N[], double *XReal, double *XImag, double *YReal, double *YImag, int Sign)
{
  fftw_plan Plan;
  fftw_iodim Dim[NumDims];
  int k, NumEl;
  for(k = 0, NumEl = 1; k < NumDims; k++)
  {
    Dim[NumDims-k-1].n = N[k];
    Dim[NumDims-k-1].is = Dim[NumDims-k-1].os = (k == 0) ? 1 : (N[k-1] * Dim[NumDims-k].is);
    NumEl *= N[k];
  }
  if(!(Plan = fftw_plan_guru_split_dft(NumDims, Dim, 0, NULL,
    XReal, XImag, YReal, YImag, FFTW_ESTIMATE))) {
    printf("Failed creating FFTW plan.\n");
    return;
  }

  if(Sign == -1)
    fftw_execute_split_dft(Plan, XReal, XImag, YReal, YImag);
  else
  {
    fftw_execute_split_dft(Plan, XImag, XReal, YImag, YReal);
    DivideArray(YReal, NumEl, NumEl);
    DivideArray(YImag, NumEl, NumEl);
  }

  fftw_destroy_plan(Plan);
  return; 
}

int main()
{
  int W = 201, H = 201;
  std::vector<double> xr(W*H);
  std::vector<double> xi(W*H);
  std::vector<double> yr(W*H);
  std::vector<double> yi(W*H);

  for (int i = 0; i < W*H; i++) {
    xr[i] = xi[i] = yr[i] = yi[i] = 5.0;
  }

  int dims[2] = { W, H };
  FFTNDSplit(2, dims, &xr[0], &xi[0], &yr[0], &yi[0], 1);
  return 0;
}

valgrind 的第一个错误输出:

==20663== Invalid write of size 8
==20663==    at 0x4E4A0EB: fftw_cpy2d (in /usr/local/lib/libfftw3.so.3.4.4)
==20663==    by 0x4E4A36C: dotile_buf (in /usr/local/lib/libfftw3.so.3.4.4)
==20663==    by 0x4E4F4CA: fftw_tile2d (in /usr/local/lib/libfftw3.so.3.4.4)
==20663==    by 0x4E4F4CA: fftw_tile2d (in /usr/local/lib/libfftw3.so.3.4.4)
==20663==    by 0x4E4A51F: fftw_cpy2d_tiledbuf (in /usr/local/lib/libfftw3.so.3.4.4)
==20663==    by 0x4E8E651: copy (in /usr/local/lib/libfftw3.so.3.4.4)
==20663==    by 0x4E877C6: apply (in /usr/local/lib/libfftw3.so.3.4.4)
==20663==    by 0x4E518F8: apply (in /usr/local/lib/libfftw3.so.3.4.4)
==20663==    by 0x400CDF: FFTNDSplit(int, int const*, double*, double*, double*, double*, int) (in /.../fftw3_split_test)
==20663==    by 0x400FDA: main (in /.../fftw3_split_test)
==20663==  Address 0x5e56b80 is 0 bytes inside a block of size 32 free'd
==20663==    at 0x4C2BCD7: free (vg_replace_malloc.c:473)
==20663==    by 0x4E4F368: fftw_tensor_destroy2 (in /usr/local/lib/libfftw3.so.3.4.4)
==20663==    by 0x4E562BB: fftw_mkproblem_dft_d (in /usr/local/lib/libfftw3.so.3.4.4)
==20663==    by 0x4F1B4CB: fftw_plan_guru_split_dft (in /usr/local/lib/libfftw3.so.3.4.4)
==20663==    by 0x400C7D: FFTNDSplit(int, int const*, double*, double*, double*, double*, int) (in /.../fftw3_split_test)
==20663==    by 0x400FDA: main (in /.../fftw3_split_test)

Answer 1

问题是您计划在数组 (Real, Imag) 上执行，但您正在执行 (Imag, Real)，试图对正向和反向变换使用相同的计划。根据手册第 4.6 节，这是不允许的:

• 对于拆分数组，实部和虚部之间的分离 部分“ii-ri”和“io-ro”与输入相同 和创建计划时的输出数组。 (这个条件是 交错数组自动满足。)

为逆创建一个单独的计划将解决问题。