mike.dl
mike.dl

Reputation: 67

Problems using Eigen c++ library in Rcpp for Galerkin Matrix

I'm trying working out the following C++ code in RStudio.

// [[Rcpp::depends(RcppEigen)]]
#include <RcppEigen.h>

using namespace Rcpp;

// [[Rcpp::export]]
#include <iostream>
#include <cmath>

using Eigen::Dense;
using Eigen::SparseLU;
using Eigen::Sparse;
using Eigen::SparseMatrix;

using namespace std;

Eigen::SparseMatrix<double> getGalMat(int N) {

  // Note: N hast to bigger or equal to 5!!!
  assert(N >= 5);

  typedef Eigen::SparseMatrix<double>::Index index_t;
  typedef Eigen::Triplet<double> triplet_t;
  std::vector<triplet_t> triplets;

  // reserve "minimal" vector size (the number of non-zero entries)
  triplets.reserve(5*N - 6);

  // N minus 1
  int Nm = N - 1;

  // set the (off-) diagonals
  double diag_m2 = + 16;
  double diag_m1 = - 64;
  double diag    = + 96;
  double diag_p1 = - 64;
  double diag_p2 = + 16;

  // set first and last 2 rows by hand
  // A(1  ,:)
  triplets.push_back({0   ,0   ,diag   }); // A(1,1)
  triplets.push_back({0   ,1   ,diag_p1}); // A(1,2)
  triplets.push_back({0   ,2   ,diag_p2}); // A(1,3)
  // A(2  ,:)
  triplets.push_back({1   ,0   ,diag_m1}); // A(2,1)
  triplets.push_back({1   ,1   ,diag   }); // A(2,2)
  triplets.push_back({1   ,2   ,diag_p1}); // A(2,3)
  triplets.push_back({1   ,3   ,diag_p2}); // A(2,4)
  // A(N-1,:)
  triplets.push_back({Nm-1,Nm-3,diag_m2}); // A(N-1,N-3)
  triplets.push_back({Nm-1,Nm-2,diag_m1}); // A(N-1,N-2)
  triplets.push_back({Nm-1,Nm-1,diag   }); // A(N-1,N-1)
  triplets.push_back({Nm-1,Nm  ,diag_p1}); // A(N-1,N  )
  // A(N  ,:)
  triplets.push_back({Nm  ,Nm-2,diag_m2}); // A(N,N-2)
  triplets.push_back({Nm  ,Nm-1,diag_m1}); // A(N,N-1)
  triplets.push_back({Nm  ,Nm  ,diag   }); // A(N,N  )
  // loop over remaining rows
  for (int i = 2; i < Nm-1; i++) {
    triplets.push_back({i,i-2,diag_m2}); // A(i,i-2)
    triplets.push_back({i,i-1,diag_m1}); // A(i,i-1)
    triplets.push_back({i,i  ,diag   }); // A(i,i  )
    triplets.push_back({i,i+1,diag_p1}); // A(i,i+1)
    triplets.push_back({i,i+2,diag_p2}); // A(i,i+2)
  }

  // let EIGEN build the sparse matrix from our triplets
  Eigen::SparseMatrix<double> spMat(N,N);
  spMat.setFromTriplets(triplets.begin(), triplets.end());

  // return
  return spMat;

}

When running it I have a long list of error from row 41 onwords, corresponding to:

triplets.push_back({0   ,0   ,diag   }); // A(1,1)

getting the error:'extended initializer lists only available with -std=c++0x or -std=gnu++0x [enabled by default]'

Would anyone have suggestions?

Thanks in advance for your help.

Upvotes: 0

Views: 194

Answers (2)

Avi Ginsburg
Avi Ginsburg

Reputation: 10596

As you don't seem to be able to enable C++11, the alternative is to make the code C++98 compatible. Instead of:

triplets.push_back({0   ,0   ,diag   }); // A(1,1)

explicitly use the constructor:

triplets.push_back(triplet_t(0   ,0   ,diag   )); // A(1,1)

Upvotes: 0

coatless
coatless

Reputation: 20746

The script is using syntax that requires a higher C++ standard than c++98. To fix this, there are two ways to approach this.

  1. Enable C++11 as R favors either compilation under C++98 or C++11 using the C++11 Rcpp plugin. This is preferred way to go about fixing this issue. Add the following to the top of your cpp code:

    // [[Rcpp::plugins(cpp11)]]
    
  2. This fix is non-portable (doesn't move with the code file) and is session-based so it must be either repeated in each new session or set in the .Rprofile. This sets a compilation flag commonly found in ~/.R/Makevars

    Sys.setenv("PKG_CXXFLAGS"="-std=c++11")
    

Upvotes: 2

Related Questions