稀疏矩阵的填充速度非常慢,Eigen 中没有内存增益

Very slow filling of sparse matrix and no memory gain in Eigen

我正在尝试将 Eigen tutorial 中给出的示例实现为伪代码。据我了解,它说明了填充稀疏矩阵的推荐方法,前提是每列的非零条目数已知。

伪代码在Header“Filling a sparse Matrix”下找到,写成如下:

1: SparseMatrix<double> mat(rows,cols);         // default is column major
2: mat.reserve(VectorXi::Constant(cols,6));
3: for each i,j such that v_ij != 0
4:   mat.insert(i,j) = v_ij;                    // alternative: mat.coeffRef(i,j) += v_ij;
5: mat.makeCompressed();                        // optional

下面是我尝试将其转换为 C 的尝试。我已经(希望)这样写 vee() ,它会在每列中创建 2500 个非零元素。因此 2500 应对应示例中的 6。我也将它设置为 3000 来测试 make.Compressed。

不幸的是,我不了解该程序的行为。它在几秒钟内执行 i=0...3000,然后卡住了几分钟。然后它达到 6000 并再次卡住几分钟。为什么会这样以及如何获得更好的性能?

另外,内存占用也很奇怪。您可以看到,有时在接近尾声时,Eigen 使用的内存明显多于 GSL.The 中相应密集矩阵所需的内存,使用的内存也波动很大。步长大于 100MB

我这样编译 运行:

ludi@ludi-M17xR4:~/Desktop/tests$ g++ -o eigenfill.x eigenfill.cc -L/usr/local/lib -lgsl -lgslcblas && ./eigenfill.x

#include <iostream>
#include <stdio.h>
#include <stdlib.h>     
#include <Eigen/Sparse>
#include <gsl/gsl_matrix.h>
#define rows 1e4
#define cols 1e4



/*-- DECLARATIONES --*/
int FillMatrix(Eigen::SparseMatrix<double> mat);
double vee(int i, int j);

int main()
{
printf("---> Watch gsl matrix memory usage!\n");    
gsl_matrix *testmat = gsl_matrix_calloc(rows, cols);
sleep(20);
gsl_matrix_free(testmat);
printf("---> Watch eigen matrix memory usage!\n");  

Eigen::SparseMatrix<double> mat(rows,cols);         // default is column major
FillMatrix(mat);
printf("------------------------DONE");
return(0);
}

/*-- --*/
int FillMatrix(Eigen::SparseMatrix<double> mat)
{
int i, j;
Eigen::VectorXd Vres;
mat.reserve(Eigen::VectorXi::Constant(cols,3000));
    for(i=0;i<rows;i++)
    {
        if(i%500==0){printf("i= %i\n", i);}
    for(j=0;j<cols;j++)
        {
        if (vee(i,j) != 0){mat.insert(i,j) = vee(i,j);    /*alternative: mat.coeffRef(i,j) += v_ij;*/ }
        }

    }
printf("--->starting compression"); 
mat.makeCompressed();

return(0);
}
/*-- --*/
double vee(int i, int j)
{
    double result = 0.0;

    if(j%4 == 0){result =1.0;}

    return result;
}
/*-- --*/

编辑

一个答案提醒我需要使用地址,因为 FillMatrix() 的局部变量在 运行 结束时将消失。我尝试了以下,但无法编译:

#include <iostream>
#include <stdio.h>
#include <stdlib.h>     
#include <Eigen/Sparse>
#include <gsl/gsl_matrix.h>
#define rows 1e4
#define cols 1e4



/*-- DECLARATIONES --*/
int FillMatrix(Eigen::SparseMatrix<double> & mat);
double vee(int i, int j);

int main()
{
printf("---> Watch gsl matrix memory usage!\n");    
gsl_matrix *testmat = gsl_matrix_calloc(rows, cols);
sleep(20);
gsl_matrix_free(testmat);
printf("---> Watch eigen matrix memory usage!\n");  

Eigen::SparseMatrix<double> mat(rows,cols);         // default is column major
FillMatrix(& mat);
printf("------------------------>DONE\n");
return(0);
}

/*-- --*/
int FillMatrix(Eigen::SparseMatrix<double> &mat)
{
int i, j;
Eigen::VectorXd Vres;
mat.reserve(Eigen::VectorXi::Constant(cols,3000));
    for(i=0;i<rows;i++)
    {
        if(i%500==0){printf("i= %i\n", i);}
    for(j=0;j<cols;j++)
        {
        if (vee(i,j) != 0){mat.insert(i,j) = vee(i,j);    /*alternative: mat.coeffRef(i,j) += v_ij;*/ }
        }

    }
printf("--->starting compression\n");   
mat.makeCompressed();

return(0);
}
/*-- --*/
double vee(int i, int j)
{
    double result = 0.0;

    if(i%4 == 0){result =1.0;}

    return result;
}
/*-- --*/

错误信息是:

ludi@ludi-M17xR4:~/Desktop/tests$ g++ -o eigenfill.x eigenfill.cc -L/usr/local/lib -lgsl -lgslcblas && ./eigenfill.x
eigenfill.cc: In function ‘int main()’:
eigenfill.cc:24:17: error: invalid initialization of non-const reference of type ‘Eigen::SparseMatrix<double>&’ from an rvalue of type ‘Eigen::SparseMatrix<double>*’
 FillMatrix(& mat);
                 ^
eigenfill.cc:12:5: error: in passing argument 1 of ‘int FillMatrix(Eigen::SparseMatrix<double>&)’
 int FillMatrix(Eigen::SparseMatrix<double> & mat);
     ^
ludi@ludi-M17xR

编辑 如果我写它编译:

FillMatrix(mat);

而不是

FillMatrix(&mat);

我不明白。最后一个不应该是正确的吗?

这是因为您将四栏中的一栏填满了。看起来像下面这样:

1 0 0 0 1 0 0 0 1 ...
1 0 0 0 1 0 0 0 1 ...
1 0 0 0 1 0 0 0 1 ...
1 0 0 0 1 0 0 0 1 ...
...

j%4 == 0 替换为 i%4 == 0 即可。内存使用量将是:2500*cols*(sizeof(double)+sizeof(int)) + cols*sizeof(int).

您还必须修复 FillMatrix 的原型,以便填充真实矩阵而不是临时矩阵:

int FillMatrix(Eigen::SparseMatrix<double> &mat);

最后,如果您的目标确实是使用具有大约 1/4 非零值的稀疏矩阵,那么密集表示可能会更快。