为什么矩阵乘法算法中循环的顺序会影响性能？
为什么矩阵乘法算法中循环的顺序会影响性能？
我有两个用于计算两个矩阵乘积的函数：\n
void MultiplyMatrices_1(int **a, int **b, int **c, int n){
for (int i = 0; i < n; i++)
for (int j = 0; j < n; j++)
for (int k = 0; k < n; k++)
c[i][j] = c[i][j] + a[i][k]*b[k][j];
}
void MultiplyMatrices_2(int **a, int **b, int **c, int n){
for (int i = 0; i < n; i++)
for (int k = 0; k < n; k++)
for (int j = 0; j < n; j++)
c[i][j] = c[i][j] + a[i][k]*b[k][j];
}
\n我使用`gprof`运行并分析了两个可执行文件，除了这两个函数外，其他代码完全相同。在大小为2048 x 2048的矩阵上，第二个函数的执行速度明显（大约快了5倍）。有任何想法为什么会出现这种情况？
为什么矩阵乘法算法中循环的顺序会影响性能？
除了内存局部性之外，编译器优化也是一个关键因素。对于向量和矩阵运算来说，一个关键的优化是循环展开。
在这个内部循环中，可以看到i和j都没有变化。这意味着它可以被重新编写为：
for (int k = 0; k < n; k+=4) {
int * aik = &a[i][k];
c[i][j] +=
+ aik[0]*b[k][j]
+ aik[1]*b[k+1][j]
+ aik[2]*b[k+2][j]
+ aik[3]*b[k+3][j];
}
可以看到，循环次数减少了四倍，访问c[i][j]的次数也减少了四倍。而且，a[i][k]在内存中连续访问。内存访问和乘法可以在CPU中并行进行。
但是，如果n不是4、6或8的倍数（或者编译器决定将其展开为其他倍数），编译器会自动处理这个问题。
为了加速这个解决方案，你可以尝试先转置矩阵b。这需要一些额外的工作和编码，但这样对b的访问也是连续的（因为你交换了[k]和[j]）。
另一个可以提高性能的方法是使用多线程进行乘法运算。在4核CPU上，这可以提高性能约3倍。
最后，你还可以考虑使用float或double类型。你可能会认为int类型会更快，但事实并非总是如此，因为浮点运算可以在硬件和编译器中进行更重的优化。
第二个示例中，c[i][j]在每次迭代中都在变化，这使得优化变得更加困难。
这个问题的出现原因可能是内存局部性。当重新排列循环时，内部循环所需的内存更接近，可以被缓存，而在效率低下的版本中，需要从整个数据集中访问内存。
测试这个假设的方法是在两段代码上运行缓存调试器（如
解决这个问题的方法是重新排列循环，使得内部循环所需的内存更接近，可以被缓存。通过这种方式，可以减少缓存未命中，从而提高性能。
以下是一个示例代码，演示了如何重新排列循环以提高内存局部性和性能：
#include
const int N = 1000; // matrix size
int main() {
int A[N][N];
int B[N][N];
int C[N][N];
// Original loop order
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];
}
}
}
// Reordered loop for better memory locality
for (int i = 0; i < N; i++) {
for (int k = 0; k < N; k++) {
for (int j = 0; j < N; j++) {
C[i][j] += A[i][k] * B[k][j];
}
}
}
return 0;
}
通过重新排列循环，将内部循环所需的内存访问集中在一起，可以减少缓存未命中并提高性能。在实际应用中，可以通过使用缓存友好的算法和数据结构，以及优化内存访问模式来进一步改善性能。
从这段内容中，我们可以得出以下结论：
问题的原因：
- 计算机的内存层次结构中，不同类型的内存具有不同的性能特征。寄存器是最快的内存，但容量有限；主内存容量大，但访问速度较慢。为了提高性能，计算机通常在处理器和主内存之间添加多级缓存，这些缓存速度介于寄存器和主内存之间。如果从缓存中读取数据，速度会比从主内存中读取数据快很多。因此，如果能够保持对缓存的连续访问，性能将会更好。
- 大多数程序的内存访问模式是，如果从内存中读取一个字节，程序往往会读取附近的多个值。为了避免每次都需要从主内存中读取这些附近的值，缓存通常是以块的形式进行读取。当读取一个值时，缓存会将该值周围的一块内存值（通常为1KB至1MB）加载到缓存中。这样，如果程序再次读取附近的值，这些值已经在缓存中，就不需要从主内存中读取。
解决方法：
- 将内层循环的顺序调整为j循环内，可以最大程度上利用局部性原理，减少缓存未命中的情况，从而提高性能。
- 可以进一步使用更复杂的循环优化技术，如分块（blocking）技术，将数组划分为可以在缓存中保持更长时间的子区域，然后在这些块上进行多个操作，以计算整体结果。
文章整理如下：
我认为你所看到的是计算机内存层次结构中局部性原理的影响。通常，计算机内存被分为不同的类型，具有不同的性能特征（这通常被称为内存层次结构）。最快的内存位于处理器的寄存器中，可以在一个时钟周期内访问和读取。然而，通常只有少数几个这样的寄存器（通常不超过1KB）。计算机的主内存则非常大（例如8GB），但访问速度较慢。为了提高性能，计算机通常在处理器和主内存之间物理上构造了多级缓存。这些缓存比寄存器慢但比主内存快得多，因此如果在缓存中查找某个值，速度要比访问主内存快得多（通常快5-25倍）。在访问内存时，处理器首先检查内存缓存中的值，然后再返回主内存读取该值。如果您始终在缓存中访问值，性能将比随机访问内存的情况要好得多。
大多数程序的内存访问模式是，如果从内存中读取一个字节，程序往往会从该内存区域附近读取多个不同的值。因此，这些缓存通常被设计成，当您从内存中读取一个值时，会将该值周围一块内存（通常在1KB至1MB之间）的值也加载到缓存中。这样，如果程序再次读取附近的值，这些值已经在缓存中，就不需要从主内存中读取。
现在，最后一个细节 - 在C/C++中，数组是按行主序（row-major order）存储的，这意味着矩阵的每一行的值都存储在一起。因此，在内存中，数组的存储方式是首先存储第一行，然后是第二行，然后是第三行，以此类推。
有了这些信息，我们来看看你的代码。第一个版本如下：
for (int i = 0; i < n; i++)
for (int j = 0; j < n; j++)
for (int k = 0; k < n; k++)
c[i][j] = c[i][j] + a[i][k]*b[k][j];
现在，让我们看看内层循环的那一行代码。在每次迭代中，k的值在增加。这意味着当运行内层循环时，每次迭代很可能在加载b[k][j]的值时发生缓存未命中。原因是，由于矩阵是按行主序存储的，每次增加k时，都会跳过整个矩阵的一整行，并且跳过缓存中已经加载的值。然而，查找c[i][j]（因为i和j是相同的）时不会发生未命中，a[i][k]的未命中可能性也较低，因为这些值是按行主序存储的，如果在上一次迭代中已经将a[i][k]缓存，那么这次迭代读取的a[i][k]值从相邻的内存位置读取。因此，内层循环的每次迭代很可能会发生一次缓存未命中。
但是考虑下面这个第二个版本的代码：
for (int i = 0; i < n; i++)
for (int k = 0; k < n; k++)
for (int j = 0; j < n; j++)
c[i][j] = c[i][j] + a[i][k]*b[k][j];
现在，由于每次迭代都增加了j，我们来思考一下在内层循环的那一行代码中，你可能会有多少次缓存未命中。由于值是按行主序存储的，c[i][j]的值很可能在缓存中，因为上一次迭代中的c[i][j]值很可能也被缓存并准备好被读取。同样，b[k][j]可能也被缓存了，而且由于i和k没有改变，a[i][k]也很可能被缓存。这意味着在内层循环的每次迭代中，你很可能没有缓存未命中。
这意味着第二个版本的代码在循环的每次迭代中都不太可能发生缓存未命中，而第一个版本几乎肯定会发生缓存未命中。因此，第二个循环比第一个循环更快，正如你所看到的。
有趣的是，许多编译器开始原型支持检测到第二个版本的代码比第一个版本更快。有些编译器会尝试自动重写代码以最大程度地发挥并行性。如果你有一本《紫龙书》，第11章讨论了这些编译器的工作原理。
此外，您还可以使用更复杂的循环优化技术来进一步优化此循环的性能。例如，可以使用分块（blocking）技术，将数组划分为可以在缓存中保持更长时间的子区域，然后在这些块上进行多个操作，以计算整体结果。
希望这可以帮助你！