[Rcpp-devel] NumericVector Double mismatch when indexing an array

Goldfeld, Keith Keith.Goldfeld at nyumc.org
Tue Sep 25 21:19:28 CEST 2012


Yes - I agree that iterating is not the best way to go, but I am using this for an agent based simulation where iteration is really my only option.  I won't be generating more than a few thousand records, so the time shouldn't be much of a factor.  Of course, in R it is painfully slow - so that is why I am going the Rcpp route.

-----Original Message-----
From: dmbates at gmail.com [mailto:dmbates at gmail.com] On Behalf Of Douglas Bates
Sent: Tuesday, September 25, 2012 3:13 PM
To: Rodney Sparapani
Cc: Goldfeld, Keith; rcpp-devel at lists.r-forge.r-project.org
Subject: Re: [Rcpp-devel] NumericVector Double mismatch when indexing an array

On Tue, Sep 25, 2012 at 1:53 PM, Rodney Sparapani <rsparapa at mcw.edu> wrote:
> On 09/25/2012 01:37 PM, Goldfeld, Keith wrote:
>
>>>  code<- 'Rcpp::RNGScope scope;
>>
>>
>> +                    Rcpp::NumericVector rn(5);
>>
>> +                    for (int i=0; i<  5; i++) {
>>
>> +                        rn(i) = rnorm(1,0,1);
>>
>> +                    }
>>
>> +                    return Rcpp::wrap(rn);
>>
>> +  '
>>
>>>
>>
>>>  fun<- cxxfunction(body=code, plugin="Rcpp")
>
>
> You just need an explicit type conversion like so:
>
> funk2 <- cxxfunction(body='RNGScope scope;
>           NumericVector rn(5);
>           for (int i=0; i < 5; i++) rn[i] = as<double>(rnorm(1,0,1));
>           return wrap(rn);',
>           includes="using namespace Rcpp;\n", plugin="Rcpp")
>
> funk2()

That works but you wouldn't want to use it as a general model because you are creating vectors of length 1 then dereferencing the first element in each of these vectors and copying it to another vector.
Obviously with n=5 the difference in execution time will be minimal but with n=5 million it won't be.  Using the Rcpp sugar function rnorm will be the easiest approach unless, like me, you get a little queasy when using Rcpp sugar functions just because it is so hard to pick them apart and see what is actually happening.

I would probably end up going back to the R API and calling norm_rand or Rf_rnorm.  That latter returns a double from two double arguments, mu and sigma.


More information about the Rcpp-devel mailing list