On 30 May 2013 at 20:46, Colin Rundel wrote:
| Hi everyone,
| 
| I have recently run into a strange situation with RcppArmadillo and was 
hoping someone might have some insight. I am working with a moderately sized 
(820x820) covariance matrix in C++ which causes arms::chol to fail. In trying 
to diagnose the issue I have written out the offending matrix using the 
arma_binary format and then read it back into R using a small helper function. 
Once in R, I am able to use the R chol function to decompose the matrix without 
issue (no pivoting).
| 
| This makes no sense to me, as it is my understanding that both R and 
armadillo should be using the system lapack library and hence both should call 
the same dpotrf. To make matters worse, I've also tried the eigen decomposition 
of the same matrix, which results in anything from a crash / system hang / 
garbage when using armadillo (eig_sym) but the same matrix read into R will 
happily return a valid result from the eigen function. Other, usually smaller, 
matrices work fine on the same system.

Not necessarily. You need to follow the actual codepath from R and Arma. 
 
| One further wrinkle, is that this also seems to be somewhat system dependent 
as it occurs on my Ubuntu box (13.04, R 3.0.1, RcppArmadillo 0.3.820) but not 
my Macbook (10.8.3, R 3.0.1, RcppArmadillo 0.3.820). This leads me to believe 
that the issue may be with lapack/blas but I see the same result when using 
base lapack, atlas (3.8.4-9ubuntu1), or openblas (0.2.6-1~exp1).

Well that can easily be the issue, and in that case you should probably
follow up with respective BLAS maintainers.  

The very nice things is that you plug BLAS in and out at essentially zero
cost (apart from the download and apt-get mechanics, say 15 secs each).  So
try

   reference blas (packages libblas3 and liblapack)
   atlas          (package libatlas3-base or a tuned variant)
   open-blas      (package libopenblas-base)  

See if one or more of these behave.

Dirk
 
| I have attached a minimal working example that replicates the issue for my 
system, and I have also included a link to the matrix file as it is somewhat 
large.
| 
| xuntyped binary data, test.R         [Click mouse-2 to save to a file]
| xuntyped binary data, util.cpp       [Click mouse-2 to save to a file]
| 
| https://www.dropbox.com/s/6i1c8xp2dsthq4y/tmp.dat
| 
| -Colin
| 
| -----
| 
| Colin Rundel
| Postdoctoral Associate
| Duke University, Department of Statistical Science
| [email protected]
| www.stat.duke.edu/~cr173/
| 
| 
| ----------------------------------------------------------------------
| _______________________________________________
| Rcpp-devel mailing list
| [email protected]
| https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel
-- 
Dirk Eddelbuettel | [email protected] | http://dirk.eddelbuettel.com
_______________________________________________
Rcpp-devel mailing list
[email protected]
https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel

Reply via email to