matrix3.c

  for (bi=0; bi<size; bi+=block_size){
    for (bj=0; bj<size; bj+=block_size){
      //Handle cases where size%block_size!=0
      max_j = (bi+block_size<=size)?block_size:size%block_size;
      max_k = (bj+block_size<=size)?block_size:size%block_size;
      for (i=0; i<block_size; i++){
	for (j=0; j<max_j; j++){
	  for (k=0; k<max_k; k++){
	    result[i*size+(bi+j)] += A[i*size+(bj+k)] * B[(bj+k)*size+(bi+j)];
	  }
	}
      }
    }
  }