c++ - 块矩阵乘法

Question

我想执行块矩阵乘法（将矩阵划分为多个 sxs 矩阵并乘以相应的块）。我按照轩尼诗建筑书的示例代码编写了代码：

for(int jj=0;jj<=(n/s);jj += s){
            for(int kk=1;kk<=(n/s);kk += s){
                    for(int i=1;i<=(n/s);i++){
                            for(int j = jj; j<=((jj+s-1)>(n/s)?(n/s):(jj+s-1)); j++){
                                    temp = 0;
                                    for(int k = kk; k<=((kk+s-1)>(n/s)?(n/s):(kk+s-1)); k++){
                                            temp += b[i][k]*a[k][j];
                                    }
                                    c[j][i] += temp;
                            }
                    }
            }
    }  

这里，nxn 是原始矩阵的大小。a, b 矩阵大小相同。我将 a,b 矩阵划分为大小为 sxs 的块。在我的程序中，我将块大小设置为 4。我将 a、b 的所有元素设置为 5、a 和 n = 1000。但是，我的结果中的值是错误的。我在这里做错了吗？从过去的 2 个小时开始就一直坚持这一点。如果可能的话，你们能帮忙吗？书中的参考代码是这样的：

for (jj = 0; jj <= size; jj += N) {
    for (kk = 1; kk <= size; kk += N) {
        for (i = 1; i <= size; i++) {
            for (j = jj; j <= findMin(jj+N-1, size); j++) {
                temp = 0;
                for (k = kk; k <= findMin(kk+N-1, size); k++) {
                    temp += B[i][k] * A[j][k];
                }
                C[j][i] += temp;
            }
        }
     }
}  

这里，s=N 和 size = n/s

Answer 1

for(int jj=0;jj<N;jj+= s){
        for(int kk=0;kk<N;kk+= s){
                for(int i=0;i<N;i++){
                        for(int j = jj; j<((jj+s)>N?N:(jj+s)); j++){
                                temp = 0;
                                for(int k = kk; k<((kk+s)>N?N:(kk+s)); k++){
                                        temp += a[i][k]*b[k][j];
                                }
                                c[i][j] += temp;
                        }
                }
        }
}

AxB 大小为 N

Answer 2

错误在这一行。你有

 temp += b[i][k]*a[k][j];

你应该有

 temp += b[i][k]*a[j][k];

反而。

如果您可以将这部分放在一个函数中而不是这一行中，那就更好了：

((jj+s-1)>(n/s)?(n/s):(jj+s-1));

Answer 3

乍一看，我很惊讶地看到 0 和 1 起始索引和 <= 循环终止测试。带有 fortran 或 matlab 代码的书籍有时会假设基于 1 的索引，而 c/c++ 使用基于 0 的索引。

您还可以分别实现和/或测试内部两个 for 循环，因为它们将用于单块矩阵乘法。

我会将 findMin 函数分开，而不是在循环测试中内联它。