用Eigen和c++做了一个大量矩阵乘积的列
Using Eigen and C++ to do a colsum of massive matrix product
我试图计算列(N * p),其中N是一个稀疏的,1M × 2500矩阵,p是一个密集的2500 × 1.5M矩阵。我正在使用Eigen c++库和英特尔的MKL库。问题是矩阵N*P实际上不能存在于内存中,它太大了(~ 10tb)。我的问题是,Eigen是否能够通过惰性求值和并行性的某种组合来处理这种计算?这里说Eigen不会不必要地生成临时矩阵:http://eigen.tuxfamily.org/dox-devel/TopicLazyEvaluation.html
但是Eigen知道如何在内存中计算N * p吗?IE:它必须做类似于colsum(N * P_1) ++ colsum(N * P_2) ++…++ colsum(N * P_n),其中P按列分成N个不同的子矩阵,"++"是串联。
我正在使用128 GB RAM。
我尝试了一下,但最终出现了一个糟糕的malloc(我只在Win8上运行8GB)。我设置了我的main()
,并使用了我编写的非内联colsum
函数。
int main(int argc, char *argv[])
{
Eigen::MatrixXd dense = Eigen::MatrixXd::Random(1000, 100000);
Eigen::SparseMatrix<double> sparse(100000, 1000);
typedef Triplet<int> Trip;
std::vector<Trip> trps(dense.rows());
for(int i = 0; i < dense.rows(); i++)
{
trps[i] = Trip(20*i, i, 2);
}
sparse.setFromTriplets(trps.begin(), trps.end());
VectorXd res = colsum(sparse, dense);
std::cout << res;
std::cin >> argc;
return 0;
}
这个尝试很简单:
__declspec(noinline) VectorXd
colsum(const Eigen::SparseMatrix<double> &sparse, const Eigen::MatrixXd &dense)
{
return (sparse * dense).colwise().sum();
}
有一个坏的malloc。所以看起来你必须自己手动拆分它(除非其他人有更好的解决方案)。
编辑
我改进了函数一点,但得到相同的坏malloc:
__declspec(noinline) VectorXd
colsum(const Eigen::SparseMatrix<double> &sparse, const Eigen::MatrixXd &dense)
{
return (sparse * dense).topRows(4).colwise().sum();
}
编辑2
另一种选择是使稀疏矩阵密集并强制执行惰性求值。我不认为它可以用稀疏矩阵(哦,好吧)。
__declspec(noinline) VectorXd
colsum(const Eigen::SparseMatrix<double> &sparse, const Eigen::MatrixXd &dense)
{
Eigen::MatrixXd denseSparse(sparse);
return denseSparse.lazyProduct(dense).colwise().sum();
}
这并没有给我坏的malloc,但计算了很多毫无意义的0*x_i
表达式。
回答您的问题:特别是当涉及到产品时,Eigen经常将部分表达式求值为临时值。在某些情况下,这可以优化但尚未实现,在某些情况下,这实际上是实现它的最有效方式。
然而,在您的情况下,您可以简单地计算N
的列(1 x 2500向量)并将其乘以P
。也许未来版本的Eigen将能够自己进行这种优化,但大多数情况下,在让计算机完成其余工作之前,自己进行特定问题的优化是一个好主意。
Btw:恐怕sparse.colwise()
还没有实现,所以你必须手动计算。如果你很懒,你可以代替计算Eigen::RowVectorXd Nsum = Eigen::RowVectorXd::Ones(N.rows())*P;
(我没有检查它,但这实际上可能得到优化接近最优的代码,与最新版本的Eigen)。
- 如何创建一个CMake变量,除非显式重写,否则使用默认值
- 删除一个线程上有数百万个字符串的大型哈希映射会影响另一个线程的性能
- 为什么两个不同的未命名名称空间可以共存于一个cpp文件中
- 运行同一解决方案的另一个项目的项目
- 挂起和取消挂起一个文件DLL
- 用C++中的一个变量定义一个常量
- 函数向量_指针有不同的原型,我可以构建一个吗
- 在c++中用vector填充一个简单的动态数组
- 如何在选项卡视图Qt中设置一个新项目,并保存以前的项目
- 预处理器:插入结构名称中的前一个行号
- 我在c++代码中生成了一个运行时#3异常
- 我想将一个对T类型的非常量左值引用绑定到一个T类型的临时值
- 从链接列表c++中删除一个项目
- 告诉一个 const char 数组,除了编译时 C 样式的字符串外,它不以 '