title: Wrong results cufft 3D in-place
body: <p>I write because I'm facing problems with the cufft 3D transform in-place, while I have no problems for the out-of-place version. I tried to follow Robert Crovella's answer <a href="https://stackoverflow.com/questions/24809179/cufft-double-to-complex">here</a> but I'm not obtaining the correct results when I make a FFT+IFT.
  This is my code:</p>
  
  <pre><code>#include &lt;stdio.h&gt;
  #include &lt;stdlib.h&gt;
  #include &lt;cuda_runtime.h&gt;
  #include &lt;complex.h&gt;
  #include &lt;cuComplex.h&gt;
  #include &lt;cufft.h&gt;
  
  // Main function
  int main(int argc, char **argv){
      int N = 4; 
      double *in = NULL, *d_in = NULL;
      cuDoubleComplex *out = NULL, *d_out = NULL;
      cufftHandle plan_r2c, plan_c2r;
  
      unsigned int out_mem_size = sizeof(cuDoubleComplex) * N*N*(N/2 + 1);
      unsigned int in_mem_size = out_mem_size;
  
      in  = (double *) malloc (in_mem_size);
      out  = (cuDoubleComplex *)in;
  
      cudaMalloc((void **)&amp;d_in, in_mem_size);
      d_out = (cuDoubleComplex *)d_in;
  
      cufftPlan3d(&amp;plan_r2c, N, N, N, CUFFT_D2Z);
      cufftPlan3d(&amp;plan_c2r, N, N, N, CUFFT_Z2D);
  
      memset(in, 0, in_mem_size);
      unsigned int idx;
      for (int z = 0; z &lt; N; z++){
          for (int y = 0; y &lt; N; y++){
              for (int x = 0; x &lt; N; x++){
                  idx = z + N * ( y + x * N);
                  in[idx] = idx;
              }
          }
      }
      printf("\nStart: \n");
      for (int z = 0; z &lt; N; z++){
          printf("plane = %d ----------------------------\n", z);
          for (int x = 0; x &lt; N; x++){
              for (int y = 0; y &lt; N; y++){
                  idx = z + N * ( y + x * N);
                  printf("%.3f \t", in[idx]);
              }
              printf("\n");
          }
      }
      cudaMemcpy(d_in, in, in_mem_size, cudaMemcpyHostToDevice);
  
      cufftExecD2Z(plan_r2c, (cufftDoubleReal *)d_in, (cufftDoubleComplex *)d_out);
      cufftExecZ2D(plan_c2r, (cufftDoubleComplex *)d_out, (cufftDoubleReal *)d_in);
  
      memset(in, 0, in_mem_size);
      CU_ERR_CHECK( cudaMemcpy(in, d_in, in_mem_size, cudaMemcpyDeviceToHost) );
  
      printf("\nAfter FFT+IFT: \n");
      for (int z = 0; z &lt; N; z++){
          printf("plane = %d ----------------------------\n", z);
          for (int x = 0; x &lt; N; x++){
              for (int y = 0; y &lt; N; y++){
                  idx = z + N * ( y + x * N);
                  // Normalisation
                  in[idx] /= (N*N*N);
                  printf("%.3f \t", in[idx]);
              }
              printf("\n");
          }
      }
  
      return 0;
  }
  </code></pre>
  
  <p>The program outputs the following data:</p>
  
  <p>Starting file </p>
  
  <p>plane = 0 ----------------------------</p>
  
  <p>0.000   4.000   8.000   12.000<br>
  16.000  20.000  24.000  28.000<br>
  32.000  36.000  40.000  44.000<br>
  48.000  52.000  56.000  60.000  </p>
  
  <p>plane = 1 ----------------------------</p>
  
  <p>1.000   5.000   9.000   13.000<br>
  17.000  21.000  25.000  29.000<br>
  33.000  37.000  41.000  45.000<br>
  49.000  53.000  57.000  61.000  </p>
  
  <p>plane = 2 ----------------------------</p>
  
  <p>2.000   6.000   10.000  14.000<br>
  18.000  22.000  26.000  30.000<br>
  34.000  38.000  42.000  46.000<br>
  50.000  54.000  58.000  62.000  </p>
  
  <p>plane = 3 ----------------------------</p>
  
  <p>3.000   7.000   11.000  15.000<br>
  19.000  23.000  27.000  31.000<br>
  35.000  39.000  43.000  47.000<br>
  51.000  55.000  59.000  63.000  </p>
  
  <p>After FFT+IFT </p>
  
  <p>plane = 0 ----------------------------</p>
  
  <p>-0.000  -0.344  8.000   12.000<br>
  -0.031  20.000  24.000  -0.031<br>
  32.000  36.000  0.031   44.000<br>
  48.000  -0.094  56.000  60.000  </p>
  
  <p>plane = 1 ----------------------------</p>
  
  <p>1.000   -0.000  9.000   13.000<br>
  -0.000  21.000  25.000  0.125<br>
  33.000  37.000  0.000   45.000<br>
  49.000  0.000   57.000  61.000  </p>
  
  <p>plane = 2 ----------------------------</p>
  
  <p>2.000   6.000   -0.000  14.000<br>
  18.000  0.000   26.000  30.000<br>
  0.000   38.000  42.000  -0.000<br>
  50.000  54.000  -0.000  62.000  </p>
  
  <p>plane = 3 ----------------------------</p>
  
  <p>3.000   7.000   0.031   15.000<br>
  19.000  -0.031  27.000  31.000<br>
  -0.031  39.000  43.000  0.031<br>
  51.000  55.000  0.031   63.000  </p>
  
  <p>I even tried to pad the data this way:</p>
  
  <pre><code>// With padding
      unsigned int idx;
      for (int x = 0; x &lt; N; x++){
          for (int y = 0; y &lt; N; y++){
              for (int z = 0; z &lt; 2*(N/2+1); z++){
                  idx = z + N * ( y + x * N);
                  if (z &lt; 4) in[idx] = idx;
                  else in[idx] = 0;
              }
          }
      }
  </code></pre>
  
  <p>What am I doing wrong?</p>