为什么稀疏密集乘法比密集稀疏乘法更快？

Posted 2021-03-27

我很好奇为什么用密集矩阵乘以稀疏矩阵需要不同于反向矩阵的时间。算法有明显的不同吗？

这是matlab 2018a中的一个例子：

a=sprand(M,M,0.01);
b=rand(M);
tic;ref1=a*b;t_axb=toc
tic;ref2=b*a;t_bxa=toc

以下是使用1个线程的Eigen 3和C ++的示例：

//prepare acol=MxM ColMajor Eigen sparse matrix with 0.01 density
...
Map<Matrix<double,M,M,ColMajor> > bcol (PR, M, M );
double tic,toc;

tic=getHighResolutionTime();
result=acol*bcol;
toc=getHighResolutionTime();
printf("
acol*bcol time: %f seconds", (toc - tic));

tic=getHighResolutionTime();
result=bcol*acol;
toc=getHighResolutionTime();
printf("
bcol*acol time: %f seconds
", (toc - tic));

当M = 4000时，结果如下：

t_axb =
    0.6877
t_bxa =
    0.4803

acol*bcol time: 0.937590 seconds
bcol*acol time: 0.532622 seconds

当M = 10000时，结果是

t_axb =
   11.5649
t_bxa =
    9.7872

acol*bcol time: 20.140380 seconds
bcol*acol time: 8.061626 seconds

在这两种情况下，对于Matlab和Eigen，稀疏密集产品比密集稀疏产品慢。我很好奇

1. 为什么会这样？稀疏密度算法与密集稀疏算法有明显区别吗？ FLOP的数量是一样的，对吧？
2. 为什么特征匹配或超过Matlab的密集稀疏性能而不是稀疏密集产品？性能上的微小差异是正常的，但考虑到两者都是高度优化的库，因此~1.4-1.8的差异似乎很奇怪。我根据文档编译了所有优化的eigen。即-fPIC -fomit-frame-pointer -O3 -DNDEBUG -fopenmp -march=native

答案

您可以通过比较稀疏矩阵时间向量乘积的列主要与行主要存储来观察相同的差异：y = A * x。如果A是行主要（相当于y的每个系数），那么A的每一行可以并行处理而没有任何开销（没有通信，没有额外的临时，没有额外的操作）。相反，如果A是列主要的多线程不能免费提供，并且在大多数情况下，开销大于增益。

即使没有多线程，您也会发现内存访问模式非常不同：

• 行主要：对x的多个随机只读访问，y的每个系数只写一个。
• 专业：x的每个系数被读取一次，但我们得到多个随机读写访问目的地y。

因此，即使没有多线程，情况自然也有利于行专业。