mt/dc_matmul.c

   1 #include "stdlib.h"
   2
   3 #include "util.h"
   4
   5 #include "dataset.h"
   6
   7 #define REG_I 8
   8 #define REG_J 2
   9 #define BLOCK_I 32
  10 #define BLOCK_J 16
  11 #define BLOCK_K 16
  12 #define LDA 32
  13 #define NCORES 2
  14 #define MIN(X,Y) (X < Y ? X : Y)
  15
  16 void __attribute__((noinline)) matmul(const int coreid, const int ncores, const int lda,  const data_t A[], const data_t B[], data_t C[] )
  17 {
  18
  19    // ***************************** //
  20    // **** ADD YOUR CODE HERE ***** //
  21    // ***************************** //
  22    //
  23    // feel free to make a separate function for MI and MSI versions.
  24
  25   int i, j, k, ri, rj, ii, jj, kk;
  26   data_t *Aj, *Cj, *Bi;
  27   data_t c[REG_I][REG_J], a[REG_J], b[REG_I];
  28   size_t start = coreid * (LDA / NCORES), end = (coreid == NCORES - 1 ? LDA : (coreid + 1) * (LDA / NCORES));
  29
  30   /* if (coreid > 0) { */
  31   /*   return; */
  32   /* } */
  33   /* start = 0, end = lda; */
  34   if (ncores == NCORES && lda == LDA) {
  35     for (jj = start; jj < end; jj += BLOCK_J) {
  36       int kk_start= (coreid == 0 ? 0 : LDA/2) ,kk_end = (coreid == 0 ? LDA/2 : LDA);
  37       for (kk = kk_start; kk < kk_end; kk += BLOCK_K) {
  38         //  for (ii = 0; ii < LDA; ii += BLOCK_I)
  39         for (j = jj; j < MIN(end, jj + BLOCK_J); j += REG_J) {
  40           Aj = A + j*LDA;
  41           Cj = C + j*LDA;
  42           for (i = 0; i < LDA/*, ii + BLOCK_I)*/; i += REG_I) {
  43             /* Load C in register blocks. */
  44             Bi = B + i;
  45             for (ri = 0; ri < REG_I; ri++) {
  46               for (rj = 0; rj < REG_J; rj++) {
  47                 c[ri][rj] = Cj[i + ri + ( rj)*LDA];
  48               }
  49             }
  50
  51
  52             for (k = kk; k < MIN(LDA, kk + BLOCK_K); k++) {
  53               for (ri = 0; ri < REG_I; ri++) {
  54                 b[ri] = Bi[k*LDA  + ri];
  55               }
  56               /* Compute C in register blocks. */
  57               for (rj = 0; rj < REG_J; rj++) {
  58                 a[rj] = Aj[(rj)*LDA + k];
  59                 for (ri = 0; ri < REG_I; ri++) {
  60                   c[ri][rj] += a[rj] * b[ri];
  61                 }
  62               }
  63             }
  64
  65             /* store C in register blocks. */
  66             for (ri = 0; ri < REG_I; ri++) {
  67               for (rj = 0; rj < REG_J; rj++) {
  68                 Cj[i + ri + ( rj)*LDA] = c[ri][rj];
  69               }
  70             }
  71           }
  72         }
  73         /* barrier(nc); */
  74
  75         /* kk_start= (coreid == 1 ? 0 : LDA/2); */
  76         /* kk_end = (coreid == 1 ? LDA/2 : LDA); */
  77         /* for (kk = kk_start; kk < kk_end; kk += BLOCK_K) { */
  78         /* //  for (ii = 0; ii < LDA; ii += BLOCK_I) */
  79         /* for (j = jj; j < MIN(end, jj + BLOCK_J); j += REG_J) { */
  80         /*   Aj = A + j*LDA; */
  81         /*   Cj = C + j*LDA; */
  82         /*   for (i = 0; i < LDA/\*, ii + BLOCK_I)*\/; i += REG_I) { */
  83         /*     /\* Load C in register blocks. *\/ */
  84         /*     Bi = B + i; */
  85         /*     for (ri = 0; ri < REG_I; ri++) { */
  86         /*       for (rj = 0; rj < REG_J; rj++) { */
  87         /*      c[ri][rj] = Cj[i + ri + ( rj)*LDA]; */
  88         /*       } */
  89         /*     } */
  90
  91
  92         /*     for (k = kk; k < MIN(LDA, kk + BLOCK_K); k++) { */
  93         /*       for (ri = 0; ri < REG_I; ri++) { */
  94         /*      b[ri] = Bi[k*LDA  + ri]; */
  95         /*       } */
  96         /*       /\* Compute C in register blocks. *\/ */
  97         /*       for (rj = 0; rj < REG_J; rj++) { */
  98         /*      a[rj] = Aj[(rj)*LDA + k]; */
  99         /*      for (ri = 0; ri < REG_I; ri++) { */
 100         /*        c[ri][rj] += a[rj] * b[ri]; */
 101         /*      } */
 102         /*       } */
 103         /*     } */
 104
 105             /* store C in register blocks. */
 106         /*     for (ri = 0; ri < REG_I; ri++) { */
 107     /*        for (rj = 0; rj < REG_J; rj++) { */
 108     /*          Cj[i + ri + ( rj)*LDA] = c[ri][rj]; */
 109     /*        } */
 110     /*      } */
 111     /*   } */
 112     /* } */
 113       }
 114     }
 115
 116
 117     //barrier(nc);
 118     for (jj = start; jj < end; jj += BLOCK_J) {
 119       int kk_start= (coreid != 0 ? 0 : LDA/2), kk_end = (coreid != 0 ? LDA/2 : LDA);
 120       for (kk = kk_start; kk < kk_end; kk += BLOCK_K) {
 121         //  for (ii = 0; ii < LDA; ii += BLOCK_I)
 122         for (j = jj; j < MIN(end, jj + BLOCK_J); j += REG_J) {
 123           Aj = A + j*LDA;
 124           Cj = C + j*LDA;
 125           for (i = 0; i < LDA/*, ii + BLOCK_I)*/; i += REG_I) {
 126             /* Load C in register blocks. */
 127             Bi = B + i;
 128             for (ri = 0; ri < REG_I; ri++) {
 129               for (rj = 0; rj < REG_J; rj++) {
 130                 c[ri][rj] = Cj[i + ri + ( rj)*LDA];
 131               }
 132             }
 133
 134
 135             for (k = kk; k < MIN(LDA, kk + BLOCK_K); k++) {
 136               for (ri = 0; ri < REG_I; ri++) {
 137                 b[ri] = Bi[k*LDA  + ri];
 138               }
 139               /* Compute C in register blocks. */
 140               for (rj = 0; rj < REG_J; rj++) {
 141                 a[rj] = Aj[(rj)*LDA + k];
 142                 for (ri = 0; ri < REG_I; ri++) {
 143                   c[ri][rj] += a[rj] * b[ri];
 144                 }
 145               }
 146             }
 147
 148               /* store C in register blocks. */
 149             for (ri = 0; ri < REG_I; ri++) {
 150               for (rj = 0; rj < REG_J; rj++) {
 151                 Cj[i + ri + ( rj)*LDA] = c[ri][rj];
 152               }
 153             }
 154           }
 155         }
 156       }
 157     }
 158     /* We only care about performance for 32x32 matrices and 2 cores. Otherwise just naive mat_mul */
 159 } else {
 160       if (coreid > 0)
 161         return;
 162
 163       for ( i = 0; i < lda; i++ )
 164         for ( j = 0; j < lda; j++ )
 165           for ( k = 0; k < lda; k++ )
 166             C[i + j*lda] += A[j*lda + k] * B[k*lda + i];
 167     }
 168   }