mt/az_matmul/az_matmul.c

   1 //**************************************************************************
   2 // Multi-threaded Matrix Multiply benchmark
   3 //--------------------------------------------------------------------------
   4 // TA     : Christopher Celio
   5 // Student:
   6 //
   7 //
   8 // This benchmark multiplies two 2-D arrays together and writes the results to
   9 // a third vector. The input data (and reference data) should be generated
  10 // using the matmul_gendata.pl perl script and dumped to a file named
  11 // dataset.h.
  12
  13
  14 // print out arrays, etc.
  15 //#define DEBUG
  16
  17 //--------------------------------------------------------------------------
  18 // Includes
  19
  20 #include <string.h>
  21 #include <stdlib.h>
  22 #include <stdio.h>
  23
  24
  25 //--------------------------------------------------------------------------
  26 // Input/Reference Data
  27
  28 typedef float data_t;
  29 #include "dataset.h"
  30
  31
  32 //--------------------------------------------------------------------------
  33 // Basic Utilities and Multi-thread Support
  34
  35 __thread unsigned long coreid;
  36 unsigned long ncores;
  37
  38 #include "util.h"
  39
  40 #define stringify_1(s) #s
  41 #define stringify(s) stringify_1(s)
  42 #define stats(code) do { \
  43     unsigned long _c = -rdcycle(), _i = -rdinstret(); \
  44     code; \
  45     _c += rdcycle(), _i += rdinstret(); \
  46     if (coreid == 0) \
  47       printf("%s: %ld cycles, %ld.%ld cycles/iter, %ld.%ld CPI\n", \
  48              stringify(code), _c, _c/DIM_SIZE/DIM_SIZE/DIM_SIZE, 10*_c/DIM_SIZE/DIM_SIZE/DIM_SIZE%10, _c/_i, 10*_c/_i%10); \
  49   } while(0)
  50
  51
  52 //--------------------------------------------------------------------------
  53 // Helper functions
  54
  55 void printArray( char name[], int n, data_t arr[] )
  56 {
  57    int i;
  58    if (coreid != 0)
  59       return;
  60
  61    printf( " %10s :", name );
  62    for ( i = 0; i < n; i++ )
  63       printf( " %3ld ", (long) arr[i] );
  64    printf( "\n" );
  65 }
  66
  67 void __attribute__((noinline)) verify(size_t n, const data_t* test, const data_t* correct)
  68 {
  69    if (coreid != 0)
  70       return;
  71
  72    size_t i;
  73    for (i = 0; i < n; i++)
  74    {
  75       if (test[i] != correct[i])
  76       {
  77          printf("FAILED test[%d]= %3ld, correct[%d]= %3ld\n",
  78             i, (long)test[i], i, (long)correct[i]);
  79          exit(-1);
  80       }
  81    }
  82
  83    return;
  84 }
  85
  86 //--------------------------------------------------------------------------
  87 // matmul function
  88
  89 // single-thread, naive version
  90 void __attribute__((noinline)) matmul_naive(const int lda,  const data_t A[], const data_t B[], data_t C[] )
  91 {
  92    int i, j, k;
  93
  94    if (coreid > 0)
  95       return;
  96
  97    for ( i = 0; i < lda; i++ )
  98       for ( j = 0; j < lda; j++ )
  99       {
 100          for ( k = 0; k < lda; k++ )
 101          {
 102             C[i + j*lda] += A[j*lda + k] * B[k*lda + i];
 103          }
 104       }
 105
 106 }
 107
 108
 109 data_t ffmul(data_t a, data_t b) {
 110   data_t result = 0;
 111
 112   for (int i=0; i < b; i++) {
 113     result += a;
 114   }
 115
 116   return result;
 117 }
 118
 119
 120 //void __attribute__((noinline)) matmul(const int lda,  const data_t A[], const data_t B[], data_t C[] )
 121 //{
 122 //
 123 //   // ***************************** //
 124 //   // **** ADD YOUR CODE HERE ***** //
 125 //   // ***************************** //
 126 //   //
 127 //   // feel free to make a separate function for MI and MSI versions.
 128 //
 129 //    static __thread int i, j, k;
 130 //    static __thread int jlda, ilda;
 131 //    static __thread data_t tempA1, tempA2, tempA3, tempA4, tempA5, tempA6, tempA7, tempA8;
 132 //    static __thread int start, end;
 133 //
 134 //    start = coreid*(lda>>1);
 135 //    end = (coreid+1)*(lda>>1);
 136 //
 137 //    for (j=start; j < end; j+=1) {
 138 //      jlda = j * lda;
 139 //      for ( i=0; i < lda; i+=1 ) {
 140 //        ilda = i*lda;
 141 //        tempA1 = A[i   + jlda];
 142 //        //tempA2 = A[i+1 + jlda];
 143 //        //tempA3 = A[i+2 + jlda];
 144 //        //tempA4 = A[i+3 + jlda];
 145 //        //tempA5 = A[i+4 + jlda];
 146 //        //tempA6 = A[i+5 + jlda];
 147 //        //tempA7 = A[i+6 + jlda];
 148 //        //tempA8 = A[i+7 + jlda];
 149 //        //tempC1 = C[i + j*lda];
 150 //        //tempC2 = C[i+1 + j*lda];
 151 //        for(k=0; k < lda; k+=1) {
 152 //          //C[k   + jlda] += tempA1 * B[k   + i*lda] + tempA2 * B[k   + (i+1)*lda] + tempA3 * B[k   + (i+2)*lda] + tempA4 * B[k   + (i+3)*lda] +
 153 //          //                 tempA5 * B[k   + (i+4)*lda] + tempA6 * B[k   + (i+5)*lda] + tempA7 * B[k   + (i+6)*lda] + tempA8 * B[k   + (i+7)*lda];
 154 //
 155 //          C[k   + jlda] += tempA1* B[k   + i*lda];// + ffmul(tempA2,B[k   + (i+1)*lda]) + tempA3 * B[k   + (i+2)*lda] + tempA4 * B[k   + (i+3)*lda] +
 156 //          //                 tempA5 * B[k   + (i+4)*lda] + tempA6 * B[k   + (i+5)*lda] + tempA7 * B[k   + (i+6)*lda] + tempA8 * B[k   + (i+7)*lda];
 157 //          //
 158 //          //C[k+1 + jlda] += tempA1 * B[k+1 + i*lda] + tempA2 * B[k+1 + (i+1)*lda] + tempA3 * B[k+1 + (i+2)*lda] + tempA4 * B[k+1 + (i+3)*lda] +
 159 //          //                 tempA5 * B[k+1   + (i+4)*lda] + tempA6 * B[k+1   + (i+5)*lda] + tempA7 * B[k+1   + (i+6)*lda] + tempA8 * B[k+1   + (i+7)*lda];
 160 //          //
 161 //          //C[k+2 + jlda] += tempA1 * B[k+2 + i*lda] + tempA2 * B[k+2 + (i+1)*lda] + tempA3 * B[k+2 + (i+2)*lda] + tempA4 * B[k+2 + (i+3)*lda] +
 162 //          //                 tempA5 * B[k+2   + (i+4)*lda] + tempA6 * B[k+2   + (i+5)*lda] + tempA7 * B[k+2   + (i+6)*lda] + tempA8 * B[k+2   + (i+7)*lda];
 163 //          //
 164 //          //C[k+3 + jlda] += tempA1 * B[k+3 + i*lda] + tempA2 * B[k+3 + (i+1)*lda] + tempA3 * B[k+3 + (i+2)*lda] + tempA4 * B[k+3 + (i+3)*lda] +
 165 //          //                 tempA5 * B[k+3   + (i+4)*lda] + tempA6 * B[k+3   + (i+5)*lda] + tempA7 * B[k+3   + (i+6)*lda] + tempA8 * B[k+3   + (i+7)*lda];
 166 //          //
 167 //          //C[k+4 + jlda] += tempA1 * B[k+4 + i*lda] + tempA2 * B[k+4 + (i+1)*lda] + tempA3 * B[k+4 + (i+2)*lda] + tempA4 * B[k+4 + (i+3)*lda] +
 168 //          //                 tempA5 * B[k+4   + (i+4)*lda] + tempA6 * B[k+4   + (i+5)*lda] + tempA7 * B[k+4   + (i+6)*lda] + tempA8 * B[k+4   + (i+7)*lda];
 169 //          //
 170 //          //C[k+5 + jlda] += tempA1 * B[k+5 + i*lda] + tempA2 * B[k+5 + (i+1)*lda] + tempA3 * B[k+5 + (i+2)*lda] + tempA4 * B[k+5 + (i+3)*lda] +
 171 //          //                 tempA5 * B[k+5   + (i+4)*lda] + tempA6 * B[k+5   + (i+5)*lda] + tempA7 * B[k+5   + (i+6)*lda] + tempA8 * B[k+5   + (i+7)*lda];
 172 //          //
 173 //          //C[k+6 + jlda] += tempA1 * B[k+6 + i*lda] + tempA2 * B[k+6 + (i+1)*lda] + tempA3 * B[k+6 + (i+2)*lda] + tempA4 * B[k+6 + (i+3)*lda] +
 174 //          //                 tempA5 * B[k+6   + (i+4)*lda] + tempA6 * B[k+6   + (i+5)*lda] + tempA7 * B[k+6   + (i+6)*lda] + tempA8 * B[k+6   + (i+7)*lda];
 175 //          //
 176 //          //C[k+7 + jlda] += tempA1 * B[k+7 + i*lda] + tempA2 * B[k+7 + (i+1)*lda] + tempA3 * B[k+7 + (i+2)*lda] + tempA4 * B[k+7 + (i+3)*lda] +
 177 //          //                 tempA5 * B[k+7   + (i+4)*lda] + tempA6 * B[k+7   + (i+5)*lda] + tempA7 * B[k+7   + (i+6)*lda] + tempA8 * B[k+7   + (i+7)*lda];
 178 //
 179 //
 180 //        }
 181 //      }
 182 //    }
 183 //}
 184
 185
 186 void __attribute__((noinline)) matmul(const int lda,  const data_t A[], const data_t B[], data_t C[] )
 187 {
 188
 189    // ***************************** //
 190    // **** ADD YOUR CODE HERE ***** //
 191    // ***************************** //
 192    //
 193    // feel free to make a separate function for MI and MSI versions.
 194
 195     static __thread int i, j, k;
 196     static __thread data_t tempA0, tempA1, tempA2, tempA3, tempA4, tempA5, tempA6, tempA7;
 197     static __thread data_t tempC0, tempC1, tempC2, tempC3, tempC4, tempC5, tempC6, tempC7; //tempC8, tempC9, tempC10, tempC11, tempC12, tempC13, tempC14, tempC15;
 198
 199     static __thread int start, end, jStride, jToRow, jToCol, iToRow;
 200
 201     start = coreid << 9;
 202     end = (coreid+1) << 9;
 203     jStride = 8;
 204
 205     for (j=start; j < end; j+=jStride) {
 206       jToRow = (j>>5)<<5;
 207       jToCol = j%32;
 208       tempC0  = 0;
 209       tempC1  = 0;
 210       tempC2  = 0;
 211       tempC3  = 0;
 212       tempC4  = 0;
 213       tempC5  = 0;
 214       tempC6  = 0;
 215       tempC7  = 0;
 216       //tempC8  = 0;
 217       //tempC9  = 0;
 218       //tempC10 = 0;
 219       //tempC11 = 0;
 220       //tempC12 = 0;
 221       //tempC13 = 0;
 222       //tempC14 = 0;
 223       //tempC15 = 0;
 224
 225       for ( i=0; i < lda; i+=2 ) {
 226         iToRow = i << 5;
 227
 228         tempA0 = A[i   + jToRow];
 229         tempA1 = A[i+1 + jToRow];
 230         //tempA2 = A[i+2 + jToRow];
 231         //tempA3 = A[i+3 + jToRow];
 232         //tempA4 = A[i+4 + jToRow];
 233         //tempA5 = A[i+5 + jToRow];
 234         //tempA6 = A[i+6 + jToRow];
 235         //tempA7 = A[i+7 + jToRow];
 236
 237         tempC0  += tempA0 * B[(jToCol   ) + (iToRow)];
 238         tempC1  += tempA0 * B[(jToCol+1 ) + (iToRow)];
 239         tempC2  += tempA0 * B[(jToCol+2 ) + (iToRow)];
 240         tempC3  += tempA0 * B[(jToCol+3 ) + (iToRow)];
 241         tempC4  += tempA0 * B[(jToCol+4 ) + (iToRow)];
 242         tempC5  += tempA0 * B[(jToCol+5 ) + (iToRow)];
 243         tempC6  += tempA0 * B[(jToCol+6 ) + (iToRow)];
 244         tempC7  += tempA0 * B[(jToCol+7 ) + (iToRow)];
 245         //tempC8  += tempA0 * B[(jToCol+8 ) + (iToRow)];
 246         //tempC9  += tempA0 * B[(jToCol+9 ) + (iToRow)];
 247         //tempC10 += tempA0 * B[(jToCol+10) + (iToRow)];
 248         //tempC11 += tempA0 * B[(jToCol+11) + (iToRow)];
 249         //tempC12 += tempA0 * B[(jToCol+12) + (iToRow)];
 250         //tempC13 += tempA0 * B[(jToCol+13) + (iToRow)];
 251         //tempC14 += tempA0 * B[(jToCol+14) + (iToRow)];
 252         //tempC15 += tempA0 * B[(jToCol+15) + (iToRow)];
 253
 254         iToRow += 32;
 255         tempC0  += tempA1 * B[(jToCol   ) + (iToRow)];
 256         tempC1  += tempA1 * B[(jToCol+1 ) + (iToRow)];
 257         tempC2  += tempA1 * B[(jToCol+2 ) + (iToRow)];
 258         tempC3  += tempA1 * B[(jToCol+3 ) + (iToRow)];
 259         tempC4  += tempA1 * B[(jToCol+4 ) + (iToRow)];
 260         tempC5  += tempA1 * B[(jToCol+5 ) + (iToRow)];
 261         tempC6  += tempA1 * B[(jToCol+6 ) + (iToRow)];
 262         tempC7  += tempA1 * B[(jToCol+7 ) + (iToRow)];
 263         //tempC8  += tempA1 * B[(jToCol+8 ) + (iToRow+32)];
 264         //tempC9  += tempA1 * B[(jToCol+9 ) + (iToRow+32)];
 265         //tempC10 += tempA1 * B[(jToCol+10) + (iToRow+32)];
 266         //tempC11 += tempA1 * B[(jToCol+11) + (iToRow+32)];
 267         //tempC12 += tempA1 * B[(jToCol+12) + (iToRow+32)];
 268         //tempC13 += tempA1 * B[(jToCol+13) + (iToRow+32)];
 269         //tempC14 += tempA1 * B[(jToCol+14) + (iToRow+32)];
 270         //tempC15 += tempA1 * B[(jToCol+15) + (iToRow+32)];
 271
 272         //iToRow += 32;
 273         //tempC0  += tempA2 * B[(jToCol   ) + (iToRow)];
 274         //tempC1  += tempA2 * B[(jToCol+1 ) + (iToRow)];
 275         //tempC2  += tempA2 * B[(jToCol+2 ) + (iToRow)];
 276         //tempC3  += tempA2 * B[(jToCol+3 ) + (iToRow)];
 277         //tempC4  += tempA2 * B[(jToCol+4 ) + (iToRow)];
 278         //tempC5  += tempA2 * B[(jToCol+5 ) + (iToRow)];
 279         //tempC6  += tempA2 * B[(jToCol+6 ) + (iToRow)];
 280         //tempC7  += tempA2 * B[(jToCol+7 ) + (iToRow)];
 281         //tempC8  += tempA2 * B[(jToCol+8 ) + (iToRow)];
 282         //tempC9  += tempA2 * B[(jToCol+9 ) + (iToRow)];
 283         //tempC10 += tempA2 * B[(jToCol+10) + (iToRow)];
 284         //tempC11 += tempA2 * B[(jToCol+11) + (iToRow)];
 285         //tempC12 += tempA2 * B[(jToCol+12) + (iToRow)];
 286         //tempC13 += tempA2 * B[(jToCol+13) + (iToRow)];
 287         //tempC14 += tempA2 * B[(jToCol+14) + (iToRow)];
 288         //tempC15 += tempA2 * B[(jToCol+15) + (iToRow)];
 289
 290         //iToRow += 32;
 291         //tempC0  += tempA3 * B[(jToCol   ) + (iToRow)];
 292         //tempC1  += tempA3 * B[(jToCol+1 ) + (iToRow)];
 293         //tempC2  += tempA3 * B[(jToCol+2 ) + (iToRow)];
 294         //tempC3  += tempA3 * B[(jToCol+3 ) + (iToRow)];
 295         //tempC4  += tempA3 * B[(jToCol+4 ) + (iToRow)];
 296         //tempC5  += tempA3 * B[(jToCol+5 ) + (iToRow)];
 297         //tempC6  += tempA3 * B[(jToCol+6 ) + (iToRow)];
 298         //tempC7  += tempA3 * B[(jToCol+7 ) + (iToRow)];
 299         //tempC8  += tempA3 * B[(jToCol+8 ) + (iToRow)];
 300         //tempC9  += tempA3 * B[(jToCol+9 ) + (iToRow)];
 301         //tempC10 += tempA3 * B[(jToCol+10) + (iToRow)];
 302         //tempC11 += tempA3 * B[(jToCol+11) + (iToRow)];
 303         //tempC12 += tempA3 * B[(jToCol+12) + (iToRow)];
 304         //tempC13 += tempA3 * B[(jToCol+13) + (iToRow)];
 305         //tempC14 += tempA3 * B[(jToCol+14) + (iToRow)];
 306         //tempC15 += tempA3 * B[(jToCol+15) + (iToRow)];
 307
 308         //iToRow += 32;
 309         //tempC0 += tempA4 * B[(jToCol   ) + (iToRow)];
 310         //tempC1 += tempA4 * B[(jToCol+1 ) + (iToRow)];
 311         //tempC2 += tempA4 * B[(jToCol+2 ) + (iToRow)];
 312         //tempC3 += tempA4 * B[(jToCol+3 ) + (iToRow)];
 313         //tempC4 += tempA4 * B[(jToCol+4 ) + (iToRow)];
 314         //tempC5 += tempA4 * B[(jToCol+5 ) + (iToRow)];
 315         //tempC6 += tempA4 * B[(jToCol+6 ) + (iToRow)];
 316         //tempC7 += tempA4 * B[(jToCol+7 ) + (iToRow)];
 317         //
 318         //iToRow += 32;
 319         //tempC0 += tempA5 * B[(jToCol   ) + (iToRow)];
 320         //tempC1 += tempA5 * B[(jToCol+1 ) + (iToRow)];
 321         //tempC2 += tempA5 * B[(jToCol+2 ) + (iToRow)];
 322         //tempC3 += tempA5 * B[(jToCol+3 ) + (iToRow)];
 323         //tempC4 += tempA5 * B[(jToCol+4 ) + (iToRow)];
 324         //tempC5 += tempA5 * B[(jToCol+5 ) + (iToRow)];
 325         //tempC6 += tempA5 * B[(jToCol+6 ) + (iToRow)];
 326         //tempC7 += tempA5 * B[(jToCol+7 ) + (iToRow)];
 327         //
 328         //iToRow += 32;
 329         //tempC0 += tempA6 * B[(jToCol   ) + (iToRow)];
 330         //tempC1 += tempA6 * B[(jToCol+1 ) + (iToRow)];
 331         //tempC2 += tempA6 * B[(jToCol+2 ) + (iToRow)];
 332         //tempC3 += tempA6 * B[(jToCol+3 ) + (iToRow)];
 333         //tempC4 += tempA6 * B[(jToCol+4 ) + (iToRow)];
 334         //tempC5 += tempA6 * B[(jToCol+5 ) + (iToRow)];
 335         //tempC6 += tempA6 * B[(jToCol+6 ) + (iToRow)];
 336         //tempC7 += tempA6 * B[(jToCol+7 ) + (iToRow)];
 337         //
 338         //iToRow += 32;
 339         //tempC0 += tempA7 * B[(jToCol   ) + (iToRow)];
 340         //tempC1 += tempA7 * B[(jToCol+1 ) + (iToRow)];
 341         //tempC2 += tempA7 * B[(jToCol+2 ) + (iToRow)];
 342         //tempC3 += tempA7 * B[(jToCol+3 ) + (iToRow)];
 343         //tempC4 += tempA7 * B[(jToCol+4 ) + (iToRow)];
 344         //tempC5 += tempA7 * B[(jToCol+5 ) + (iToRow)];
 345         //tempC6 += tempA7 * B[(jToCol+6 ) + (iToRow)];
 346         //tempC7 += tempA7 * B[(jToCol+7 ) + (iToRow)];
 347
 348       }
 349       C[j     ] = tempC0;
 350       C[j + 1 ] = tempC1;
 351       C[j + 2 ] = tempC2;
 352       C[j + 3 ] = tempC3;
 353       C[j + 4 ] = tempC4;
 354       C[j + 5 ] = tempC5;
 355       C[j + 6 ] = tempC6;
 356       C[j + 7 ] = tempC7;
 357       //C[j + 8 ] = tempC8 ;
 358       //C[j + 9 ] = tempC9 ;
 359       //C[j + 10] = tempC10;
 360       //C[j + 11] = tempC11;
 361       //C[j + 12] = tempC12;
 362       //C[j + 13] = tempC13;
 363       //C[j + 14] = tempC14;
 364       //C[j + 15] = tempC15;
 365     }
 366 }
 367
 368
 369
 370 //--------------------------------------------------------------------------
 371 // Main
 372 //
 373 // all threads start executing thread_entry(). Use their "coreid" to
 374 // differentiate between threads (each thread is running on a separate core).
 375
 376 void thread_entry(int cid, int nc)
 377 {
 378    coreid = cid;
 379    ncores = nc;
 380
 381    // static allocates data in the binary, which is visible to both threads
 382    static data_t results_data[ARRAY_SIZE];
 383
 384
 385    //// Execute the provided, naive matmul
 386    //barrier();
 387    //stats(matmul_naive(DIM_SIZE, input1_data, input2_data, results_data); barrier());
 388
 389    //
 390    //// verify
 391    //verify(ARRAY_SIZE, results_data, verify_data);
 392    //
 393    //// clear results from the first trial
 394    //size_t i;
 395    //if (coreid == 0)
 396    //   for (i=0; i < ARRAY_SIZE; i++)
 397    //      results_data[i] = 0;
 398    //barrier();
 399
 400
 401    // Execute your faster matmul
 402    barrier();
 403    stats(matmul(DIM_SIZE, input1_data, input2_data, results_data); barrier());
 404
 405 #ifdef DEBUG
 406    printArray("results:", ARRAY_SIZE, results_data);
 407    printArray("verify :", ARRAY_SIZE, verify_data);
 408 #endif
 409
 410    // verify
 411    verify(ARRAY_SIZE, results_data, verify_data);
 412    barrier();
 413
 414    exit(0);
 415 }
 416