将犰狳与 FFTW 接口

Interfacing Armadillo with FFTW

我正在尝试在 Armadillo 包的矩阵中使用 FFTW。我需要在二维矩阵上执行 N 个独立的 FFT。按照 FFTW 和其他在线资源的手册,我有以下代码:

#include <armadillo>
#include <iostream>
#include "fftw3.h"

using namespace std;
using namespace arma;

fftw_plan fplan;

int main(int argc, char *argv[])
{
  cx_mat psi;
  int NCOL=2, NROW=4;
  int frank=1;
  int howmany=NCOL;
  int n[]={NROW};         
  int idist=NROW, odist=NROW; 
  int istride=1, ostride=1;
  int *inembed=n, *onembed=n;

  psi.resize(NROW, NCOL);
  psi(0,0)=cx_double(1,1); psi(0,1)=cx_double(1,1);
  psi(1,0)=cx_double(2,1); psi(1,1)=cx_double(1,2);
  psi(2,0)=cx_double(3,1); psi(2,1)=cx_double(1,3);
  psi(3,0)=cx_double(4,1); psi(3,1)=cx_double(1,4);
  cout << psi << endl; // the output is correct at here

  fftw_complex* in = reinterpret_cast<fftw_complex*>(psi.memptr());

  fplan = fftw_plan_many_dft(frank, n, howmany,
                             in, inembed, istride, idist,
                             in, onembed, ostride, odist,
                             FFTW_FORWARD, FFTW_MEASURE);

  cout << endl << "after fftw plan: " << endl << psi << endl; // all zeros

  // fftw_execute(fplan); // output zeros again if execute
  cx_double* A_mem=psi.memptr();
  cout << endl << "by reference: " << endl << *A_mem << " " << *(A_mem+1) << endl;

  return 0;
}

代码编译没有任何错误。但是,如果我 运行 代码,在 fftw_plan_many_dft 之后,它会输出全零。如果我删除 fftw_plan_many_dft,输出是正确的。我什至不执行 fftw 计划,那么为什么它只通过设置计划来清理我的数据?

使用 FFTW_MEASURE 标志,缓冲区实际上用于调整 FFT 算法。您应该先创建计划,然后初始化输入数据。

请注意 FFTW_ESTIMATE 不会 表现出这种行为,因此您也可以只更改当前代码中的此标志来解决问题,但随后您的代码将可能 运行 慢。