Rcpp call F77 BLAS/LAPACK (eigen_sym_cpp.cpp)

#include <Rcpp.h>
#include <R_ext/Lapack.h>
#include <R_ext/BLAS.h>

//[[Rcpp::export]]
Rcpp::List eigen_sym_cpp(Rcpp::NumericMatrix X, int num_eig = -1, bool eigenvalues_only = false, double tol = 1.5e-8)
{
  Rcpp::NumericMatrix A = Rcpp::clone(X); // perform deep copy of input
  char jobz = eigenvalues_only?'N':'V', range = (num_eig == -1)?'A':'I', uplo = 'U';
  int N = A.nrow(), lda = std::max(1, N), il = 1, iu = (num_eig == -1)?N:num_eig;
  int m = (range == 'A')?N:(iu-il+1), ldz = std::max(1, N), lwork = std::max(1, 26*N), liwork = std::max(1, 10*N), info = 0;
  double abstol = tol, vl = 0.0, vu = 0.0;
  Rcpp::IntegerVector isuppz(2 * std::max(1, m)), iwork(std::max(1, lwork));
  Rcpp::NumericVector work(std::max(1, lwork)), W(N);
  Rcpp::NumericMatrix Z(ldz, std::max(1, m));
  F77_CALL(dsyevr)(&jobz, &range, &uplo, &N, A.begin(), &lda, &vl, &vu, &il, &iu, &abstol,
           &m, W.begin(), Z.begin(), &ldz, isuppz.begin(), work.begin(), &lwork, iwork.begin(), &liwork, &info);
  return Rcpp::List::create(Rcpp::Named("info") = info,
                            Rcpp::Named("vectors") = Z,
                            Rcpp::Named("values") = W);
}