解释来自Agner Fog在c++中优化软件，它减少了数据如何访问和存储在缓存中。

有关术语和详细信息，请参阅维基上关于缓存的条目，我将在这里缩小范围。

缓存被组织成集和行。一次只使用一个集合，其中包含的任何行都可以使用。一行能镜像的内存乘以行数就得到了缓存大小。

对于一个特定的内存地址，我们可以用下面的公式计算出哪个集合应该镜像它:

set = ( address / lineSize ) % numberOfsets

理想情况下，这种公式在集合上给出了一个均匀的分布，因为每个内存地址都有可能被读取(我说的是理想情况)。

很明显，重叠是可能发生的。在缓存丢失的情况下，在缓存中读取内存并替换旧的值。请记住，每个集合都有一些行，其中最近最少使用的行将被新读取的内存覆盖。

我将试着遵循Agner的例子:

Assume each set has 4 lines, each holding 64 bytes. We first attempt to read the address 0x2710, which goes in set 28. And then we also attempt to read addresses 0x2F00, 0x3700, 0x3F00 and 0x4700. All of these belong to the same set. Before reading 0x4700, all lines in the set would have been occupied. Reading that memory evicts an existing line in the set, the line that initially was holding 0x2710. The problem lies in the fact that we read addresses that are (for this example) 0x800 apart. This is the critical stride (again, for this example).

关键步幅也可以计算:

criticalStride = numberOfSets * lineSize

间隔为criticalStride或一个倍数的变量争夺相同的缓存线。

这是理论部分。接下来是解释(也是Agner，为了避免出错，我一直在仔细阅读):

假设矩阵为64x64(记住，效果根据缓存不同)，缓存为8kb，每组4行*行大小为64字节。每行可以容纳矩阵中的8个元素(64位整型)。

关键步幅为2048字节，对应矩阵的4行(在内存中是连续的)。

假设我们正在处理第28行。我们试图用这一行的元素和第28列的元素交换。一行的前8个元素组成了一条缓存线，但它们将在第28列中进入8个不同的缓存线。记住，关键步幅是间隔4行(一列中的4个连续元素)。

当列中到达元素16时(每集4行，4行间隔=麻烦)，ex-0元素将从缓存中被清除。当我们到达列的末尾时，之前的所有缓存行将丢失，需要在访问下一个元素时重新加载(整行将被覆盖)。

如果大小不是关键步幅的倍数，就会破坏这种完美的灾难场景，因为我们不再处理垂直方向上关键步幅的元素，因此缓存重新加载的数量会严重减少。

另一个免责声明——我只是理解了这个解释，希望我理解了，但我可能错了。不管怎样，我在等待神秘主义的回应(或确认)。：）

Luchian解释了为什么会发生这种行为，但我认为展示这个问题的一种可能的解决方案，同时展示一些缓参无关算法是一个好主意。

你的算法基本上是:

for (int i = 0; i < N; i++) 
   for (int j = 0; j < N; j++) 
        A[j][i] = A[i][j];

这对于现代CPU来说简直太可怕了。一种解决方案是了解您的缓存系统的细节，并调整算法以避免这些问题。只要你知道这些细节，工作就很好。不是特别便携。

我们能做得更好吗?是的，我们可以:解决这个问题的一般方法是缓存无关算法，顾名思义，避免依赖于特定的缓存大小[1]

解决方案是这样的:

void recursiveTranspose(int i0, int i1, int j0, int j1) {
    int di = i1 - i0, dj = j1 - j0;
    const int LEAFSIZE = 32; // well ok caching still affects this one here
    if (di >= dj && di > LEAFSIZE) {
        int im = (i0 + i1) / 2;
        recursiveTranspose(i0, im, j0, j1);
        recursiveTranspose(im, i1, j0, j1);
    } else if (dj > LEAFSIZE) {
        int jm = (j0 + j1) / 2;
        recursiveTranspose(i0, i1, j0, jm);
        recursiveTranspose(i0, i1, jm, j1);
    } else {
    for (int i = i0; i < i1; i++ )
        for (int j = j0; j < j1; j++ )
            mat[j][i] = mat[i][j];
    }
}

