API tools faq
paste
Advertisement
bhavyagoyal

Best

bhavyagoyal
Nov 16th, 2015
389
0
Never
Add comment
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
C++ 12.14 KB | None | 0 0
raw download clone embed print report
  1. #include <stdio.h>
  2. #include <string.h>
  3. #include <stdlib.h>
  4. #include <sys/time.h>
  5. #include <time.h>
  6. #include <math.h>
  7. #include "evaluate.h"
  8. #include <omp.h>
  9.  
  10. void initialize(double** A, double** B, double** C, int* n);
  11. void* tune_matrix_multiply();
  12. void matrix_multiply(double* A, double* B, double* C, int n, void* tunning_params);
  13. void getUserID(char userID[20],char userName[30]);
  14.  
  15.  
  16. int main(int argc, char *argv[]){
  17.     srand(time(NULL));
  18.     void* tunning_params = tune_matrix_multiply();
  19.  
  20.     int n;
  21.     double *A, *B, *C, tolerance = pow(10,-7);
  22.     char userID[20];
  23.     char userName[30];
  24.     getUserID(userID,userName);
  25.  
  26.     if (debug_run(argc, argv,userID,userName)) {  /* This is debug run for students to test their code. Do not edit this. */
  27.         initialize(&A, &B, &C, &n);
  28.         double etOpenMPst, etOpenMPsp;
  29.         etOpenMPst = omp_get_wtime();
  30.         matrix_multiply(A, B, C, n, tunning_params); /* Multiply B and C to store the results into A <- B * C */
  31.         etOpenMPsp = omp_get_wtime();
  32.         double total_flops = myCalcFlops(etOpenMPst, etOpenMPsp, n);
  33.         printf("Flops : %lf\n",total_flops);
  34.  
  35.         int check = check_Correctness(A, B, C, n, tolerance);
  36.         if(check) printf("Correct : true \n");
  37.         else printf("Correct : false \n");
  38.         free(A);
  39.         free(B);
  40.         free(C);
  41.        
  42.     } else{
  43.         /* This is an evaluation run
  44.          * and it will be used for evaluation of student answers
  45.          * so dont wory about the code in the else block.
  46.         */
  47.         int sizes[4];
  48.  
  49.         find_matrix_sizes(argc, argv, sizes);
  50.        
  51.         double limit = 100.0, flops[4]= {0.0,0.0,0.0,0.0};
  52.         int s, correctness[4] = {0,0,0,0};
  53.         for(s=0;s<4;s++){
  54.             int n0=sizes[s];
  55.             //allocate memory of A, B, C
  56.             A = (double *) malloc(n0*n0*sizeof(double));
  57.             B = (double *) malloc(n0*n0*sizeof(double));
  58.             C = (double *) malloc(n0*n0*sizeof(double));
  59.             int j, i, k;
  60.  
  61.             for(j=0;j<5;j++){
  62.                 // Randomly initializing elements B & C from 0-limit
  63.                 /* B is a Lower Triangular Matrix and C is a n/16 banded Matrix */
  64.                 for(i=0;i<n0;i++){
  65.                     for(k=0;k<n0;k++){
  66.                         *(A + n0*i + k) = 0.0;
  67.                        
  68.                         if (k>i) *(B + n0*i + k) = 0.0;
  69.                         else *(B + n0*i + k) = (double)rand()/(double)(RAND_MAX/limit);
  70.                        
  71.                         if (k>i+n0/16-1 || i>k+n0/16-1) *(C + n0*i + k) = 0.0;
  72.                         else *(C + n0*i + k) = (double)rand()/(double)(RAND_MAX/limit);
  73.                     }
  74.                 }
  75.                
  76.                 if(j==0){
  77.                     double etOpenMPst, etOpenMPsp;
  78.                     etOpenMPst = omp_get_wtime();
  79.                     matrix_multiply(A, B, C, n0, tunning_params);
  80.                     etOpenMPsp = omp_get_wtime();
  81.                 }
  82.                 else{
  83.                     double etOpenMPst, etOpenMPsp;
  84.                     etOpenMPst = omp_get_wtime();
  85.                     matrix_multiply(A, B, C, n0, tunning_params);
  86.                     etOpenMPsp = omp_get_wtime();
  87.  
  88.                     flops[s] += myCalcFlops(etOpenMPst,etOpenMPsp, n0);
  89.                    
  90.                     int check = check_Correctness(A, B, C, n0, tolerance);
  91.                     if(check==1){
  92.                         correctness[s]++;
  93.                     }
  94.                 }
  95.             }
  96.             flops[s] = flops[s]/4.0;
  97.             free(A);
  98.             free(B);
  99.             free(C);
  100.         }
  101.         double overall_flop_score = (flops[0] + 2 * flops[1] + 3 * flops[2] + 4 * flops[3])/10;
  102.         int q,overall_correctness = 1;
  103.         for(q=0;q<4;q++){
  104.             if(correctness[q] != 4){
  105.                 overall_correctness = 0;
  106.             }
  107.         }
  108.             printf("UserName:%s\tuserID:%s\n",userName,userID);
  109.         printf("Size\tcorrectness\tflops\n");
  110.         printf("small \t\t%d\t%lf\n",correctness[0],flops[0]);
  111.         printf("medium \t\t%d\t%lf\n",correctness[1],flops[1]);
  112.         printf("large \t\t%d\t%lf\n",correctness[2],flops[2]);
  113.         printf("extralarge \t%d\t%lf\n",correctness[3],flops[3]);
  114.         printf("\n");
  115.         printf("Score = %lf\n",overall_correctness * overall_flop_score);
  116.        
  117.         FILE* outputfile=fopen("output.csv","w");
  118.         fprintf(outputfile, "\"%s\",\"%s\",%d,%lf,%d,%lf,%d,%lf,%d,%lf,%lf",userID,userName,correctness[0],flops[0],correctness[1],flops[1],correctness[2],flops[2],correctness[3],flops[3],overall_correctness * overall_flop_score);
  119.         fclose(outputfile);
  120.         FILE* outputfile2=fopen("output2.csv","w");
  121.         fprintf(outputfile, "\"%s\",%d,%lf,%d,%lf,%d,%lf,%d,%lf,%lf",userName,correctness[0],flops[0],correctness[1],flops[1],correctness[2],flops[2],correctness[3],flops[3],overall_correctness * overall_flop_score);
  122.         fclose(outputfile2);
  123.     }
  124.     return 0;
  125. }
  126. //Add only code in the functions
  127. void getUserID(char userID[20],char userName[30]){
  128.     //Enter your moodle userID here i.e. your entry number
  129.     strcpy(userID,"jca142468");
  130.     //Enter your name here
  131.     strcpy(userName,"Rohit Katariya");
  132. }
  133.  
  134. void initialize(double** A, double** B, double** C, int* n){
  135.     /*
  136.     Initialize matrix A , B & C where each of them is a square matrix of size n*n
  137.     Try 4 sizes of matrix:
  138.     a) small - should fit in the L1 cache
  139.     b) medium - should fit in the L2 cache
  140.     b) large - should fit in the L3 cache
  141.     b) extra-large - should fit in the memory
  142.     */
  143.     *n=1024;
  144.     *A = (double *) malloc(*n**n*sizeof(double));
  145.     *B = (double *) malloc(*n**n*sizeof(double));
  146.     *C = (double *) malloc(*n**n*sizeof(double));
  147.     int i, k;
  148.     for(i=0;i<*n;i++){
  149.         for(k=0;k<*n;k++){
  150.             *(*A + *n*i + k) = 0.0;
  151.                        
  152.             if(k>i) *(*B + *n*i + k) = 0.0;
  153.             else *(*B + *n*i + k) = (double)rand()/(double)(RAND_MAX/100.0);
  154.                        
  155.             if(k>i+*n/16-1 || i>k+*n/16-1) *(*C + *n*i + k) = 0.0;
  156.             else *(*C + *n*i + k) = (double)rand()/(double)(RAND_MAX/100.0);
  157.         }
  158.                    
  159.     }
  160. }
  161. #include <assert.h>
  162.  
  163. void* tune_matrix_multiply(){
  164.     /* Write code to test various methods to optimize the matrix multiplication
  165.        and store these parameters in a data structure(for ex, array of values ). This
  166.        data structure can be returned by this function and will be passed to the
  167.        matrix multiply function so that it can efficiently multiply the matrices.
  168.     */
  169.    
  170.    
  171. }
  172.  
  173. #define MAXT(A,B) ({ __typeof__ (A) a = (A); __typeof__ (B) b = (B); a > b ? a : b; })
  174. #define MINT(A,B) ({ __typeof__ (A) a = (A); __typeof__ (B) b = (B); a > b ? b : a; })
  175.  
  176. #define maX(A,B) ({ int a = (A); int b = (B); a > b ? a : b; })
  177. #define miN(A,B) ({ int a = (A); int b = (B); a > b ? b : a; })
  178. #define MAX(A,B) ({ ((A) > (B)) ? A : B; })
  179. #define MIN(A,B) ({ ((A) > (B)) ? B : A; })
  180.  
  181. void matrix_multiply(double* A, double* B, double* C, int n, void* tunning_params){
  182.     // Multiply matrix B & C and store the results in A
  183.  
  184.     int i,j,k,procsize,threadsize;
  185.     // double *pA, *pB, *pC;
  186.     // int chunksize;
  187.     // chunksize = (n==256)?2:1;
  188.     // omp_set_num_threads(2);
  189.     // Change omp_proc_bind
  190.     procsize = omp_get_num_procs();
  191.     threadsize = omp_get_max_threads();
  192.  
  193.     if(n&32) {
  194.         // #pragma omp for collapse(2) schedule(static, 2)
  195.         // for(i=0;i<n;i++){
  196.         //     for(j=0;j<n;j++){
  197.         //         A[n*i+j] = 0;
  198.         //         int lim1 = MAX(0,j-n/16+1);
  199.         //         int lim2 = MIN(n, j+n/16);
  200.         //         register double temp=0;
  201.         //         for(k=lim1;k<lim2;k++){
  202.         //             temp += B[n*i+k]*C[n*k+j];
  203.         //         }
  204.         //         A[n*i+j] = temp;
  205.         //     }
  206.         // }
  207.         for(i = 0 ; i < n ; i++){
  208.             int l, p ;
  209.             for(j = 0 , l = MIN(i+2,n) , p = i<<5 ; j < i + 1 ; j++){
  210.                 A[p+j] += B[p + j] * C[(j<<5) + j];
  211.                 if(i && j < i)  A[p + j] += B[p + j + 1] * C[((j+1)<<5) + j];
  212.                 if(j + 1 && j + 1 < l)  A[p + j + 1] += B[p + j]*C[(j<<5) + j + 1];
  213.             }
  214.         }
  215.         return;
  216.         // for(i=0;i<n;i++){
  217.         //     double *pA = A + (n*i);
  218.         //     for(k=0;k< (i+1);k++){
  219.         //         double *pC = C + n*k;
  220.         //         int lim2 = MIN(n, k+(n>>4));
  221.         //         int lim1 = MAX(0,k-(n>>4)+1);
  222.         //         // register double temp = B[n*i+k];
  223.         //         for(j=lim1;j<lim2;j++){
  224.         //             *(pA+j) += B[n*i+k]**(pC+j);
  225.         //         }
  226.         //     }
  227.         // }
  228.         // return;
  229.     }
  230.  
  231.     if(n&64) {
  232.         int l , q , p;
  233.         for(i = 0 ; i < n ; i++){
  234.             for(j = 0 , l = MIN(i+2,n) ,q = MIN(i+3,n) , p = i<<6 ; j < i + 1 ; j++){
  235.                 A[p+j] += B[p + j] * C[(j<<6) + j];
  236.                 if(i && j < i)  A[p + j] += B[p + j + 1] * C[((j+1)<<6) + j];
  237.                 if(i-1 && j < i-1)  A[p + j] += B[p + j + 2] * C[((j+2)<<6) + j];
  238.                 if(i-2 && j < i-2)  A[p + j] += B[p + j + 3] * C[((j+3)<<6) + j];
  239.                 if(j + 1 && j + 1 < l)  A[p + j + 1] += B[p + j]*C[(j<<6) + j + 1];
  240.                 if(j + 2 && j + 2 < q)  A[p + j + 2] += B[p + j]*C[(j<<6) + j + 2];
  241.                 if(j + 3 && j + 3 < MIN(i+4,n))  A[p + j + 3] += B[p + j]*C[(j<<6) + j + 3];
  242.             }
  243.         }
  244.         return;
  245.         // #pragma omp parallel for num_threads(2) default(none) private(k,j) shared(n,A,B,C)
  246.         // for(i=0;i<n;i++){
  247.         //     double *pA = A + (n*i);
  248.         //     for(k=0;k< (i+1);k++){
  249.         //         double *pC = C + n*k;
  250.         //         int lim2 = MIN(n, k+(n>>4));
  251.         //         int lim1 = MAX(0,k-(n>>4)+1);
  252.         //         // register double temp = B[n*i+k];
  253.         //         for(j=lim1;j<lim2;j++){
  254.         //             *(pA+j) += B[n*i+k]**(pC+j);
  255.         //         }
  256.         //     }
  257.         // }
  258.         // return;
  259.     }
  260.  
  261.         // #pragma omp for schedule(guided)
  262.         // for(i=0;i<n;i++){
  263.         //     for(j=i+1;j<n;j++){
  264.         //         double temp=B[n*i+j];
  265.         //         B[n*i+j] = B[n*j+i];
  266.         //         B[n*j+i] = temp;
  267.         //     }
  268.         // }
  269.  
  270.  
  271.     if(n&256){
  272.         // #pragma omp parallel default(none) private(i,k,j) shared(n,A,B,C)
  273.         // #pragma omp for schedule(dynamic, 4)
  274.         #pragma omp parallel for default(none) schedule(dynamic,4) private(k,j) shared(n,A,B,C) num_threads(7)
  275.         for(i=n-1 ; i>=0 ; i--){
  276.             double *pA = A+(n*i);
  277.             for(k=0;k<(i+1);k++){
  278.                 double *pC = C+n*k;
  279.                 // double *pB = B+n*i+k
  280.                 int lim2 = MIN(n, k+(n>>4));
  281.                 int lim1 = MAX(0,k-(n>>4)+1);
  282.                 register double temp = B[n*i+k];
  283.                 for(j=lim1;j<lim2;j++){
  284.                     *(pA+j) += temp**(pC+j);
  285.                 }
  286.             }
  287.         }
  288.         return;
  289.     }
  290.    
  291.  
  292.             // double *pA = A+(n*i);
  293.             // double *pB = B+(n*i);
  294.  
  295.         #pragma omp parallel for default(none) schedule(dynamic) private(i,k,j) shared(n,A,B,C) num_threads(procsize+1)
  296.         for(i=0; i<n; i++){
  297.             double *pA = A+(n*i);
  298.             for(k=0;k<(i+1);k++){
  299.                 double *pC = C+n*k;
  300.                 int lim2 = MIN(n, k+(n>>4));
  301.                 int lim1 = MAX(0, k-(n>>4)+1);
  302.                 register double temp=B[n*i+k];
  303.                 for(j=lim1;j<lim2;j++){
  304.                     *(pA+j) += temp**(pC+j);
  305.                 }
  306.             }
  307.             // int ii = i+n/2;
  308.             // pA = A + (n*(ii));
  309.             // for(k=0;k<(ii+1);k++){
  310.             //     double *pC = C+n*k;
  311.             //     int lim2 = MIN(n, k+(n>>4));
  312.             //     int lim1 = MAX(0, k-(n>>4)+1);
  313.             //     register double temp=B[n*ii+k];
  314.             //     for(j=lim1;j<lim2;j++){
  315.             //         *(pA+j) += temp**(pC+j);
  316.             //     }
  317.             // }
  318.         }
  319.         return;
  320.  
  321.  
  322.  
  323.  
  324.     // Basic Parallel
  325.     // int i=0,j=0,k=0;
  326.     // // int temp=0;
  327.     // // printf("Number of procs %d\n",omp_get_num_procs());
  328.     // // omp_set_num_threads(omp_get_num_procs());
  329.     // #pragma omp parallel for private(j)
  330.     //     for(i=0;i<n;i++){
  331.     //         for(k=0;k<n;k++){
  332.     //             for(j=0;j<n;j++){
  333.     //                 A[n*i+j] += B[n*i+k]*C[n*k+j];
  334.     //             }
  335.     //         }
  336.     //     }
  337.  
  338.  
  339.  
  340. /*
  341.  
  342. Version 2 - 1.65 Gflops
  343.     int i=0,j=0,k=0;
  344.     int chunk=10;
  345.     // int temp=0;
  346.     // printf("Number of procs %d\n",omp_get_num_procs());
  347.     // omp_set_num_threads(omp_get_num_procs());
  348.     #pragma omp parallel default(none) shared(A,B,C,n,chunk) private(i,j,k)
  349.     {
  350.         #pragma omp for collapse(2)
  351.         for(i=0;i<n;i++){
  352.             for(j=0;j<n;j++){
  353.                 // A[n*i+j] = 0;
  354.                 double temp=0;
  355.                 for(k=0;k<n;k++){
  356.                     temp += B[n*i+k]*C[n*k+j];
  357.                 }
  358.                 A[n*i+j]=temp;
  359.             }
  360.         }
  361.     }
  362. */
  363. }
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement
Public Pastes
Advertisement
We use cookies for various purposes including analytics. By continuing to use Pastebin, you agree to our use of cookies as described in the Cookies Policy.  OK, I Understand
Not a member of Pastebin yet?
Sign Up, it unlocks many cool features!