Using SuperLU sparse solver with RcppArmadillo

Question

I am trying to use the SparseLU solver from armadillo (http://arma.sourceforge.net/docs.html#spsolve) through RcppArmadillo:

#define ARMA_USE_SUPERLU
// [Rcpp::depends(RcppArmadillo)]
#include <RcppArmadillo.h>

// [[Rcpp::export]]
arma::vec sp_solver(arma::sp_mat K, arma::vec x) {
  arma::superlu_opts opts;
  opts.symmetric = true;
  arma::vec res;
  arma::spsolve(res, K, x, "superlu", opts);
  return res;
}

/*** R
library(Matrix)
K <- sparseMatrix(i = c(1, 2, 1), j = c(1, 2, 2), x = c(1, 1, 0.5), symmetric = TRUE)
x <- runif(2)
sp_solver(K, x)
*/

I get the error undefined reference to 'superlu_free'. I guess I'm missing some library linking. Any idea how to solve this issue?

---

I'm on Windows 10.

Answer 1

RcppArmadillo is super-convenient and I use it myself all the time. Because all Rcpp* code is going to be called from R, we can assume some functionality to be present.

This includes LAPACK and BLAS, and explains why we can use RcppArmadillo "without linking" even when the Armadillo documentation clearly states you need LAPACK and BLAS. Why? Well because R already has LAPACK and BLAS.

(As an aside, this lead to considerable early issues if and only if R was built with its own subset of LAPACK as certain complex valued routines were not available. Baptiste was hit quite a bit by this as his packages need(ed) these. Over the years, Brian Ripley was most helpful in adding those missing routines to R's LAPACK. And one never had issues when R was built with an external LAPACK and BLAS as is common for e.g. the Debian/Ubuntu package I maintain.)

Here you want SuperLU. This is optional, and it is your job to ensure this gets linked. In short, it does not just work magically. And it is hard to automate as it involves linking which gets us platform-dependendy and installation requirements we cannot control for easily.

But the question is not new, and there is in fact an entire Rcpp Gallery post on the topic.

Edit: And with one-liner from that post that is suitable for my system, your code works just fine too (once I correct the Rcpp::depends to have the double brackets it needs:

R> Sys.setenv("PKG_LIBS"="-lsuperlu")
R> sourceCpp("answer.cpp")

R> library(Matrix)

R> K <- sparseMatrix(i = c(1, 2, 1), j = c(1, 2, 2), x = c(1, 1, 0.5), symmetric = TRUE)

R> x <- runif(2)

R> sp_solver(K, x)
         [,1]
[1,] 0.264929
[2,] 0.546050
R> 

Corrected Code

#define ARMA_USE_SUPERLU
// [[Rcpp::depends(RcppArmadillo)]]
#include <RcppArmadillo.h>

// [[Rcpp::export]]
arma::vec sp_solver(arma::sp_mat K, arma::vec x) {
  arma::superlu_opts opts;
  opts.symmetric = true;
  arma::vec res;
  arma::spsolve(res, K, x, "superlu", opts);
  return res;
}

/*** R
library(Matrix)
K <- sparseMatrix(i = c(1, 2, 1), j = c(1, 2, 2), x = c(1, 1, 0.5), symmetric = TRUE)
x <- runif(2)
sp_solver(K, x)
*/