Parallel Programming for FPGAs 矩阵乘法 阅读笔记

621 阅读3分钟

第七章 矩阵乘法

矩阵乘法

void matrixmul(int A[N][M], int B[M][P], int AB[N][P]) {
  #pragma HLS ARRAY_RESHAPE variable=A complete dim=2
  #pragma HLS ARRAY_RESHAPE variable=B complete dim=1
  /* for each row and column of AB */
  row: for(int i = 0; i < N; ++i) {
    col: for(int j = 0; j < P; ++j) {
      #pragma HLS PIPELINE II=1
      /* compute (AB)i,j */
      int ABij = 0;
    product: for(int k = 0; k < M; ++k) {
        ABij += A[i][k] * B[k][j];
      }
      AB[i][j] = ABij;
    }
  }
}

pipeline

directive放在函数 面积(MACs) 速度(interval)
最顶层 N × M × P 1
row M × P N
最内层 1 N × M × P(II=1:possible but difficult because of the recurrence involved in the accumulation of variable ABij.)

array_partition

每个时钟周期执行大量的运算,就需要为它们提供运算的对象,保存每个运算出处的结果。

当数组的分块越来越小,每次存储的访问顺序在编译时就可以确定。所以数组分块是最简单和高效的方法用于提升在单位时间内,对内存的访问。

array_partition可以:

  • 将存储存储空间的地址分为不同的区间
  • 也可以将多个存储空间合并以一个。

这种变化存储空间数据位宽变大用于保存数组,但是不能改变总存储的比特数目。这其中的区别如图7.2所示。

图7.2: A[N][M]

  • 左边: 采用原始的方式,保存N*M个元素。
  • 中间: 使用了array_partition directive进行了变形,结果是使用了M个存储单元,每个保存N个元素。
  • 右边: 使用array_shape directive进行了变形,结果保存在一个存储单元,带有N个入口,每个入口保存原始矩阵的M个元素

分块维度的确定: 关键看哪个维度需要一个时钟内对该维度的多个元素进行读写

note: 不可分割的维度通常相邻操作的元素那个维度是相同的,即通过相同的变量进行访问。这里如果pipeline之后,A的列、B的行是通过常量访问的,故应当对这个维度进行partition。

注意这里的分块是有效的因为分块的维度(A 的列和 B 的行)是通过常量访问的。另外,不可分块的维度是通过相同变量访问的。当对多维矩阵进行分块是,这是一条很用用的经验用于确定哪个维度可以进行分块。

块矩阵乘法

streaming architecture

  • 减少了输入数据和输出数据对存储的需求

  • 头文件
#ifndef _BLOCK_MM_H_
#define _BLOCK_MM_H_
#include "hls_stream.h"
#include <iostream>
#include <iomanip>
#include <vector>
using namespace std;

typedef int DTYPE;
const int SIZE = 8;
const int BLOCK_SIZE = 4;

typedef struct { DTYPE a[BLOCK_SIZE]; } blockvec;

typedef struct { DTYPE out[BLOCK_SIZE][BLOCK_SIZE]; } blockmat;

void blockmatmul(hls::stream<blockvec> &Arows, hls::stream<blockvec> &Bcols,
                                 blockmat & ABpartial, DTYPE iteration);
#endif
  • 源文件
#include "block_mm.h"
void blockmatmul(hls::stream<blockvec> &Arows, hls::stream<blockvec> &Bcols,
        blockmat &ABpartial, int it) {
  #pragma HLS DATAFLOW
  int counter = it % (SIZE/BLOCK_SIZE);
  static DTYPE A[BLOCK_SIZE][SIZE];
  if(counter == 0){ //only load the A rows when necessary
    loadA: for(int i = 0; i < SIZE; i++) {
      blockvec tempA = Arows.read();
      for(int j = 0; j < BLOCK_SIZE; j++) {
        #pragma HLS PIPELINE II=1
        A[j][i] = tempA.a[j];
      }
    }
  }
  DTYPE AB[BLOCK_SIZE][BLOCK_SIZE] = { 0 };
  partialsum: for(int k=0; k < SIZE; k++) {
    blockvec tempB = Bcols.read();
    for(int i = 0; i < BLOCK_SIZE; i++) {
      for(int j = 0; j < BLOCK_SIZE; j++) {
        AB[i][j] = AB[i][j] +  A[i][k] * tempB.a[j];
      }
    }
  }
  writeoutput: for(int i = 0; i < BLOCK_SIZE; i++) {
    for(int j = 0; j < BLOCK_SIZE; j++) {
      ABpartial.out[i][j] = AB[i][j];
    }
  }
}
  • 测试文件
// The beginning of the testbench is shown in the previous figure
int main() {
  int row, col, it = 0;
  for(int it1 = 0; it1 < SIZE; it1 = it1 + BLOCK_SIZE) {
    for(int it2 = 0; it2 < SIZE; it2 = it2 + BLOCK_SIZE) {
      row = it1; //row + BLOCK_SIZE * factor_row;
      col = it2; //col + BLOCK_SIZE * factor_col;

      for(int k = 0; k < SIZE; k++) {
        for(int i = 0; i < BLOCK_SIZE; i++) {
          if(it % (SIZE/BLOCK_SIZE) == 0) strm_matrix1_element.a[i] = A[row+i][k];
          strm_matrix2_element.a[i] = B[k][col+i];
        }
        if(it % (SIZE/BLOCK_SIZE) == 0) strm_matrix1.write(strm_matrix1_element);
        strm_matrix2.write(strm_matrix2_element);
      }
      blockmatmul(strm_matrix1, strm_matrix2, block_out, it);

      for(int i = 0; i < BLOCK_SIZE; i++)
        for(int j = 0; j < BLOCK_SIZE; j++)
          matrix_hwout[row+i][col+j] = block_out.out[i][j];
      it = it + 1;
    }
  }

  matmatmul_sw(A, B, matrix_swout);

  for(int i = 0; i<SIZE; i++)
    for(int j = 0; j<SIZE; j++)
      if(matrix_swout[i][j] != matrix_hwout[i][j]) { fail=1; }

  if(fail==1) cout << "failed" << endl;
  else cout << "passed" << endl;

  return 0;
}