On 14 June 2012 at 10:14, Julian Smith wrote:
| That's useful to know about what svd_econ() exactly does. I will give that a
| shot and report back.
It doesn't buy anything -- same example as before but now with svd_econ().
R> require(RcppArmadillo)
Loading required package: RcppArmadillo
Loading required package: Rcpp
R> require(inline)
Loading required package: inline
R> require(rbenchmark)
Loading required package: rbenchmark
R>
R> arma.code <- 'using namespace arma; NumericMatrix Xr(Xs); int n = Xr.nrow(),
k = Xr.ncol(); mat X(Xr.begin(), n, k, false); mat U; vec s; mat V; svd_econ(U,
s, V, X); return wrap(s);'
R> rcppsvd <- cxxfunction(signature(Xs="numeric"),
body=arma.code,plugin="RcppArmadillo")
R>
R> N <- 1000
R> A <- matrix(rnorm(N^2), N)
R>
R> res <- benchmark(rcppsvd(A), svd(A), replications=10)
R> print(res)
test replications elapsed relative user.self sys.self user.child
sys.child
1 rcppsvd(A) 10 108.502 9.45304 113.343 1.664 0
0
2 svd(A) 10 11.478 1.00000 18.946 3.788 0
0
R>
Still as bad as before.
Dirk
|
| Has anyone looked at integrating something like SLEPc, Anasazi(via Trilinos)
or
| ARPACK++ into rcpp? These would be some really cool tools to have available.
|
| http://www.grycap.upv.es/slepc/description/summary.htm
|
| http://trilinos.sandia.gov/packages/anasazi/
|
| http://www.ime.unicamp.br/~chico/arpack++/
|
| Someone actually wrote a wrapper for ARPACK++ for Eigen
|
| https://github.com/beam2d/arpaca/blob/master/README.md
| [cleardot]
|
|
| On Thu, Jun 14, 2012 at 9:04 AM, Conrad Sand <conrads...@gmail.com> wrote:
|
|
| On Jun 15, 2012 12:11 AM, "Dirk Eddelbuettel" <e...@debian.org> wrote:
| > Thanks for that earlier hint re 'thin' and 'full' SVDs.
|
| armadillo has the standard svd() and the thin version too: svd_econ().
|
| > Conrad, any interest in switching to dgesdd?
|
| yes, but probably not simply changing svd() directly to dgedd(). lapack
| gave the function a different name for a reason: the results may be
| slightly different. if armadillo starts giving different results all of a
| sudden, there would be a lot of displeased people. a big no no, given that
| armadillo is used for critical stuff.
|
| I'll probably add an option to svd() to optionally use dgedd().
|
| I've done something very similar for eig_sym() in armadillo 3.2, where an
| alternative faster algorithm for eigen decomposition can be optionally
| used.
|
|
| > Dirk, at useR and across the room from Doug
| >
| > | > On Wed, Jun 13, 2012 at 4:07 PM, Douglas Bates <ba...@stat.wisc.edu>
| wrote:
| > | >>
| > | >> On Wed, Jun 13, 2012 at 5:16 PM, Dirk Eddelbuettel
<e...@debian.org>
| wrote:
| > | >> >
| > | >> > On 13 June 2012 at 15:05, Julian Smith wrote:
| > | >> > | I agree that RcppEigen is a little bit faster, but ease of use
| is
| > | >> > important to
| > | >> > | me, so I feel like RcppArmadillo might win out in my
| application.
| > | >> >
| > | >> > Yup, that my personal view too.
| > | >> >
| > | >> > | | RcppArmadillo will use the very same LAPACK and BLAS libs
your
| R
| > | >> > session
| > | >> > | | uses. So MKL, OpenBlas, ... are all options. Eigen actually
| has its
| > | >> > own
| > | >> > | code
| > | >> > | | outperforming LAPACK, so it doesn't as much there.
| > | >> > |
| > | >> > | Why do you think R outperforms RcppArmadillo in this example
| below?
| > | >> > Anyway to
| > | >> > | speed this up?
| > | >> >
| > | >> > That is odd. "I guess it shouldn't." I shall take another look --
| as I
| > | >> > understand it both should go to the same underlying Lapack
| routine. I
| > | >> > may
| > | >> > have to consult with Conrad on this.
| > | >> >
| > | >> > Thanks for posting a full and reproducible example!
| > | >> >
| > | >> > Dirk
| > | >> >
| > | >> > | require(RcppArmadillo)
| > | >> > | require(inline)
| > | >> > |
| > | >> > | arma.code <- '
| > | >> > | using namespace arma;
| > | >> > | NumericMatrix Xr(Xs);
| > | >> > | int n = Xr.nrow(), k = Xr.ncol();
| > | >> > | mat X(Xr.begin(), n, k, false);
| > | >> > | mat U;
| > | >> > | vec s;
| > | >> > | mat V;
| > | >> > | svd(U, s, V, X);
| > | >> > | return wrap(s);
| > | >> > | '
| > | >>
| > | >> Because the arma code is evaluating the singular vectors (U and V)
| as
| > | >> well as the singular values (S) whereas the R code is only
| evaluating
| > | >> the singular values. There is considerably more effort required to
| > | >> evaluate the singular vectors in addition to the singular values.
| > | >>
| > | >> > | rcppsvd <- cxxfunction(signature(Xs="numeric"),
| > | >> > | arma.code,
| > | >> > | plugin=
| "RcppArmadillo")
| > | >> > |
| > | >> > | A<-matrix(rnorm(5000^2), 5000)
| > | >> > |
| > | >> > | > system.time(rcppsvd(A))
| > | >> > | user system elapsed
| > | >> > | 1992.406 4.862 1988.737
| > | >> > |
| > | >> > | > system.time(svd(A))
| > | >> > | user system elapsed
| > | >> > | 652.496 2.641 652.614
| > | >> > |
| > | >> > | On Wed, Jun 13, 2012 at 11:43 AM, Dirk Eddelbuettel <
| e...@debian.org>
| > | >> > wrote:
| > | >> > |
| > | >> > |
| > | >> > | On 13 June 2012 at 10:57, Julian Smith wrote:
| > | >> > | | I've been toying with both RcppArmadillo and RcppEigen
| the past
| > | >> > few days
| > | >> > | and
| > | >> > | | don't know which library to continue using. RcppEigen
| seems
| > | >> > really slick,
| > | >> > | but
| > | >> > | | appears to be lacking some of the decompositions I want
| and
| > | >> > isn't nearly
| > | >> > | as
| > | >> > | | fast to code. RcppArmadillo seems about as fast, easier
| to code
| > | >> > up etc.
| > | >> > | What
| > | >> > | | are some of the advantages/disadvantages of both?
| > | >> > |
| > | >> > | That's pretty close. I have been a fan of [Rcpp]
| Armadillo which I
| > | >> > find
| > | >> > | easier to get my head around. Doug, however, moved from
| > | >> > [Rcpp]Armadillo
| > | >> > | to
| > | >> > | [Rcpp]Eigen as it has some things he needs. Eigen should
| have a
| > | >> > "larger"
| > | >> > | API
| > | >> > | than Armadillo, but I find the code and docs harder to
| navigate.
| > | >> > |
| > | >> > | And you should find Eigen to be a little faster. Andreas
| Alfons
| > | >> > went as far
| > | >> > | as building 'robustHD' using RcppArmadillo with a drop-in
| for
| > | >> > RcppEigen (in
| > | >> > | package 'sparseLTSEigen'; both package names from memmory
| and I
| > | >> > may have
| > | >> > | mistyped). He reported a performance gain of around 25%
| for his
| > | >> > problem
| > | >> > | sets. On the 'fastLm' benchmark, we find the fast
| Eigen-based
| > | >> > | decompositions
| > | >> > | to be much faster than Armadillo.
| > | >> > |
| > | >> > | | Can you call LAPACK or BLAS from either? Is there a
| wrapper in
| > | >> > RcppEigen
| > | >> > | to
| > | >> > | | call LAPACK functions? Want some other decomposition
| methods,
| > | >> > dont like
| > | >> > | the
| > | >> > | | JacobiSVD method in Eigen.
| > | >> > |
| > | >> > | You need to differentiate between the Eigen and Armadillo
| docs
| > | >> > _for their
| > | >> > | libraries_ and what happens when you access the Rcpp*
| variant from
| > | >> > R.
| > | >> > |
| > | >> > | RcppArmadillo will use the very same LAPACK and BLAS libs
| your R
| > | >> > session
| > | >> > | uses. So MKL, OpenBlas, ... are all options. Eigen
| actually has
| > | >> > its own
| > | >> > | code
| > | >> > | outperforming LAPACK, so it doesn't as much there.
| > | >> > |
| > | >> > | Hope this helps, Dirk (at useR!)
| > | >> > |
| > | >> > | |
| > | >> > | |
| > | >> >
| ----------------------------------------------------------------------
| > | >> > | | _______________________________________________
| > | >> > | | Rcpp-devel mailing list
| > | >> > | | Rcpp-devel@lists.r-forge.r-project.org
| > | >> > | |
| > | >> > https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/
| rcpp-devel
| > | >> > | --
| > | >> > | Dirk Eddelbuettel | e...@debian.org | http://
| dirk.eddelbuettel.com
| > | >> > |
| > | >> > |
| > | >> >
| > | >> > --
| > | >> > Dirk Eddelbuettel | e...@debian.org |
http://dirk.eddelbuettel.com
| > | >> > _______________________________________________
| > | >> > Rcpp-devel mailing list
| > | >> > Rcpp-devel@lists.r-forge.r-project.org
| > | >> > https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/
| rcpp-devel
| > | >
| > | >
| >
| > --
| > Dirk Eddelbuettel | e...@debian.org | http://dirk.eddelbuettel.com
|
|
|
| ----------------------------------------------------------------------
| _______________________________________________
| Rcpp-devel mailing list
| Rcpp-devel@lists.r-forge.r-project.org
| https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel
--
Dirk Eddelbuettel | e...@debian.org | http://dirk.eddelbuettel.com
_______________________________________________
Rcpp-devel mailing list
Rcpp-devel@lists.r-forge.r-project.org
https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel
