0

CuFFT Lib で double から std::complex への FFT を作成したいと考えています。私のコードは次のようになります

#include <complex>
#include <iostream>
#include <cufft.h>
#include <cuda_runtime_api.h>

typedef std::complex<double> Complex;
using namespace std;

int main(){
  int n = 100;
  double* in;
  Complex* out;
  in = (double*) malloc(sizeof(double) * n);
  out = (Complex*) malloc(sizeof(Complex) * n/2+1);
  for(int i=0; i<n; i++){
     in[i] = 1;
  }

  cufftHandle plan;
  plan = cufftPlan1d(&plan, n, CUFFT_D2Z, 1);
  unsigned int mem_size = sizeof(double)*n;
  cufftDoubleReal *d_in;
  cufftDoubleComplex *d_out;
  cudaMalloc((void **)&d_in, mem_size);
  cudaMalloc((void **)&d_out, mem_size);
  cudaMemcpy(d_in, in, mem_size, cudaMemcpyHostToDevice);
  cudaMemcpy(d_out, out, mem_size, cudaMemcpyHostToDevice);
  int succes = cufftExecD2Z(plan,(cufftDoubleReal *) d_in,(cufftDoubleComplex *) d_out);
  cout << succes << endl;
  cudaMemcpy(out, d_out, mem_size, cudaMemcpyDeviceToHost);

  for(int i=0; i<n/2; i++){
     cout << "out: " << i << " "  << out[i].real() << " " <<  out[i].imag() << endl;
  }
  return 0;
}

しかし、変換された値は 1 0 0 0 0 .... または正規化なしで 100 0 0 0 0 .... である必要があると思うので、これは間違っているように思えますが、 0 0 0 0 0 を取得するだけです。 ..

さらに、cufftExecD2Z が適切に機能する場合は、それが可能である必要がありますが、正しく行う方法がわかりません。誰でも助けることができますか？

4

1 に答える 1

1

コードにさまざまなエラーがあります。おそらく、カフトのドキュメントとサンプル コードを確認する必要があります。

  1. すべての API 戻り値に対して適切な cuda エラー チェックと適切な cufft エラー チェックを行う必要があります。

  2. 関数の戻り値はcufftPlan1d計画に入りません。

    plan = cufftPlan1d(&plan, n, CUFFT_D2Z, 1);

    関数自体がプランを設定します (その&planため、関数に渡すのです)。戻り値をプランに代入すると、関数によって設定されたプランが台無しになります。

  3. 出力が size になる可能性があることを正しく識別しましたが((N/2)+1)、ホスト側でも適切にスペースを割り当てませんでした:

    out = (Complex*) malloc(sizeof(Complex) * n/2+1);

    またはデバイス側で：

    unsigned int mem_size = sizeof(double)*n;
...
cudaMalloc((void **)&d_out, mem_size);

次のコードでは、上記の問題のいくつかが修正されており、目的の結果 (100、0、0、...) を得るのに十分です。

#include <complex>
#include <iostream>
#include <cufft.h>
#include <cuda_runtime_api.h>

#define cudaCheckErrors(msg) \
    do { \
        cudaError_t __err = cudaGetLastError(); \
        if (__err != cudaSuccess) { \
            fprintf(stderr, "Fatal error: %s (%s at %s:%d)\n", \
                msg, cudaGetErrorString(__err), \
                __FILE__, __LINE__); \
            fprintf(stderr, "*** FAILED - ABORTING\n"); \
            exit(1); \
        } \
    } while (0)


typedef std::complex<double> Complex;
using namespace std;

int main(){
  int n = 100;
  double* in;
  Complex* out;
#ifdef IN_PLACE
  in = (double*) malloc(sizeof(Complex) * (n/2+1));
  out = (Complex*)in;
#else
  in = (double*) malloc(sizeof(double) * n);
  out = (Complex*) malloc(sizeof(Complex) * (n/2+1));
#endif
  for(int i=0; i<n; i++){
     in[i] = 1;
  }

  cufftHandle plan;
  cufftResult res = cufftPlan1d(&plan, n, CUFFT_D2Z, 1);
  if (res != CUFFT_SUCCESS)  {cout << "cufft plan error: " << res << endl; return 1;}
  cufftDoubleReal *d_in;
  cufftDoubleComplex *d_out;
  unsigned int out_mem_size = (n/2 + 1)*sizeof(cufftDoubleComplex);
#ifdef IN_PLACE
  unsigned int in_mem_size = out_mem_size;
  cudaMalloc((void **)&d_in, in_mem_size);
  d_out = (cufftDoubleComplex *)d_in;
#else
  unsigned int in_mem_size = sizeof(cufftDoubleReal)*n;
  cudaMalloc((void **)&d_in, in_mem_size);
  cudaMalloc((void **)&d_out, out_mem_size);
#endif
  cudaCheckErrors("cuda malloc fail");
  cudaMemcpy(d_in, in, in_mem_size, cudaMemcpyHostToDevice);
  cudaCheckErrors("cuda memcpy H2D fail");
  res = cufftExecD2Z(plan,d_in, d_out);
  if (res != CUFFT_SUCCESS)  {cout << "cufft exec error: " << res << endl; return 1;}
  cudaMemcpy(out, d_out, out_mem_size, cudaMemcpyDeviceToHost);
  cudaCheckErrors("cuda memcpy D2H fail");

  for(int i=0; i<n/2; i++){
     cout << "out: " << i << " "  << out[i].real() << " " <<  out[i].imag() << endl;
  }
  return 0;
}

実際のケースから複雑なケースへのインプレース変換に必要なものについては、ドキュメントを参照してください。上記のコードを再コンパイルし-DIN_PLACEて、インプレース変換の動作を確認し、必要なコードを変更することができます。

于 2014-07-17T17:11:08.553 に答える