Some basic features of 'MUMPS' (Multifrontal Massively Parallel
sparse direct Solver) are wrapped in a class whose methods can be used
for sequentially solving a sparse linear system (symmetric or not)
with one or many right hand sides (dense or sparse).
There is a possibility to do separately symbolic analysis,
LU (or LDL^t) factorization and system solving.
The ordering libraries 'PORD' and 'METIS' are included.
'MUMPS' method was first described in Amestoy et al. (2001)
MUMPS stands for "a MUltifrontal Massively Parallel sparse direct Solver" see more on their official site https://mumps-solver.org/. rmumps bundles MUMPS 5.9.1 (sequential, double precision) with the orderings AMD, AMF, QAMD, PORD and METIS. Currently, it is one of the most competitive direct solvers for sparse matrices. On my CPU (Xenon E5-2609 v2 @ 2.50GHz) I have a speedup ranging from 3 to 16 compared to the default solver from Matrix package. In addition, the precision of the solution is equal or even better than with Matrix.
Example of use:
# prepare random sparse matrix
library(Matrix)
library(rmumps)
n=2000
a=Matrix(0, n, n)
set.seed(7)
ij=sample(1:(n*n), 15*n)
a[ij]=runif(ij)
diag(a)=0
diag(a)=-rowSums(a)
a[1,1]=a[1,1]-1
am=Rmumps$new(a)
b=as.double(a%*%(1:n)) # rhs for an exact solution vector 1:n
# following time includes symbolic analysis, LU factorization and system solving
system.time(x<-am$solve(b))
bb=2*b
# this second time should be much shorter
# as symbolic analysis and LU factorization are already done
system.time(xx<-am$solve(bb))
# compare to Matrix corresponding times
system.time(xm<-solve(a, b))
system.time(xxm<-solve(a, bb))
# compare to Matrix precision
range(x-1:n) # mumps
range(xm-1:n) # Matrix
# matrix inversion
system.time(aminv <- am$inv())
system.time(ainv <- solve(a)) # the same in Matrix
# clean up by hand to avoid possible interference between gc() and
# Rcpp object destructor after unloading this namespace
rm(am)