[Rcpp-devel] R.e. First foray into Rcpp and comparison of speed with R

Christian Gunning xian at unm.edu
Fri Sep 2 01:51:02 CEST 2011


On Thu, Sep 1, 2011 at 9:19 AM, Dirk Eddelbuettel <edd at debian.org> wrote:
>
>
> I think (but feel to free me wrong by looking at headers and code :-) that
> operator() and operator[] do the same thing.

I see a huge difference in this set of
examples ( which are painfully repetitive on many levels)
that are identical except for [] vs. (). This is true whether
or not ret appears on the RHS. Am I missing something here?

library(inline);
library(Rcpp);
library(rbenchmark);
src1<-'
int nn=as<int>(n);
NumericVector ret(nn);
for (int i=0; i<nn; i++) {
    for (int j=0; j<nn; j++) {
        ret[i] = ret[i]+1;
        // alternately,
        // ret[i] = 1;
    }
};
return(ret);
'

src1a<-'
int nn=as<int>(n);
NumericVector ret(nn);
for (int i=0; i<nn; i++) {
    for (int j=0; j<nn; j++) {
        ret(i) = ret(i)+1;
        // alternately,
        // ret(i) = 1;
    }
};
return(ret);
'


mysample<-cxxfunction( signature(n="numeric"), src1, plugin='Rcpp')
mysamplea<-cxxfunction( signature(n="numeric"), src1a, plugin='Rcpp')
benchmark(
    bracket = mysample( 1e4),
    parens = mysamplea( 1e4),
    replications = 10
)

     test replications elapsed relative user.self sys.self user.child sys.child
1 bracket           10   1.132 1.000000      1.13        0          0         0
2  parens           10  10.190 9.001767     10.18        0          0

> There is a noNA() wrapper for Rcpp sugar to push performance -- NA checking
> is implemented on access 'because that is how R does' (and our first task to
> reproduce numbers you'd get at the R prompt) but if you know what you are
> doing and are aware of possible pitfalls you can skip this.

Ah, so that's for sugar.  Thanks.
So I *can* do linear access on a NumericMatrix, ala
myMatrix[ computed_index ], which goes in column-major order.
Curiously, I don't see any timing difference here.  Now I'm
confused...

src2<-'
int index, nn=as<int>(n);
NumericMatrix ret(nn, nn);
for (int i=0; i<nn; i++) {
    for (int j=0; j<nn; j++) {
        index = i+j*nn;
        ret[index] = ret[index] + index;
    }
};
return(ret);
'

src2a<-'
int index, nn=as<int>(n);
NumericMatrix ret(nn, nn);
for (int i=0; i<nn; i++) {
    for (int j=0; j<nn; j++) {
        index = i+j*nn;
        ret(i, j) = ret(i, j) + index;
    }
};
return(ret);
'

matfil<-cxxfunction( signature(n="numeric"), src2, plugin='Rcpp')
matfila<-cxxfunction( signature(n="numeric"), src2a, plugin='Rcpp')

benchmark(
    bracket = matfil( 1e4),
    parens = matfila( 1e4),
    replications = 10
)

     test replications elapsed relative user.self sys.self user.child sys.child
1 bracket           10  17.495 1.000000     13.14     4.24          0         0
2  parens           10  17.514 1.001086     13.10     4.31          0         0

-xian

-- 
A man, a plan, a cat, a ham, a yak, a yam, a hat, a canal – Panama!


More information about the Rcpp-devel mailing list