在医学图像重建实施中改善局部性并减少缓存污染

Question

我正在为我的大学做一项与医学用图像重建算法相关的研究.

我陷入了长达3周的困境,我需要提高以下代码的性能:

for (lor=lor0[mypid]; lor <= lor1[mypid]; lor++)
{
  LOR_X = P.symmLOR[lor].x;
  LOR_Y = P.symmLOR[lor].y;
  LOR_XY = P.symmLOR[lor].xy;
  lor_z = P.symmLOR[lor].z;
  LOR_Z_X = P.symmLOR[lor_z].x;
  LOR_Z_Y = P.symmLOR[lor_z].y;
  LOR_Z_XY = P.symmLOR[lor_z].xy;  

  s0 = P.a2r[lor];
  s1 = P.a2r[lor+1];

  for (s=s0; s < s1; s++)
  {
    pixel     = P.a2b[s];
    v         = P.a2p[s]; 

    b[lor]    += v * x[pixel];

    p          = P.symm_Xpixel[pixel];
    b[LOR_X]  += v * x[p];

    p          = P.symm_Ypixel[pixel];
    b[LOR_Y]  += v * x[p];

    p          = P.symm_XYpixel[pixel];
    b[LOR_XY] += v * x[p];


    // do Z symmetry.
    pixel_z    = P.symm_Zpixel[pixel];
    b[lor_z]  += v * x[pixel_z];


    p          = P.symm_Xpixel[pixel_z];
    b[LOR_Z_X]  += v * x[p];


    p          = P.symm_Ypixel[pixel_z];
    b[LOR_Z_Y]  += v * x[p];

    p          = P.symm_XYpixel[pixel_z];
    b[LOR_Z_XY] += v * x[p];

   }

}

对于任何想知道的人来说,代码实现了MLEM转发函数,所有变量都是FLOAT.

经过几次测试后,我注意到代码的这一部分出现了很大的延迟.(你知道,90 - 10规则).

后来,我使用Papi(http://cl.cs.utk.edu/papi/)来测量L1D缓存未命中.正如我所想的那样,Papi确认由于更多的未命中而导致性能下降,特别是对于随机访问b矢量(大小很大).

在互联网上阅读信息我只知道到目前为止提高性能的两个选项:改善数据局部性或减少数据污染.

为了做第一个改进,我将尝试将代码更改为高速缓存感知,就像Ulrich Drepper关于每个程序员应该知道的内存(www.akkadia.org/drepper/cpumemory.pdf)A所提出的那样. 1矩阵乘法.

我相信阻塞SpMV(稀疏矩阵向量乘法)将改善性能.

另一方面,每次程序试图访问b vector时,我们都会遇到所谓的缓存污染.

有没有办法在不使用Cache的情况下使用SIMD指令从b向量加载值？

此外,可以使用像void _mm_stream_ps(float*p,__ m128 a)这样的函数在向量b上存储一个浮点值而不会污染缓存？

我不能使用_mm_stream_ps因为总是存储4个浮点数,但是对b矢量的访问显然是随机的.

我希望在我的困境中明白.

更多信息:v是具有CRS格式的Sparse Matrix商店的列值.我意识到,如果我尝试将CRS格式更改为其他格式,可以进行其他优化,但是,正如我之前所说,我已经进行了几个月的测试,我知道性能下降与向量b上的随机访问有关.从400.000.000 L1D未命中当我不存储在矢量b时,我可以转到100~未命中.

谢谢.

Answer 1

我想说，首先尝试为你的编译器提供一些帮助。

• const像循环之前一样声明外循环的边界。
• 将所有可能的变量（所有LOR_..）声明为局部变量，例如： float LOR_X = P.symmLOR[lor].x;或 size_t s0 = P.a2r[lor];
• 如果您碰巧拥有符合 C99 标准的现代编译器，这尤其适用于循环变量：for (size_t s=s0; s < s1; s++)
• 单独加载和存储向量b 。您在那里访问的项目的位置不依赖于s。因此，创建一个局部变量来保存在内循环之前处理的所有不同情况的实际值 ，在内循环内更新这些局部变量，并在内循环之后存储结果。
• 也许将你的内部循环分成几个。索引计算相对便宜，然后您的系统可能会更好地识别对某些向量的流式访问。
• 查看编译器生成的汇编程序并识别内部循环的代码。尝试一下将尽可能多的“远”负载和存储移出该循环。

编辑：在重新阅读 gravitron 的答案和您的评论之后，这里重要的事情实际上是尽可能将变量声明为本地变量，并检查编译器是否成功将缓存缺失加载并存储在内循环之外。