特征:如何加速a += coeffs * coeffs.转置()
Eigen: how to speed up a += coeffs * coeffs.transpose()
我需要计算许多(大约400k)小线性最小二乘问题的解。每个问题包含10-300个方程,只有7个变量。为了解决这些问题,我使用了特征库。直接求解太费时间了,我把每个问题都转换成求解7x7的线性方程组通过手工求导。
我收到了很好的加速,但我想再次提高性能。
我使用vagrind来分析我的程序,我发现自代价最高的操作是特征矩阵的算子+=。此操作需要调用a.l ldlt().solve(b);
我用这个算子来组合每个方程组的A矩阵和B向量
//I cal these code to solve each problem
const int nVars = 7;
//i really need double precision
Eigen::Matrix<double, nVars, nVars> a = Eigen::Matrix<double, nVars, nVars>::Zero();
Eigen::Matrix<double, nVars, 1> b = Eigen::Matrix<double, nVars, 1>::Zero();
Eigen::Matrix<double, nVars, 1> equationCoeffs;
//............................
//Somewhere in big cycle.
//equationCoeffs and z are updated on each iteration
a += equationCoeffs * equationCoeffs.transpose();
b += equationCoeffs * z;
其中z是某个标量
所以我的问题是:我怎样才能提高这些操作的性能?
PS抱歉我的英语很差
你可以尝试一次分配一个足够大的矩阵(例如300 x 7)来存储所有系数,然后让Eigen优化的矩阵-矩阵乘积核为你完成这项工作,而不是手工形成普通方程的矩阵和向量分量:
Matrix<double,Dynamic,nbVars> D(300,nbVars);
VectorXd f(300);
for(...)
{
int nb_equations = ...;
for(i=0..nb_equations-1)
{
D.row(i) = equationCoeffs;
f(i) = z;
}
a = D.topRows(nb_equations).transpose() * D.topRows(nb_equations);
b = D.topRows(nb_equations).transpose() * f.head(nb_equations);
// solve ax=b
}
您可以对矩阵D
同时使用列为主存储和行为主存储,看看哪一个是最好的。
另一种可能的方法是将a
, equationCoeffs
和b
声明为8x8或8x1矩阵或向量,以确保equationCoeffs(7)==0
。这样可以最大化SIMD的使用。然后在调用LDLT时使用a.topLeftCorners<7,7>()
, b.head<7>()
。您甚至可以将此策略与前一个策略结合使用。
最后,如果您的CPU支持AVX或FMA,您可以使用devel分支并使用-mavx
或-mfma
进行编译,以获得显着的加速
如果您可以使用g++5.1,您可能会想看看OpenMP(http://openmp.org/mp-documents/OpenMP4.0.0.Examples.pdf)。g++ 5.1(或C的gcc5.1)也有一些对OpenACC的基本支持,您也可以尝试一下。未来应该会有更多的OpenACC实现。
另外,如果你有访问英特尔编译器(icc, icpc),它加快了我的代码,甚至只是通过使用它。
如果你可以使用nvidia的nvcc,你可以使用thrust库(http://docs.nvidia.com/cuda/thrust/#axzz3g8xJPGHe),他们的github上也有很多示例代码(https://github.com/thrust/thrust)。然而,使用推力并不是那么简单,需要一些真正的思考。
编辑:Thrust也需要Nvidia的GPU。对于AMD卡,我相信有一个名为ArrayFire的库,它看起来非常类似于Thrust(我还没有尝试过那个)
我有一个单一的问题Ax=b与480k浮动变量。矩阵A是稀疏的,用特征值BiCGSTAB求解需要4.8秒。
我以前也用过ViennaCL,所以我试着在那里解决同样的问题,只花了1.2秒。扩展的增加实现了
- 转置矩阵:交换元素不会更改值
- 使用 Eigen 3 库编写一个带有转置作为参数的函数
- 犰狳(C++)中的快速阵列置换(广义张量转置)
- 为什么转置这个 std::vector<std::vector<std::string> > 这么慢?
- 转置结构容器
- 我的转置矩阵代码有什么问题?
- 数组的转置和乘法
- 在C++中使用矢量转置 2D 矩阵
- 特征:块转置
- MKL矩形矩阵Inplace转置:不使用多个核心
- 并行转置不同的矩阵
- 关于次级对角线的转置(翻转)矩阵
- 输出是从您输入的矩阵中打印出矩阵的转置,但我的代码只是打印出您输入的第一个矩阵
- 如何在阵列火中避免翻转和转置的memcpy?
- CUDA矩阵与共享内存转置
- C 阵列的复合物共轭转置
- 转置期间的动态内存分配
- 转置的一维矢量的平均矢量
- 在CUDA中具有共享MEM的非方面矩阵转置
- 特征:如何加速a += coeffs * coeffs.转置()