Facebook
From Morose Pig, 6 Years ago, written in Plain Text.
Embed
Download Paste or View Raw
Hits: 264
  1. // classic.cpp : "Textbook" implementation of matrix multiply
  2.  
  3. // Author:  Paul J. Drongowski
  4. // Address: Boston Design Center
  5. //          Advanced Micro Devices, Inc.
  6. //          Boxborough, MA 01719
  7. // Date:    20 October 2005
  8. //
  9. // Copyright (c) 2005 Advanced Micro Devices, Inc.
  10.  
  11. // Celem tego programu jest prezentacja pomiaru i analizy
  12. //efektywnosci programu za pomocą  CodeAnalyst(tm).
  13. // Implementacja mnożenia macierzy jest realizowana za pomoca typowego
  14. // algorytmu podręcznikowego.
  15. //  Efektywność tego podejscia jest niska poprzez
  16. // nieefektywną  kolejnosć odwołań do elementów macierzy.
  17. #include <stdio.h>
  18. #include <time.h>
  19. #include <windows.h>
  20. #include "omp.h"
  21. #include <string>
  22.  
  23. #define USE_MULTIPLE_THREADS true
  24. #define MAXTHREADS 128
  25. int NumThreads;
  26. double start;
  27.  
  28. static const int ROWS = 512;     // liczba wierszy macierzy
  29. static const int COLUMNS = 512;  // lizba kolumn macierzy
  30. static const int R = 128;
  31.  
  32. float matrix_a[ROWS][COLUMNS];    // lewy operand
  33. float matrix_b[ROWS][COLUMNS];    // prawy operand
  34. float matrix_r[ROWS][COLUMNS];    // wynik
  35.  
  36. float result = 0;
  37. float resultParallel = 0;
  38. float result6 = 0;
  39. float resultParallel6 = 0;
  40.  
  41. float minResult = 100;
  42. float minResultParallel = 100;
  43. float minResult6 = 100;
  44. float minResultParallel6 = 100;
  45.  
  46. float maxResult = 0;
  47. float maxResultParallel = 0;
  48. float maxResult6 = 0;
  49. float maxResultParallel6 = 0;
  50.  
  51. float tmp;
  52. int r = 256;
  53.  
  54. FILE *result_file;
  55.  
  56. void initialize_matrices()
  57. {
  58.         // zdefiniowanie zawarosci poczatkowej macierzy
  59. #pragma omp parallel for
  60.         for (int i = 0; i < ROWS; i++) {
  61.                 for (int j = 0; j < COLUMNS; j++) {
  62.                         matrix_a[i][j] = (float)rand() / RAND_MAX;
  63.                         matrix_b[i][j] = (float)rand() / RAND_MAX;
  64.                         matrix_r[i][j] = 0.0;
  65.                 }
  66.         }
  67. }
  68.  
  69. void initialize_matricesZ()
  70. {
  71.         // zdefiniowanie zawarosci poczatkowej macierzy
  72.  
  73. #pragma omp parallel for
  74.         for (int i = 0; i < ROWS; i++) {
  75.                 for (int j = 0; j < COLUMNS; j++) {
  76.                         matrix_r[i][j] = 0.0;
  77.                 }
  78.         }
  79. }
  80.  
  81. void print_result()
  82. {
  83.         // wydruk wyniku
  84.         for (int i = 0; i < ROWS; i++) {
  85.                 for (int j = 0; j < COLUMNS; j++) {
  86.                         fprintf(result_file, "%6.4f ", matrix_r[i][j]);
  87.                 }
  88.                 fprintf(result_file, "\n");
  89.         }
  90. }
  91.  
  92. void multiply_matrices_IKJ()
  93. {
  94.         for (int i = 0; i < ROWS; i++)
  95.                 for (int k = 0; k < COLUMNS; k++)
  96.                         for (int j = 0; j < COLUMNS; j++)
  97.                                 matrix_r[i][j] += matrix_a[i][k] * matrix_b[k][j];
  98. }
  99.  
  100. void multiply_matrices_IKJ_Parallel()
  101. {
  102.  
  103. #pragma omp parallel for
  104.         for (int i = 0; i < ROWS; i++)
  105.                 for (int k = 0; k < COLUMNS; k++)
  106.                         for (int j = 0; j < COLUMNS; j++)
  107.                         {
  108.                                 matrix_r[i][j] += matrix_a[i][k] * matrix_b[k][j];
  109.                         }
  110.  
  111. }
  112.  
  113. void multiply_matrices_JKI()
  114. {
  115.         // mnozenie macierzy
  116.         for (int j = 0; j < COLUMNS; j++)
  117.                 for (int k = 0; k < COLUMNS; k++)
  118.                         for (int i = 0; i < ROWS; i++)
  119.                                 matrix_r[i][j] += matrix_a[i][k] * matrix_b[k][j];
  120.  
  121. }
  122.  
  123. void multiply_matrices_JKI_Parallel()
  124. {
  125.         // mnozenie macierzy
  126. #pragma omp parallel for
  127.         for (int j = 0; j < COLUMNS; j++)
  128.                 for (int k = 0; k < COLUMNS; k++)
  129.                         for (int i = 0; i < ROWS; i++)
  130.                                 matrix_r[i][j] += matrix_a[i][k] * matrix_b[k][j];
  131.  
  132. }
  133.  
  134. void multiply_matrices_IJK_IKJ()
  135. {
  136.  
  137.         for (int i = 0; i < ROWS; i += r)
  138.                 for (int j = 0; j < COLUMNS; j += r)
  139.                         for (int k = 0; k < COLUMNS; k += r)
  140.                                 for (int ii = i; ii < i + r; ii++)  // WYNIK CZĘŚCIOWY
  141.                                         for (int kk = k; kk < k + r; kk++)
  142.                                                 for (int jj = j; jj < j + r; jj++)
  143.                                                         matrix_r[ii][jj] += matrix_a[ii][kk] * matrix_b[kk][jj];
  144. }
  145.  
  146. void multiply_matrices_IJK_IKJ_Parallel()
  147. {
  148.  
  149. #pragma omp parallel for
  150.         for (int i = 0; i < ROWS; i += r)
  151.                 for (int j = 0; j < COLUMNS; j += r)
  152.                         for (int k = 0; k < COLUMNS; k += r)
  153.                                 for (int ii = i; ii < i + r; ii++)  // WYNIK CZĘŚCIOWY
  154.                                         for (int kk = k; kk < k + r; kk++)
  155.                                                 for (int jj = j; jj < j + r; jj++)
  156.                                                         matrix_r[ii][jj] += matrix_a[ii][kk] * matrix_b[kk][jj];
  157. }
  158.  
  159. float print_elapsed_time(std::string info)
  160. {
  161.         double elapsed;
  162.         double resolution;
  163.  
  164.         // wyznaczenie i zapisanie czasu przetwarzania
  165.         elapsed = (double)clock() / CLK_TCK;
  166.         resolution = 1.0 / CLK_TCK;
  167.  
  168.         fprintf(result_file,
  169.                 "Czas wykonania programu %s: %8.4f sec (%6.4f sec rozdzielczosc pomiaru)\n",
  170.                 info.c_str(), elapsed - start, resolution);
  171.         return elapsed - start;
  172. }
  173.  
  174. void enterLine()
  175. {
  176.         fprintf(result_file,
  177.                 "\n"
  178.         );
  179.         fprintf(result_file,
  180.                 "\n"
  181.         );
  182. }
  183.  
  184. void load_IKJ()
  185. {
  186.         for (int i = 0;i < 100;i++)
  187.         {
  188.                 initialize_matricesZ();
  189.                 start = (double)clock() / CLK_TCK;
  190.                 multiply_matrices_IKJ();
  191.                 tmp = print_elapsed_time("Sekwencyjnie IKJ ");
  192.                 result += tmp;
  193.  
  194.                 if (tmp > maxResult)
  195.                 {
  196.                         maxResult = tmp;
  197.                 }
  198.  
  199.                 if (tmp < minResult)
  200.                 {
  201.                         minResult = tmp;
  202.                 }
  203.         }
  204.         enterLine();
  205. }
  206.  
  207. void loadParallel_IKJ()
  208. {
  209.         for (int i = 0; i < 100; i++)
  210.         {
  211.                 initialize_matricesZ();
  212.                 start = (double)clock() / CLK_TCK;
  213.                 multiply_matrices_IKJ_Parallel();
  214.                 tmp = print_elapsed_time("Rownolegle IKJ ");
  215.                 resultParallel += tmp;
  216.  
  217.                 if (tmp > maxResultParallel)
  218.                 {
  219.                         maxResultParallel = tmp;
  220.                 }
  221.  
  222.                 if (tmp < minResultParallel)
  223.                 {
  224.                         minResultParallel = tmp;
  225.                 }
  226.         }
  227.         enterLine();
  228. }
  229.  
  230. void load_IJK_IKJ()
  231. {
  232.         for (int i = 0; i < 100; i++)
  233.         {
  234.                 initialize_matricesZ();
  235.                 start = (double)clock() / CLK_TCK;
  236.                 multiply_matrices_IJK_IKJ();
  237.                 tmp = print_elapsed_time("Sekwencyjnie IJK_IKJ ");
  238.                 result6 += tmp;
  239.  
  240.                 if (tmp > maxResult6)
  241.                 {
  242.                         maxResult6 = tmp;
  243.                 }
  244.  
  245.                 if (tmp < minResult6)
  246.                 {
  247.                         minResult6 = tmp;
  248.                 }
  249.         }
  250.         enterLine();
  251. }
  252.  
  253. void loadParallel_IJK_IKJ()
  254. {
  255.         for (int i = 0; i < 100; i++)
  256.         {
  257.                 initialize_matricesZ();
  258.                 start = (double)clock() / CLK_TCK;
  259.                 multiply_matrices_IJK_IKJ_Parallel();
  260.                 tmp = print_elapsed_time("Rownolegle IJK_IKJ");
  261.                 resultParallel6 += tmp;
  262.  
  263.                 if (tmp > maxResultParallel6)
  264.                 {
  265.                         maxResultParallel6 = tmp;
  266.                 }
  267.  
  268.                 if (tmp < minResultParallel6)
  269.                 {
  270.                         minResultParallel6 = tmp;
  271.                 }
  272.         }
  273.         enterLine();
  274. }
  275.  
  276. void save_IKJ()
  277. {
  278.         float avgResult = result / 100;
  279.         fprintf(result_file,
  280.                 "Średnia sekwencyjnego %8.4f \n", avgResult
  281.         );
  282.         fprintf(result_file,
  283.                 "Max sekwencyjnego %8.4f \n", maxResult
  284.         );
  285.         fprintf(result_file,
  286.                 "Min sekwencyjnego %8.4f \n", minResult
  287.         );
  288.         enterLine();
  289. }
  290.  
  291. void saveParallel_IKJ()
  292. {
  293.         float avgResultParallel = resultParallel / 100;
  294.         fprintf(result_file,
  295.                 "Średnia równoległego %8.4f \n", avgResultParallel
  296.         );
  297.         fprintf(result_file,
  298.                 "Min równoległego %8.4f \n", minResultParallel
  299.         );
  300.         fprintf(result_file,
  301.                 "Max równoległego %8.4f \n", maxResultParallel
  302.         );
  303.         fclose(result_file);
  304. }
  305.  
  306. void saveParallel_IJK_IKJ()
  307. {
  308.         float avgResultParallel = resultParallel6 / 100;
  309.         fprintf(result_file,
  310.                 "Średnia równoległego IJK_IKJ   %8.4f \n", avgResultParallel
  311.         );
  312.         fprintf(result_file,
  313.                 "Min równoległego IJK_IKJ   %8.4f \n", minResultParallel6
  314.         );
  315.         fprintf(result_file,
  316.                 "Max równoległego IJK_IKJ   %8.4f \n", maxResultParallel6
  317.         );
  318.         fclose(result_file);
  319. }
  320.  
  321. void save_IJK_IKJ()
  322. {
  323.         float avgResult = result6 / 100;
  324.         fprintf(result_file,
  325.                 "Średnia sekwencyjnego IJK_IKJ %8.4f \n", avgResult
  326.         );
  327.         fprintf(result_file,
  328.                 "Max sekwencyjnego IJK_IKJ %8.4f \n", maxResult6
  329.         );
  330.         fprintf(result_file,
  331.                 "Min sekwencyjnego IJK_IKJ %8.4f \n", minResult6
  332.         );
  333.         enterLine();
  334. }
  335.  
  336. int main(int argc, char* argv[])
  337. {
  338.         omp_set_num_threads(4);
  339.  
  340.         //       start = (double) clock() / CLK_TCK ;
  341.         if ((result_file = fopen("classic.txt", "w")) == NULL) {
  342.                 fprintf(stderr, "nie mozna otworzyc pliku wyniku \n");
  343.                 perror("classic");
  344.                 return(EXIT_FAILURE);
  345.         }
  346.  
  347.  
  348.         //Determine the number of threads to use
  349.         if (USE_MULTIPLE_THREADS) {
  350.                 SYSTEM_INFO SysInfo;
  351.                 GetSystemInfo(&SysInfo);
  352.                 NumThreads = SysInfo.dwNumberOfProcessors;
  353.                 if (NumThreads > MAXTHREADS)
  354.                         NumThreads = MAXTHREADS;
  355.         }
  356.         else
  357.                 NumThreads = 1;
  358.         fprintf(result_file, "Klasyczny algorytm mnozenia macierzy, liczba watkow %d \n", NumThreads);
  359.         printf("liczba watkow  = %d\n\n", NumThreads);
  360.  
  361.         //load_IKJ();
  362.         //loadParallel_IKJ();
  363.  
  364.         //save_IKJ();
  365.         //saveParallel_IKJ();
  366.  
  367.         load_IJK_IKJ();
  368.         loadParallel_IJK_IKJ();
  369.        
  370.        
  371.         save_IJK_IKJ();
  372.         saveParallel_IJK_IKJ();
  373.  
  374.         return(0);
  375. }