我想加速倍增打击:
这是我需要改进的代码:
void multiply(int size, int **matA, int **matB, int ** matC) {
for(int i=0;i<size;i++) {
for(int j=0;j<size;j++) {
int t = matC[i][j];
for(int k=0;k<size;k++) {
t += matA[i][k] * matB[k][j];
}
matC[i][j] += t;
}
}
}
我有两个矩阵和一个大小为 5000x5000 的结果矩阵。 巨大的矩阵可能意味着它们无法完全加载到缓存中? for循环中是否出现过多页错误?我想知道如何加快乘法速度,以及如何组织数据(使用一维数组还是二维数组?)
我的答案代码是列表blow,我选择使用1d数组来模拟2d数组(每个矩阵使用
new []
一次),但我不确定使用2d数组时是否更快。
我使用临时矩阵来存储 matB
的转置矩阵,以避免 for 循环中的页面错误。
我添加 AVX2 以获得更高的性能。
使用大小为 5000x5000 的一维数组或二维数组哪个更好? 还有其他想法或技巧吗?
int** allocate(int rows, int cols) {
int ** mat;
mat = new int*[rows];
int *temp = new int[rows*cols];
for(int i = 0; i<rows; i++) {
mat[i] = temp + i * cols;
}
return mat;
}
void multiply(int size, int **matA, int **matB, int ** matC) {
int i;
int n = size*size; // total size
// column-major order
int **transMatB = allocate(size, size);
int *transArrB = transMatB[0];
//copy transposed data, maybe many page faults here.
#pragma omp parallel for
for(i = 0; i < n; i++) {
transMatB[i/size][i%size] = matB[i%size][i/size];
}
#pragma omp parallel for
for(i = 0; i < n; i ++) {
int *row = matA[i / size];
int *col = transMatB[i % size];
int temp;
#ifdef __AVX2__
temp = multiplyAndSumArrays(row, col, size);
#else
temp = 0;
for (int k = 0; k < size; k ++) {
temp += row[k] * col[k];
}
#endif
matC[i / size][i % size] += temp;
}
// remove temp transposed mastrix
delete[] transMatB[0];
delete[] transMatB; transMatB = nullptr;
}
在优化方面,矩阵-矩阵乘法可能是研究最多的内核。对于最终结果,请阅读 Goto 和 van de Geijn 的论文,引用如下。
关键在于
特别是最后一点:简而言之,3 个循环中的每一个都分为两个循环,一个在块上,一个在块内。然后你有 6 个循环(意味着 5 个!左右不同的算法)和 3 个块大小作为调整参数。上面的论文对此分析得很完整。
请注意,这并不简单!对于合理可行的解决方案,请执行递归 2x2 乘法:将每个矩阵划分为 2x2 块结构,然后递归地相乘。当块足够小以适合缓存时,您将停止递归。
这应该可以作为课堂作业来实现,从而提高成绩。您甚至可以简单地进行多线程处理。
Goto, Kazushige / Geijn, Robert A. van de
Anatomy of high-performance matrix multiplication
2008
ACM Trans. Math. Softw. , Vol. 34, No. 3
ACM: New York, NY, USA p. 1-25