缓存友好的矩阵乘法方法

Question

缓存友好的矩阵乘法方法

matrix-multiplicationcachegrind

9

我打算使用高缓存友好性的方法来计算两个矩阵的乘积（以减少cache miss的次数）。

我发现可以通过缓存友好的转置函数实现此目的。

但是，我无法找到这个算法。请问如何实现？

- Aakash Anuj

2个回答

1

@Cesar的答案是不正确的。例如，内部循环。

for (int k = 0; k < n; k++)
   s += a[i,k] * Bcolj[k];

遍历 a 的第 i 列。

以下代码应确保我们始终按行访问数据。

void multiply(const double (&a)[I][K], 
              const double (&b)[K][J], 
              double (&c)[I][J]) 
{
    for (int j=0; j<J; ++j) {
       // iterates the j-th row of c
       for (int i=0; i<I; ++i) {
         c[i][j] = 0;
       } 

       // iterates the j-th row of b
       for (int k=0; k<K; ++k) {
          double t = b[k][j];
          // iterates the j-th row of c
          // iterates the k-th row of a
          for (int i=0; i<I; ++i) {
            c[i][j] += a[i][k] * t;
          } 
       }
    }
}

- Joe C

2

你的代码也有问题。 c[i][j] 的重置完全是可选的（这取决于调用者是否将矩阵重置为零）。此外，k 的循环从1开始，但应该从零开始。 - greywolf82

@greywolf82 需要重置 c[i][j]，因为 "c[i][j] += a[i][k] * t;" 的累加需要一个初始值。"k 从0开始" 是正确的。已修复。 - Joe C

是的，我知道，但如果调用者例如对零进行了memset，则不需要循环。在您的代码中添加注释以澄清。 - greywolf82

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- Cesar · Accepted Answer

你要找的词是抖动。在谷歌上搜索抖动矩阵乘法可以得到更多结果。

标准的乘法算法c=a*b看起来像这样：

void multiply(double[,] a, double[,] b, double[,] c)
{
    for (int i = 0; i < n; i++)
        for (int j = 0; j < n; j++)
            for (int k = 0; k < n; k++)
                C[i, j] += a[i, k] * b[k, j]; 
}

基本上，在大步快速地浏览内存会对性能产生不利影响。在B[k, j]中访问k的访问模式正是如此。因此，我们可以重新排列操作，使得最内层循环仅在矩阵的第二个访问索引上进行操作。

void multiply(double[,] a, double[,] B, double[,] c)
{  
   for (i = 0; i < n; i++)
   {  
      double t = a[i, 0];
      for (int j = 0; j < n; j++)
         c[i, j] = t * b[0, j];

      for (int k = 1; k < n; k++)
      {
         double s = 0;
         for (int j = 0; j < n; j++ )
            s += a[i, k] * b[k, j];
         c[i, j] = s;
      }
   }
}

这是该页面上给出的示例。但另一种选择是事先将B [k，*]的内容复制到数组中，并在内部循环计算中使用此数组。即使涉及数据复制，这种方法通常比其他方法快得多。即使这可能看起来违反直觉，请随意尝试。

void multiply(double[,] a, double[,] b, double[,] c)
{
    double[] Bcolj = new double[n];
    for (int j = 0; j < n; j++)
    {
        for (int k = 0; k < n; k++)
            Bcolj[k] = b[k, j];

        for (int i = 0; i < n; i++)
        {
            double s = 0;
            for (int k = 0; k < n; k++)
                s += a[i,k] * Bcolj[k];
            c[j, i] = s;
        }
   }
}