Used astyle --style=kr formatted source codes.
[charm.git] / examples / charm++ / nodeHelper / fft-trans / fft_bench.cpp
1 //#include <stdio.h>
2 //#include <iostream>
3
4 #include <mpi.h>
5 #include <fftw3.h>
6 #include <cstdlib>
7 #include <limits>
8 #include <cmath>
9
10 using namespace std;
11
12 int main(int argc, char *argv[]) 
13 {
14   MPI_Init(&argc, &argv);
15
16   long N = atol(argv[1]);
17   long n = N * N;
18
19   fftw_complex *data;
20   fftw_plan p, q;
21
22   data = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * n);
23
24   int length[] = {N};
25   p = fftw_plan_many_dft(1, length, N, data, length, 1, N,
26       data, length, 1, N, FFTW_FORWARD, FFTW_ESTIMATE); 
27   q = fftw_plan_many_dft(1, length, N, data, length, 1, N,
28       data, length, 1, N, FFTW_BACKWARD, FFTW_ESTIMATE);
29
30   srand48(0);
31   for(int i = 0; i<n; i++) {
32     data[i][0] = drand48();
33     data[i][1] = drand48();
34   }
35
36   double t1, t2; 
37   t1 = MPI_Wtime(); 
38
39   fftw_execute(p);
40
41   t2 = MPI_Wtime() - t1; 
42   printf( "On %ld elements, Elapsed time is %f with %f Gflop/s\n", n, t2, 5*(double)n*log2((double)n)/(t2*1000000000)); 
43
44   fftw_execute(q);
45
46   double infNorm = 0.0;
47
48   srand48(0);
49   for(int i=0; i<n; i++) {
50     data[i][0] = data[i][0]/N - drand48();
51     data[i][1] = data[i][1]/N - drand48();
52
53     double mag = sqrt(pow(data[i][0],2) + pow(data[i][1], 2));
54     if(mag > infNorm) infNorm = mag;
55   } 
56
57   double r = infNorm / (std::numeric_limits<double>::epsilon() * log((double)n));
58
59   printf("residual = %f\n", r);
60
61   fftw_destroy_plan(p);
62   fftw_destroy_plan(q);
63   fftw_free(data);
64
65   MPI_Finalize();
66   return 0;
67 }