home *** CD-ROM | disk | FTP | other *** search
/ SGI Developer Toolbox 6.1 / SGI Developer Toolbox 6.1 - Disc 4.iso / lib / mathlib / libconv / TRY / speed1d.c < prev    next >
Encoding:
C/C++ Source or Header  |  1994-08-02  |  11.2 KB  |  499 lines
  1. #include <stdio.h>
  2. #include <sys/time.h>
  3. #include <math.h>
  4. #include "conv.h"
  5.  
  6. /* --------- The following definitions change 
  7.         according to precision required -------- */
  8.  
  9.  
  10. #ifdef DOUBLE            /* real double precision    */
  11.  
  12. typedef double this_type;
  13.  
  14. #define OPS_PER_ITER    2
  15.  
  16. void dfir1d_(), ornl_dfir1d_();
  17. void diir1d_(), ornl_diir1d_();
  18. void dcor1d_(), ornl_dcor1d_();
  19. void simple_dfir1d_(), simple_diir1d_(), simple_dcor1d_();
  20. #define SIMPLE_FIR        simple_dfir1d_
  21. #define SIMPLE_IIR        simple_diir1d_
  22. #define SIMPLE_COR        simple_dcor1d_
  23. #define ORNL_FIR        ornl_dfir1d_
  24. #define ORNL_IIR        ornl_diir1d_
  25. #define ORNL_COR        ornl_dcor1d_
  26. #define MY_COR            dcor1d_
  27. #define MY_FIR            dfir1d_
  28. #define MY_IIR            diir1d_
  29.  
  30. #define MY_INIT            dinit_
  31. #define MY_ONE            done_
  32. #define THIS_REAL
  33. this_type zero = 0.;
  34. this_type one  = 1.;
  35. #endif
  36. #ifdef SINGLE            /* real single precision    */
  37.  
  38. typedef float this_type;
  39.  
  40. #define OPS_PER_ITER    2
  41.  
  42. void sfir1d_(), ornl_sfir1d_();
  43. void siir1d_(), ornl_siir1d_();
  44. void scor1d_(), ornl_scor1d_();
  45. void simple_sfir1d_(), simple_siir1d_(), simple_scor1d_();
  46. #define SIMPLE_FIR        simple_sfir1d_
  47. #define SIMPLE_IIR        simple_siir1d_
  48. #define SIMPLE_COR        simple_scor1d_
  49. #define ORNL_FIR        ornl_sfir1d_
  50. #define ORNL_IIR        ornl_siir1d_
  51. #define ORNL_COR        ornl_scor1d_
  52. #define MY_COR            scor1d_
  53. #define MY_FIR            sfir1d_
  54. #define MY_IIR            siir1d_
  55. #define MY_INIT            sinit_
  56. #define MY_ONE            sone_
  57. #define THIS_REAL
  58. this_type zero = 0.;
  59. this_type one  = 1.;
  60. #endif
  61.  
  62. #ifdef ZOMPLEX            /* complex double precision    */
  63.  
  64. typedef struct {double real, imag;} this_type;
  65.  
  66. #define OPS_PER_ITER    8
  67.  
  68. void zfir1d_(), ornl_zfir1d_();
  69. void ziir1d_(), ornl_ziir1d_();
  70. void zcor1d_(), ornl_zcor1d_();
  71. void simple_zfir1d_(), simple_ziir1d_(), simple_zcor1d_();
  72. #define SIMPLE_FIR        simple_zfir1d_
  73. #define SIMPLE_IIR        simple_ziir1d_
  74. #define SIMPLE_COR        simple_zcor1d_
  75. #define ORNL_FIR        ornl_zfir1d_
  76. #define ORNL_IIR        ornl_ziir1d_
  77. #define ORNL_COR        ornl_zcor1d_
  78. #define MY_COR            zcor1d_
  79. #define MY_FIR            zfir1d_
  80. #define MY_IIR            ziir1d_
  81. #define MY_INIT            zinit_
  82. #define MY_ONE            zone_
  83. #define THIS_COMPLEX
  84. this_type zero = {    0., 0.};
  85. this_type one  = {    1., 0.};
  86. #endif
  87. #ifdef COMPLEX            /* complex single precision    */
  88.  
  89. typedef struct {float real, imag;} this_type;
  90.  
  91. #define OPS_PER_ITER    8
  92.  
  93. void cfir1d_(), ornl_cfir1d_();
  94. void ciir1d_(), ornl_ciir1d_();
  95. void ccor1d_(), ornl_ccor1d_();
  96. void simple_cfir1d_(), simple_ciir1d_(), simple_ccor1d_();
  97. #define SIMPLE_FIR        simple_cfir1d_
  98. #define SIMPLE_IIR        simple_ciir1d_
  99. #define SIMPLE_COR        simple_ccor1d_
  100. #define ORNL_FIR        ornl_cfir1d_
  101. #define ORNL_IIR        ornl_ciir1d_
  102. #define ORNL_COR        ornl_ccor1d_
  103. #define MY_COR            ccor1d_
  104. #define MY_FIR            cfir1d_
  105. #define MY_IIR            ciir1d_
  106. #define MY_INIT            cinit_
  107. #define MY_ONE            cone_
  108. #define THIS_COMPLEX
  109. this_type zero = {    0., 0.};
  110. this_type one  = {    1., 0.};
  111. #endif
  112.  
  113. /* ---------- The rest is the same    ---------------- */
  114.  
  115. #define MAX_SIZE        111
  116. #define MAX_STRIDE        7
  117. #define INC_STRIDE        2
  118. #define MAX_TIMES        3
  119. #define MIN_OPS            5.e+6
  120.  
  121. #define ABS(a) ( ((a)>0) ? (a) : -(a))
  122. #define MAX(a,b) (((a) < (b)) ? (b) : a)
  123.  
  124. void (*ornl_fir)(), (*my_fir)();
  125. void (*ornl_iir)(), (*my_iir)();
  126. void (*ornl_cor)(), (*my_cor)();
  127. void (*simple_fir)();
  128. void (*simple_iir)();
  129. void (*simple_cor)();
  130.  
  131. void GetArguments();
  132. double second();
  133.  
  134. int parallel;
  135. int is_speedup;
  136. int is_optimal;
  137. int n_trials;
  138. int all_run;
  139. int len = 4;
  140.  
  141.     int ldx, size, ldy, n_trials, n_times, nx, ny;
  142.     int min_size, max_size, inc_size, xsize, ysize, zsize;
  143.     this_type *vx, *vy, *vz;
  144.  
  145.     double t, x, y, z;
  146.     double total_flops;
  147.  
  148. main(argc,argv)
  149. int argc;
  150. char *argv[];
  151. {
  152.     int i, j, k;
  153.  
  154. /* ******************************************************* */
  155.     GetArguments( argc, argv);
  156. /* ******************************************************* */
  157.  
  158.     if( is_speedup) 
  159.     fprintf( stderr, "Measuring Relative Performances Lib_version/Simple_version \n");
  160.     else if( is_optimal)
  161.     fprintf( stderr, "Measuring  Performances for Lib_version (Mflops)\n");
  162.     else
  163.     fprintf( stderr, "Measuring  Performances for Simple_version (Mflops)\n");
  164.  
  165.     fprintf( stderr, 
  166.         "Size     FIR(1)    IIR(1)    COR(1)    FIR(%d)    IIR(%d)    COR(%d)\n\n",
  167.         MAX_STRIDE, MAX_STRIDE, MAX_STRIDE);
  168.     
  169.     xsize = MAX(max_size, nx) * MAX_STRIDE;
  170.     ysize = MAX(max_size, ny) * MAX_STRIDE;
  171.     zsize = xsize + ysize;
  172.  
  173.     vx = (this_type *)malloc( xsize * sizeof( this_type));
  174.     vy = (this_type *)malloc( ysize * sizeof( this_type));
  175.     vz = (this_type *)malloc( zsize * sizeof( this_type));
  176.  
  177.     if( (vx == (this_type *)0) || (vy == (this_type *)0))  {
  178.     fprintf( stderr, "Malloc problem ... Exiting");
  179.     exit( -2);
  180.     }
  181.     MY_INIT( (&xsize), vx);
  182.     MY_ONE( (&ysize), vy);
  183.  
  184.   for( size = min_size; size <= max_size ; size += inc_size ) { 
  185.     ysize = (ny == 0) ? size : ny;
  186.     xsize = (nx == 0) ? size : nx;
  187.     zsize = xsize+ysize;
  188.     printf("%4d ", size);
  189.     do_it();
  190.  
  191.     printf("\n", x);
  192.     fflush(stdout);
  193.   }
  194.  
  195.     free ( vx);
  196.     free ( vy);
  197.  
  198.     return(0);
  199. }
  200.  
  201.  
  202. do_it() 
  203. {
  204.     int ii, jj, inc, i1, i2, j1, j2, k1, k2;
  205.     int l1, l2;
  206.     i1 = 0; i2 = xsize;
  207.     j1 = 0; j2 = ysize;
  208.     k1=i1+j1; k2=i2+j2-1;
  209.     l1 = -(j1+j2-1);
  210.     l2 = j2;
  211.     
  212. /*
  213. *
  214. * Do it first with UNIT stride
  215. *
  216. */
  217.     inc = 1;
  218.  
  219. /*
  220. *   Finite Impulse Response Filter
  221. */
  222.     total_flops = xsize * ysize * OPS_PER_ITER;
  223.     n_times = MIN_OPS / total_flops;
  224.     if( n_times < 1) 
  225.     n_times = 1;
  226.     total_flops *= n_times;
  227.  
  228.     if( is_speedup || is_optimal ) { 
  229.     t = second();
  230.         for( ii = 0 ; ii < n_times ; ii++)
  231.         my_fir(vx,&inc,&i1,&i2,vy,&inc,&j1,&j2,    vz,&inc,&k1,&k2,
  232.             &one, &zero);
  233.     t = second() - t;
  234.     x = total_flops * 1.e-6 / t;
  235.     } 
  236.     if( is_speedup || !is_optimal) { 
  237.     x = second();
  238.         for( ii = 0 ; ii < n_times ; ii++)
  239.         simple_fir(&i2,&j2,vx,vy,vz);
  240.     x = second() - x;
  241.     if( is_speedup )
  242.         x = x / t;
  243.     else 
  244.         x = total_flops * 1.e-6 / x;
  245.     }
  246.     printf (" %8.3f ", x);
  247. /*
  248. *   Infinite Impulse Response Filter
  249. */
  250.     if( xsize < ysize)
  251.         total_flops = .5 * xsize * xsize * OPS_PER_ITER;
  252.     else
  253.         total_flops = ysize * ( xsize - .5 * ysize) * OPS_PER_ITER;
  254.     n_times = MIN_OPS / total_flops;
  255.     if( n_times < 1) 
  256.     n_times = 1;
  257.     total_flops *= n_times;
  258.  
  259.     if( is_speedup || is_optimal ) { 
  260.     t = second();
  261.         for( ii = 0 ; ii < n_times ; ii++)
  262.         my_iir(vx,&inc,&i1,&i2,vy,&inc,&j1,&j2,vz,&inc,&k1,&i2);
  263.     t = second() - t;
  264.     x = total_flops * 1.e-6 / t;
  265.     } 
  266.     if( is_speedup || !is_optimal) { 
  267.     x = second();
  268.         for( ii = 0 ; ii < n_times ; ii++)
  269.         simple_iir(&i2,&j2,vx,vy,vz);
  270.     x = second() - x;
  271.     if( is_speedup )
  272.         x = x / t;
  273.     else 
  274.         x = total_flops * 1.e-6 / x;
  275.     }
  276.     printf (" %8.3f ", x);
  277. /*
  278. *   CORRELATION
  279. */
  280.     total_flops = xsize * ysize * OPS_PER_ITER;
  281.     n_times = MIN_OPS / total_flops;
  282.     if( n_times < 1) 
  283.     n_times = 1;
  284.     total_flops *= n_times;
  285.  
  286.     if( is_speedup || is_optimal ) { 
  287.     t = second();
  288.         for( ii = 0 ; ii < n_times ; ii++)
  289.         my_cor(vx,&inc,&i1,&i2,vy,&inc,&l1,&l2,vz,&inc,&k1,&k2);
  290.     t = second() - t;
  291.     x = total_flops * 1.e-6 / t;
  292.     } 
  293.     if( is_speedup || !is_optimal) { 
  294.     x = second();
  295.         for( ii = 0 ; ii < n_times ; ii++)
  296.         simple_cor(&i2,&j2,vx,vy,vz);
  297.     x = second() - x;
  298.     if( is_speedup )
  299.         x = x / t;
  300.     else 
  301.         x = total_flops * 1.e-6 / x;
  302.     }
  303.     printf (" %8.3f ", x);
  304. /*
  305. *
  306. * Now do it again with NON_UNIT stride
  307. *
  308. */
  309.     inc = MAX_STRIDE;
  310.  
  311. /*
  312. *   Finite Impulse Response Filter
  313. */
  314.     total_flops = xsize * ysize * OPS_PER_ITER;
  315.     n_times = MIN_OPS / total_flops;
  316.     if( n_times < 1) 
  317.     n_times = 1;
  318.     total_flops *= n_times;
  319.  
  320.     if( is_speedup || is_optimal ) { 
  321.     t = second();
  322.         for( ii = 0 ; ii < n_times ; ii++)
  323.         my_fir(vx,&inc,&i1,&i2,vy,&inc,&j1,&j2,vz,&inc,&k1,&k2,
  324.             &one, &zero);
  325.     t = second() - t;
  326.     x = total_flops * 1.e-6 / t;
  327.     } 
  328.     if( is_speedup || !is_optimal) { 
  329.     x = second();
  330.         for( ii = 0 ; ii < n_times ; ii++)
  331.         ornl_fir(vx,&inc,&i1,&i2,vy,&inc,&j1,&j2,vz,&inc,&k1,&k2,
  332.             &one, &zero);
  333.     x = second() - x;
  334.     if( is_speedup )
  335.         x = x / t;
  336.     else 
  337.         x = total_flops * 1.e-6 / x;
  338.     }
  339.     printf (" %8.3f ", x);
  340. /*
  341. *   Infinite Impulse Response Filter
  342. */
  343.     if( xsize < ysize)
  344.         total_flops = .5 * xsize * xsize * OPS_PER_ITER;
  345.     else
  346.         total_flops = ysize * ( xsize - .5 * ysize) * OPS_PER_ITER;
  347.     n_times = MIN_OPS / total_flops;
  348.     if( n_times < 1) 
  349.     n_times = 1;
  350.     total_flops *= n_times;
  351.  
  352.     if( is_speedup || is_optimal ) { 
  353.     t = second();
  354.         for( ii = 0 ; ii < n_times ; ii++)
  355.         my_iir(vx,&inc,&i1,&i2,vy,&inc,&j1,&j2,vz,&inc,&k1,&i2);
  356.     t = second() - t;
  357.     x = total_flops * 1.e-6 / t;
  358.     } 
  359.     if( is_speedup || !is_optimal) { 
  360.     x = second();
  361.         for( ii = 0 ; ii < n_times ; ii++)
  362.         ornl_iir(vx,&inc,&i1,&i2,vy,&inc,&j1,&j2,vz,&inc,&k1,&i2);
  363.     x = second() - x;
  364.     if( is_speedup )
  365.         x = x / t;
  366.     else 
  367.         x = total_flops * 1.e-6 / x;
  368.     }
  369.     printf (" %8.3f ", x);
  370. /*
  371. *   CORRELATION
  372. */
  373.     total_flops = xsize * ysize * OPS_PER_ITER;
  374.     n_times = MIN_OPS / total_flops;
  375.     if( n_times < 1) 
  376.     n_times = 1;
  377.     total_flops *= n_times;
  378.  
  379.     if( is_speedup || is_optimal ) { 
  380.     t = second();
  381.         for( ii = 0 ; ii < n_times ; ii++)
  382.         my_cor(vx,&inc,&i1,&i2,vy,&inc,&l1,&l2,vz,&inc,&k1,&k2);
  383.     t = second() - t;
  384.     x = total_flops * 1.e-6 / t;
  385.     } 
  386.     if( is_speedup || !is_optimal) { 
  387.     x = second();
  388.         for( ii = 0 ; ii < n_times ; ii++)
  389.         ornl_cor(vx,&inc,&i1,&i2,vy,&inc,&l1,&l2,vz,&inc,&k1,&k2);
  390.     x = second() - x;
  391.     if( is_speedup )
  392.         x = x / t;
  393.     else 
  394.         x = total_flops * 1.e-6 / x;
  395.     }
  396.     printf (" %8.3f ", x);
  397. }
  398.  
  399. void GetArguments( argc, argv)
  400. int argc;
  401. char *argv[];
  402. {
  403.     int i, j, k;
  404.     int nerror = 0;
  405.  
  406.     srandom( (123*getpid()) | 0x01);
  407.  
  408. #define ON    1
  409.  
  410.     is_speedup = 0;
  411.     is_optimal = 0;
  412.     ny = 0;
  413.     nx = 0;
  414.     all_run = 0;
  415.  
  416.     min_size = 16;
  417.     max_size = 1024;
  418.     inc_size = 8;
  419.     ornl_fir = ORNL_FIR;
  420.     ornl_iir = ORNL_IIR;
  421.     ornl_cor = ORNL_COR;
  422.     simple_fir = SIMPLE_FIR;
  423.     simple_iir = SIMPLE_IIR;
  424.     simple_cor = SIMPLE_COR;
  425.     my_fir = MY_FIR;
  426.     my_iir = MY_IIR;
  427.     my_cor = MY_COR;
  428. /* ******************************************************* */
  429.     for ( i = 1 ; (i < argc) && (nerror != ON) ; i ++ ) {
  430.     if( argv[i][0] == '-') {
  431.         switch ( argv[i][1]) {
  432.  
  433.         case 'x' :
  434.         case 'X' :  
  435.                 nx = atoi( argv[++i]);
  436.                 break;
  437.         case 'y' :
  438.         case 'Y' :  
  439.                 ny = atoi( argv[++i]);
  440.                 break;
  441.         case 'z' :
  442.         case 'Z' :  
  443.                 is_speedup = 1;
  444.                 is_optimal = 0;
  445.                 break;
  446.         case 'p' :
  447.         case 'P' :  
  448.         case 's' :
  449.         case 'S' :  
  450.                 is_optimal = 1;
  451.                 is_speedup = 0;
  452.                 break;
  453.         default  :  nerror = ON;
  454.         }
  455.     }
  456.     else {
  457.         if( i+1 > argc)
  458.         nerror = ON;
  459.         else { 
  460.         min_size = atoi( argv[i]);i++;
  461.         max_size = atoi( argv[i]);i++;
  462.         inc_size = atoi( argv[i]);
  463.         }
  464.     }
  465.     }
  466.     if( nerror == ON) {
  467.     fprintf( stderr, 
  468.         "Usage : %s [-s / -z] [-x nx] [-y ny] <min max inc>\n", argv[0]);
  469.     exit(-1);
  470.     }
  471. }
  472. /* **********************************************************************
  473.  
  474.    give the elapsed wall clock time
  475.  
  476. ********************************************************************** */
  477. double second()
  478. {
  479.         struct timeval s_val;
  480.     struct timezone s_z;
  481.  
  482.     static double zero_time = 0.0;
  483.     static long zero_sec = 0;
  484.     double time;
  485.     long n_sec, n_usec;
  486.  
  487.         gettimeofday(&s_val, &s_z);
  488.  
  489.     n_sec = s_val.tv_sec;
  490.     n_usec = s_val.tv_usec;
  491.     if( zero_time == 0.0 ) {
  492.         zero_sec = n_sec;
  493.         zero_time = 1.0e-6 * (double)n_usec;
  494.     }
  495.     time = (double)(n_sec-zero_sec) + (double)n_usec * 1.0E-6 - zero_time;
  496.     return( time );
  497. }
  498.  
  499.