[eigen] Initializing a sparse matrix - surprisingly slow

[Helmut Jarausch <jarausch@xxxxxxxxxxxxxxxxxxx>]

Hi,

probably I haven't understood the docs:

#include <Eigen/Sparse>
using Eigen::SparseMatrix;
typedef SparseMatrix<double>  SpMatrix;

#define Slower_by_a_factor_of_1000

SpMatrix Laplace( int m ) {  // 2D Laplace : 5 elements per row/column  
at most
  const int  n = m*m,  n0= n-1, n_m = n-m,
             NumNonZeros= 5*n;  // number of non-zeros
  SpMatrix A(n,n);

#ifdef Slower_by_a_factor_of_1000
  A.reserve(NumNonZeros);  // 1000 times slower for n=9801 <<<<< WHY  
???
#else
  A.reserve(Eigen::VectorXi::Constant(n,5));
#endif

  double d = 4,  b = -1, c = -1;

  for (int j= 0; j < n; j++) {
    if ( j >= m ) A.insert(j-m,j)= c;
    if ( j > 0 &&  j%m != 0 ) A.insert(j-1,j)= b;
    A.insert(j,j)= d;
    if ( j < n0 && (j+1)%m != 0 ) A.insert(j+1,j)= b;
    if ( j < n_m ) A.insert(j+m,j)= c;
  }
  return A;
}

This is with the GIT-version of Eigen.

Many thanks for an explanation,
Helmut.