稍微复杂一点，但是一个简短的测试显示了一些非常有趣的东西，在我的古老的e8400与VS2010 x64版本，测试代码为MATSIZE 8192

int main() {
    LARGE_INTEGER start, end, freq;
    QueryPerformanceFrequency(&freq);
    QueryPerformanceCounter(&start);
    recursiveTranspose(0, MATSIZE, 0, MATSIZE);
    QueryPerformanceCounter(&end);
    printf("recursive: %.2fms\n", (end.QuadPart - start.QuadPart) / (double(freq.QuadPart) / 1000));

    QueryPerformanceCounter(&start);
    transpose();
    QueryPerformanceCounter(&end);
    printf("iterative: %.2fms\n", (end.QuadPart - start.QuadPart) / (double(freq.QuadPart) / 1000));
    return 0;
}

results: 
recursive: 480.58ms
iterative: 3678.46ms

Edit: About the influence of size: It is much less pronounced although still noticeable to some degree, that's because we're using the iterative solution as a leaf node instead of recursing down to 1 (the usual optimization for recursive algorithms). If we set LEAFSIZE = 1, the cache has no influence for me [8193: 1214.06; 8192: 1171.62ms, 8191: 1351.07ms - that's inside the margin of error, the fluctuations are in the 100ms area; this "benchmark" isn't something that I'd be too comfortable with if we wanted completely accurate values])

[1]这些东西的来源:好吧，如果你不能从与Leiserson和co一起工作的人那里得到一个讲座。我认为他们的论文是一个很好的起点。这些算法仍然很少被描述——CLR只有一个脚注。但这仍然是给人们惊喜的好方法。

编辑(注:我不是发布这个答案的人;我只是想补充一下): 以下是上述代码的完整c++版本:

template<class InIt, class OutIt>
void transpose(InIt const input, OutIt const output,
    size_t const rows, size_t const columns,
    size_t const r1 = 0, size_t const c1 = 0,
    size_t r2 = ~(size_t) 0, size_t c2 = ~(size_t) 0,
    size_t const leaf = 0x20)
{
    if (!~c2) { c2 = columns - c1; }
    if (!~r2) { r2 = rows - r1; }
    size_t const di = r2 - r1, dj = c2 - c1;
    if (di >= dj && di > leaf)
    {
        transpose(input, output, rows, columns, r1, c1, (r1 + r2) / 2, c2);
        transpose(input, output, rows, columns, (r1 + r2) / 2, c1, r2, c2);
    }
    else if (dj > leaf)
    {
        transpose(input, output, rows, columns, r1, c1, r2, (c1 + c2) / 2);
        transpose(input, output, rows, columns, r1, (c1 + c2) / 2, r2, c2);
    }
    else
    {
        for (ptrdiff_t i1 = (ptrdiff_t) r1, i2 = (ptrdiff_t) (i1 * columns);
            i1 < (ptrdiff_t) r2; ++i1, i2 += (ptrdiff_t) columns)
        {
            for (ptrdiff_t j1 = (ptrdiff_t) c1, j2 = (ptrdiff_t) (j1 * rows);
                j1 < (ptrdiff_t) c2; ++j1, j2 += (ptrdiff_t) rows)
            {
                output[j2 + i1] = input[i2 + j1];
            }
        }
    }
}

为了说明Luchian Grigore回答中的解释，下面是64x64和65x65矩阵两种情况下的矩阵缓存存在情况(有关数字的详细信息，请参阅上面的链接)。

下面动画中的颜色意味着:

-不在缓存中， -在缓存中， -缓存命中， -从RAM中读取， - cache Miss。

64x64机箱:

注意，几乎每一次对新行的访问都会导致缓存丢失。现在看看正常情况下，65x65矩阵的情况:

在这里，您可以看到在初始热身之后的大多数访问都是缓存命中。这就是CPU缓存的工作原理。

为上面的动画生成帧的代码可以在这里看到。