/*
 * Thrashing demo using a fast fourier transform (fftw3 library)
 * on a macbook pro intel core i7, 2.6GHz, 8GB memory, SolidState hard disk
 *
 *  in another window :
 *  $ top -o cpu -stats pid,command,cpu,time,csw,rsize,vsize,faults
 *
 *  compiling and running it :
 *  $ gcc test_fftw.c -lfftw3 -o test_fftw
 *  $ ./test_fftw
 *
 *    input & output arrays of size 1048576 needs  0.03355 GB.
 *    Setting up fft with N=1048576;executing fft...
 *    done; elapsed time =    0.135 sec.
 *    
 *    input & output arrays of size 2097152 needs  0.06711 GB.
 *    Setting up fft with N=2097152;executing fft...
 *    done; elapsed time =    0.251 sec.
 *    
 *    input & output arrays of size 4194304 needs  0.13422 GB.
 *    Setting up fft with N=4194304;executing fft...
 *    done; elapsed time =    0.553 sec.
 *    
 *    input & output arrays of size 8388608 needs  0.26844 GB.
 *    Setting up fft with N=8388608;executing fft...
 *    done; elapsed time =    1.169 sec.
 *    
 *    input & output arrays of size 16777216 needs  0.53687 GB.
 *    Setting up fft with N=16777216;executing fft...
 *    done; elapsed time =    2.384 sec.
 *    
 *    input & output arrays of size 33554432 needs  1.07374 GB.
 *    Setting up fft with N=33554432;executing fft...
 *    done; elapsed time =    4.905 sec.
 *    
 *    input & output arrays of size 67108864 needs  2.14748 GB.
 *    Setting up fft with N=67108864;executing fft...
 *    done; elapsed time =   25.126 sec.
 *    
 *    input & output arrays of size 134217728 needs  4.29497 GB.
 *    Setting up fft with N=134217728;executing fft...
 *    done; elapsed time =   72.185 sec.
 *
 * We're looking at twice the array sizes in each run.
 * For the first few, execution time also doubles at each step.
 * But in the step from 1GB to 2GB of meory, execution time
 * goes up 5 times .. and then 3 times ; in both cases
 * the page fault counts in 'top' go way up, and cpu usage drops,
 * since the data flow can't keep the processor busy.
 *
 * I'm surprised it does even this well;
 * it may be that fftw3 is being very clever.
 *
 */

#include <time.h>
#include <fftw3.h>
#include <stdio.h>

int main(){
  clock_t start, end;
  fftw_complex *in, *out;
  fftw_plan p;
  int N = 1048576;  // 2**20
  int i, loop;
  float size;

  for (loop=0; loop<8; loop++){

    start = clock();

    // printf("Each complex number is %li bytes \n", sizeof(double)*2);
    size = sizeof(double)*2*2*N / 1.0e9;
    printf("input & output arrays of size %i needs %8.5f GB.\n", 
	   N, size);
    printf("Setting up fft with N=%i;", N);
    // i'th real and complex parts are e.g. in[i][0] and in[i][1]
    in = (fftw_complex*) fftw_malloc(sizeof(fftw_complex)*N);
    for (i=0; i<N; i++){
      in[i][0] = 0.0;
      in[i][1] = 0;
    }
    in[N/2][0] = 1.0;
    in[N/2 - 1][0] = 1.0;
    in[N/2 + 1][0] = 1.0;
    out = (fftw_complex*) fftw_malloc(sizeof(fftw_complex)*N);
    // This sets up the fft without actually doing it yet.
    p = fftw_plan_dft_1d(N, in, out, FFTW_FORWARD, FFTW_ESTIMATE);
    // And this runs it.
    printf("executing fft...\n");
    fftw_execute(p);
    end = clock();
    printf("done; elapsed time = %8.3f sec.\n", ((double)(end-start))/CLOCKS_PER_SEC);
    // for (i=0; i<100; i++){
    //   printf("  %8.4f + %8.4f i \n", out[i][0], out[i][1]);
    // }

    printf("\n");
    N = N * 2;
  }

  return 0;
}