[Rd] Problem using F77_CALL(dgemm) in a package
Jason Rudy
jcrudy at gmail.com
Tue Feb 22 23:45:42 CET 2011
Thanks for the tip. That API could make my work considerably easier.
Jason
On Tue, Feb 22, 2011 at 7:37 AM, Douglas Bates <bates at stat.wisc.edu> wrote:
> On Sun, Feb 20, 2011 at 4:39 PM, Jason Rudy <jcrudy at gmail.com> wrote:
>> I've never used C++ before, so for this project I think I will stick
>> with just using the BLAS and LAPACK routines directly. Another issue
>> is that I will need to do some sparse matrix computations, for which I
>> am planning to use CSPARSE, at least to begin with. I am interested
>> by RcppArmadillo, and would consider it for future projects. If you
>> don't mind, what in your opinion are the major pros and cons of an
>> RcppArmadillo solution compared to simply using the BLAS or LAPACK
>> routines through the .C interface?
>
> You may want to consider the API exported by the Matrix package that
> allows access to CHOLMOD functions in that package's library. The
> entire CSPARSE library is also included in the Matrix package but most
> of it is not exported because the CHOLMOD functions are generally more
> effective. (Both CHOLMOD and CSPARSE are written by Tim Davis.
> CSPARSE is good code but it was written more for instructional purpose
> than as an "industrial strength" package.)
>
>> On Sun, Feb 20, 2011 at 8:11 AM, Dirk Eddelbuettel <edd at debian.org> wrote:
>>>
>>> On 20 February 2011 at 09:50, Dirk Eddelbuettel wrote:
>>> | There is of course merit in working through the barebones API but in case you
>>> | would consider a higher-level alternative, consider these few lines based on
>>> | RcppArmadillo (which end up calling dgemm() for you via R's linkage to the BLAS)
>>>
>>> PS I always forget that we have direct support in Rcpp::as<> for Armadillo
>>> matrices. The examples reduces to three lines in C++, and you never need to
>>> worry about row or column dimension, or memory allocation or deallocation:
>>>
>>> R> suppressMessages(library(inline))
>>> R> txt <- '
>>> + arma::mat Am = Rcpp::as< arma::mat >(A);
>>> + arma::mat Bm = Rcpp::as< arma::mat >(B);
>>> + return Rcpp::wrap( Am * Bm );
>>> + '
>>> R> mmult <- cxxfunction(signature(A="numeric", B="numeric"),
>>> + body=txt,
>>> + plugin="RcppArmadillo")
>>> R> A <- matrix(1:9, 3, 3)
>>> R> B <- matrix(9:1, 3, 3)
>>> R> C <- mmult(A, B)
>>> R> print(C)
>>> [,1] [,2] [,3]
>>> [1,] 90 54 18
>>> [2,] 114 69 24
>>> [3,] 138 84 30
>>> R>
>>>
>>> Matrices A and B from directly initialise Armadillo matrices, and the result
>>> can be returned directly.
>>>
>>> Hth, Dirk
>>>
>>> --
>>> Dirk Eddelbuettel | edd at debian.org | http://dirk.eddelbuettel.com
>>>
>>
>> ______________________________________________
>> R-devel at r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-devel
>>
>
More information about the R-devel
mailing list