I have the following Rcpp code
#include <RcppArmadillo.h>
#define _USE_MATH_DEFINES
#include<cmath>
#include <boost/math/special_functions/bessel.hpp>
#include <mlpack.h>
#include <mlpack/methods/kmeans.hpp>
// [[Rcpp::depends(RcppEnsmallen)]]
// [[Rcpp::depends(mlpack)]]
// [[Rcpp::plugins(cpp17)]]
// [[Rcpp::depends(RcppArmadillo)]]
// [[Rcpp::plugins(openmp)]]
// [[Rcpp::depends(BH)]]
using namespace arma;
using namespace Rcpp;
#define ARMA_USE_ARPACK 1
#define ARMA_USE_SUPERLU 1
#define tolr 1e-5
#define crossprod(x) symmatu((x).t() * (x))
#define tcrossprod(x) symmatu((x) * (x).t())
#define pi_sq 9.869604401089358
// [[Rcpp::export]]
Rcpp::List NNGP_fit(const int n){
sp_mat A = sprandu<sp_mat>(n, n, 0.1);
sp_mat B = A.t()*A;
vec eigval;
mat eigvec;
eigs_sym(eigval, eigvec, B, n-1); // find n-1 eigenvalues/eigenvectors
Rcpp::Rcout<<"flag eig"<<endl;
superlu_opts opts;
opts.symmetric = true;
opts.equilibrate =true;
vec mu_f,
y(n, fill::randn);
bool status = spsolve(mu_f,B, y, "superlu",opts);
Rcpp::Rcout<<"status= "<<status<<endl;
return Rcpp::List::create(Rcpp::Named("covmat") =wrap(B),
Rcpp::Named("eigval") =wrap(eigval),
Rcpp::Named("eigvec") =wrap(eigvec),
Rcpp::Named("y") =wrap(y),
Rcpp::Named("mu_f") =wrap(mu_f)
);
}
I am getting the following output
flag eig
Error: spsolve(): use of SuperLU must be enabled
That is, the code stops working at the spsolve
command but executes eigs_sym
.
However, removing mlpack
by deleting the following solves the problem.
#include <mlpack.h>
#include <mlpack/methods/kmeans.hpp>
// [[Rcpp::depends(RcppEnsmallen)]]
// [[Rcpp::depends(mlpack)]]
So I believe I have linked superlu
correctly.
This is only a reproducible part of the code and I want to use mlpack
in my project. How can I solve this issue? I also need functions from the boost
library for the project and kept that header as well.
When you work with SuperLU you are leaving the (very convenient) world of header-only libraries such as the ones you use here:
But sparse linear algebra needs linking. This is explained and demonstrated over by an Rcpp Gallery article from two of the RcppArmadillo authors.
In general it is often advisable to solve smaller subset problems first. You have not real mlpack
issue here: you have an issue using RcppArmadillo
with sparse decomposition and for that you to link with -lsuperlu
as the Rcpp Gallery describes and shows. I would try to get that to work, and then try to proceed to mlpack
issues.
PS 1 Your code has another issue preventing it from running: the behavior-altering #define
must come before the #include
so please move
#define ARMA_USE_ARPACK 1
#define ARMA_USE_SUPERLU 1
to the very top. It then works as the run-time has a normal coding issue:
> Rcpp::sourceCpp("soquestion20241227.cpp")
> set.seed(123)
> NNGP_fit(3)
flag eig
status= 0
Error in getClass("dgCMatrix") : “dgCMatrix” is not a defined class
Note that in order to get this far I had to enable linking with SuperLU via the envoirnment variable, and I validated against the Rcpp Gallery example code:
> Sys.setenv("PKG_LIBS"="-lsuperlu")
> Rcpp::sourceCpp("gallery20170217.cpp")
> superLuDemo()
Done.
>
>
It is important to set this before the first time you call Rcpp::sourceCpp()
as some of these things appear to be cached per session.
PS 2 Here is a reproducible success. Note that I load the Matrix
package (for sparse matrix support) and set the seed to be reproducible.
> library(Matrix)
> Sys.setenv("PKG_LIBS"="-lsuperlu")
> Rcpp::sourceCpp("soquestion20241227.cpp")
> set.seed(123)
> NNGP_fit(5)
flag eig
status= 0
$covmat
5 x 5 sparse Matrix of class "dgCMatrix"
[1,] 0.0827008 . . . 0.226699
[2,] . . . . .
[3,] . . . . .
[4,] . . . . .
[5,] 0.2266988 . . . 0.788687
$eigval
[,1]
[1,] 7.05517e-33
[2,] 6.13069e-17
[3,] 1.61746e-02
[4,] 8.55213e-01
$eigvec
[,1] [,2] [,3] [,4]
[1,] 1.97882e-16 -7.23362e-15 9.59537e-01 -2.81582e-01
[2,] 1.64671e-01 7.22721e-01 5.40204e-15 2.43186e-17
[3,] -3.27270e-02 -6.76148e-01 -5.07065e-15 -1.23420e-17
[4,] -9.85805e-01 1.43172e-01 1.05140e-15 1.94107e-18
[5,] 3.29803e-17 2.09006e-15 -2.81582e-01 -9.59537e-01
$y
[,1]
[1,] 0.0699539
[2,] 0.9767215
[3,] 2.1650415
[4,] -1.8262054
[5,] 0.5804146
$mu_f
[,1]
>