examples/fdtd-1d/fdtd-1d.c

   1 #include <stdio.h>
   2 #include <stdlib.h>
   3 #include <math.h>
   4 #include <unistd.h>
   5 #include <sys/time.h>
   6 #include <assert.h>
   7
   8 #define N 1000000
   9 #define T 10000
  10 double h[N];
  11 double e[N+1];
  12 #define coeff1 0.5
  13 #define coeff2 0.7
  14
  15 #ifdef TIME
  16 #define IF_TIME(foo) foo;
  17 #else
  18 #define IF_TIME(foo)
  19 #endif
  20
  21 void init_array()
  22 {
  23     int i, j;
  24
  25         for (j=0; j<N; j++) {
  26             h[j] = ((double)j)/N;
  27             e[j] = ((double)j)/N;
  28         }
  29 }
  30
  31 void print_array()
  32 {
  33     int i, j;
  34
  35     for (j=0; j<N; j++) {
  36             fprintf(stdout, "%lf ", h[j]);
  37             if (j%80 == 79) fprintf(stdout, "\n");
  38     }
  39     fprintf(stdout, "\n");
  40 }
  41
  42 double rtclock()
  43 {
  44     struct timezone Tzp;
  45     struct timeval Tp;
  46     int stat;
  47     stat = gettimeofday (&Tp, &Tzp);
  48     if (stat != 0) printf("Error return from gettimeofday: %d",stat);
  49     return(Tp.tv_sec + Tp.tv_usec*1.0e-6);
  50 }
  51
  52 int main()
  53 {
  54     int t, i, j, k, l;
  55
  56     double t_start, t_end;
  57
  58     init_array();
  59
  60         IF_TIME(t_start = rtclock());
  61
  62 #pragma scop
  63     for (t=1; t<=T; t++){
  64             for (i=1; i<=N-1; i++)
  65                     e[i] = e[i] - coeff1*(h[i]-h[i-1]);
  66             for (i=0; i<=N-1; i++)
  67                     h[i] = h[i] - coeff2*(e[i+1]-e[i]);
  68     }
  69 #pragma endscop
  70
  71     IF_TIME(t_end = rtclock());
  72     IF_TIME(fprintf(stderr, "%0.6lfs\n", t_end - t_start));
  73
  74     if (fopen(".test", "r")) {
  75         print_array();
  76     }
  77
  78     return 0;
  79 }