I am trying to compute eigenvalues in C++ using the Armadillo function eig_sym via RcppArmadillo. The results are not entirely the same as the output of the R function eigen():
In R:
set.seed(1) X=matrix(sample(1:25), 5) X # [,1] [,2] [,3] [,4] [,5] #[1,] 7 18 4 19 6 #[2,] 9 22 3 25 16 #[3,] 14 12 24 8 2 #[4,] 20 11 21 23 15 #[5,] 5 1 13 10 17 Xcov=cov(X) eigen(Xcov)$values #[1] 1.585160e+02 7.475128e+01 5.938207e+01 3.250609e+00 -4.293203e-15 In C++:
//————————–C++ Code ————————
#include <RcppArmadillo.h> // [[Rcpp::depends(RcppArmadillo)]] using namespace Rcpp; using namespace arma; // [[Rcpp::export]] arma::vec eigenval(arma::mat M) { arma::vec values=arma::eig_sym(M); return values; } I then compile the code in R using sourceCpp():
sourceCpp("eigenval.cpp") sort(as.vector(eigenval(Xcov)), TRUE) #1.585160e+02 7.475128e+01 5.938207e+01 3.250609e+00 1.065789e-14 Compare the above results with the R function result (repeated below), we see that the last value is different.
#[1] 1.585160e+02 7.475128e+01 5.938207e+01 3.250609e+00 -4.293203e-15 Not sure if this is the right place to ask, but I wonder if anyone has any idea about the difference.
Contents
hide
Best Answer
Numerical errors. The smallest eigenvalue is flirting with machine accuracy.