Using Armadillo with SuperLU

Qiang Kou and Dirk Eddelbuettel — written Feb 18, 2017 — source

Armadillo is a very versatile C++ library for linear algebra, brought to R via the RcppArmadillo package. It has proven to be very useful and popular, and is (as of February 2017) used by well over 300 CRAN packages as indicated by the reverse depends / linking-to on its CRAN page. Well over a dozen earlier posts on this Rcpp Gallery site also demonstrate the popularity of the package as well as different usage patterns.

Armadillo has a core focus on dense matrices, but continues to expand its capabilities for sparse matrices. Basic operation are supported directly via the templated header files, along with calls back into the default (dense) LAPACK and BLAS libraries we can access easily as R uses them.

Armadillo also supports dedicated sparse solvers via the SuperLU package. However, this requires access to the SuperLU library. Many Linux distributions ship it, see e.g. the Debian package page and the Ubuntu package page; there is also a homebrew recipe for OS X / macOS (or other systems using brew). As of this writing, the version in the current Ubuntu release is behind the version Debian. But it is the 5.2.* version that is in Debian that is also required with current Armadillo versions 7.700.* so we prepared ‘backports’ via Dirk’s PPA repo.)

Recently, a GitHub issue ticket asked how to use SuperLU along with RcppArmadillo. This post essentially repeats the main answer, which was spread out over multiple posts, in a single article.

In a nutshell, two things need to happen:

  1. One needs to define the required variable ARMA_USE_SUPERLU which has to be done before the Armadillo headers are included. One possibility (shown below) is a #define statement right in the code file.

  2. One needs to tell the linker to use the SuperLU library. This step is of course not perfectly portable, and does require that the library be installed. (A genuinely portable solution would either need to test for presence of SuperLU, or include its sources. Both aspects are beyond the scope of this post._

Sys.setenv("PKG_LIBS"="-lsuperlu")

Now that R knows about this library, we can present the code requiring it:

// Important: this definition ensures Armadillo enables SuperLU
#define ARMA_USE_SUPERLU 1

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

using namespace arma;

// [[Rcpp::export]]
void superLuDemo() {
    sp_mat A = sprandu<sp_mat>(1000, 1000, 0.1);

    vec b = randu<vec>(1000);
    mat B = randu<mat>(1000, 5);

    vec x = spsolve(A, b);  // solve one system
    mat X = spsolve(A, B);  // solve several systems

    bool status = spsolve(x, A, b);  // use default solver
    if (status == false)  { Rcpp::Rcout << "no solution" << endl; }

    spsolve(x, A, b, "lapack");   // use LAPACK  solver
    spsolve(x, A, b, "superlu");  // use SuperLU solver

    Rcpp::Rcout << "Done.\n";
}

This code snippet defines a function superLuDemo() which we can call from R:

superLuDemo()
Done.

As the data generated here is random, we did not bother printing the (dense) result vector of length 1000.

tags: armadillo  sparse 

Related Articles