eigen decomposition in Rcpp and LAPACK

Sys.setenv("PKG_LIBS" = "$(LAPACK_LIBS) $(BLAS_LIBS) $(FLIBS)")
Rcpp::sourceCpp("eigen_sym_cpp.cpp")

n_row <- 1e4L
n_col <- 1e1L
A <- matrix(rnorm(n_row * n_col), n_row, n_col)
S <- cov(A)

eig_res_cpp1 <- eigen_sym_cpp(S)
eig_res_cpp2 <- eigen_sym_cpp(S, 5L)
eig_res_cpp3 <- eigen_sym_cpp(S, -1L, TRUE)

eig_res_r <- eigen(S)

# check eigenvalues
all.equal(eig_res_cpp1$values, eig_res_r$values) # TRUE
all.equal(eig_res_cpp2$values, eig_res_r$values[1L:5L]) # TRUE
all.equal(eig_res_cpp3$values, eig_res_r$values) # TRUE

# check eigenvectors
sign_vectors_cpp1 <- colSums(eig_res_r$vectors * eig_res_cpp1$vectors)
eig_res_cpp1$vectors <- sweep(eig_res_cpp1$vectors, 2, sign_vectors_cpp1, '*')
all.equal(eig_res_cpp1$vectors, eig_res_r$vectors) # TRUE

# check eigenvectors
sign_vectors_cpp2 <- colSums(eig_res_r$vectors[ , 1L:5L] * eig_res_cpp2$vectors)
eig_res_cpp2$vectors <- sweep(eig_res_cpp2$vectors, 2, sign_vectors_cpp2, '*')
all.equal(eig_res_cpp2$vectors, eig_res_r$vectors[ , 1L:5L]) # TRUE

library(microbenchmark)
microbenchmark(
  Rcpp_BLAS = eigen_sym_cpp(S),
  R = eigen(S),
  times = 100L
)
# Unit: microseconds
#       expr      min       lq      mean    median       uq      max neval
#  Rcpp_BLAS  788.180  852.546  872.9262  860.4455  874.196 1000.293   100
#          R 1203.336 1271.505 1363.7922 1281.4520 1391.457 2432.417   100