对于默认的行主维(即 matrix[i] 表示一行)形式的二维数组,下面两段代码
for(int i = 0; i < Size; ++i)
for(int j = 0; j < Size; ++j)
for(int k = 0; k < Size; ++k)
for(int l = 0; l < Size; ++l)
matS[i][j] += (vecs[i][k] * vecs[j][l] * ovlp[k][l]);
以及
for(int i = 0; i < Size; ++i)
for(int j = 0; j < Size; ++j)
for(int k = 0; k < Size; ++k)
for(int l = 0; l < Size; ++l)
matS[i][j] += (vecs[k][i] * vecs[l][j] * ovlp[k][l]);
看似区别很小,而且如果对第一段代码中的 vecs 矩阵在执行此四重循环之前先做转置,则两段代码的结果也是相同的。不同的是,第一段代码的执行效率比第二段高出许多(以 1240*1240 大小的矩阵为例,每次对一个最外层变量 i 执行完里层的三重循环,第一段代码大约需要 6 秒,而第二段需要近 18 秒)。
导致此问题的原因是,四个指标变化速度从快到慢依次是 l, k, j 和 i,代码 1 中始终将变化更频繁的指标放到后面(即 “[i][k]” 这样),而代码 2 中则有变化更频繁的指标在变化较慢的指标之前(即“[k][i]”这样)。CPU 执行时,会将连续数列中的相邻元素一起放入寄存器,使用矢量化指令集优化时尤其如此,故代码 1 的形式更有利于 CPU 快速读取数据,代码 2 则需要更多的内存跳转。