algorithm - 缓存高效的矩阵转置程序？

Question

所以转置矩阵的明显方法是使用:

  for( int i = 0; i < n; i++ )

    for( int j = 0; j < n; j++ )

      destination[j+i*n] = source[i+j*n];

但我想要一些可以利用位置和缓存阻塞的东西。我正在查找它，但找不到可以执行此操作的代码，但有人告诉我这应该是对原始代码的非常简单的修改。有任何想法吗？

编辑:我有一个 2000x2000 矩阵，我想知道如何使用两个 for 循环更改代码，基本上将矩阵拆分为我单独转置的 block ，比如 2x2 block ，或者40x40 block ，看看哪个 block 大小最有效。

Edit2:矩阵按列主顺序存储，也就是说对于矩阵

a1 a2    
a3 a4

存储为a1 a3 a2 a4。

Answer 1

您可能需要四个循环 - 两个循环遍历 block ，然后另外两个循环执行单个 block 的转置复制。为简单起见，假设一个 block 大小除以矩阵的大小，我想是这样的，尽管我想在信封背面画一些图片以确保:

for (int i = 0; i < n; i += blocksize) {
    for (int j = 0; j < n; j += blocksize) {
        // transpose the block beginning at [i,j]
        for (int k = i; k < i + blocksize; ++k) {
            for (int l = j; l < j + blocksize; ++l) {
                dst[k + l*n] = src[l + k*n];
            }
        }
    }
}

一个重要的进一步见解是，实际上有一个缓存不经意的算法(参见 http://en.wikipedia.org/wiki/Cache-oblivious_algorithm ，它使用这个确切的问题作为示例)。 “无视缓存”的非正式定义是，您无需尝试调整任何参数(在本例中为 block 大小)即可达到良好/最佳缓存性能。这种情况下的解决方案是通过递归地将矩阵分成两半来转置，并将两半转置到目标中的正确位置。

无论缓存大小实际是多少，此递归都会利用它。我希望与您的策略相比会有一些额外的管理开销，您的策略是使用性能实验实际上直接跳到缓存真正启动的递归点，并且不再继续。另一方面，您的性能实验可能会给您一个适用于您的机器但不适用于客户机器的答案。