[Rcpp-devel] Rcpp bug with rgamma

[Romain Francois]

Le 03/10/10 23:29, Dirk Eddelbuettel a écrit :
>
> Colin,
>
> Just as one last comment, the issue at hand is described in Section 6.3 of
> Writing R Extensions:
>
>     6.3 Random number generation
>     ============================
>
>     The interface to R's internal random number generation routines is
>
>          double unif_rand();
>          double norm_rand();
>          double exp_rand();
>
>     giving one uniform, normal or exponential pseudo-random variate.
>     However, before these are used, the user must call
>
>          GetRNGstate();
>
>     and after all the required variates have been generated, call
>
>          PutRNGstate();
>
>     These essentially read in (or create) `.Random.seed' and write it out
>     after use.
>
>     [...]
>
> We had not made that sufficiently clear. What you encountered with rgamma was
> an internal algorithm using a rejection sampleing scheme calling unif_rand()
> and exp_rand() ... but creating an infinite loop at it would always end up
> with the same draws.
>
> I amended the 'Rcpp-sugar' vignette as well as the header files to point
> users to this Section 6.3 to make sure GetRNGstate() and PutRNGstate() are
> used.
>
> Cheers,  Dirk

In addition, we have the RNGScope class, whose constructor calls 
GetRNGstate and destruvctor calls PutRNGstate, so that you can do :

fx <- cxxfunction( , '
	RNGScope scope ;
	NumericVector x = rgamma( 10, 1, 1 ) ;
	return x ;
', plugin = "Rcpp" )

fx()

Romain


-- 
Romain Francois
Professional R Enthusiast
+33(0) 6 28 91 30 30
http://romainfrancois.blog.free.fr
|- http://bit.ly/cCmbgg : Rcpp 0.8.6
|- http://bit.ly/bzoWrs : Rcpp svn revision 2000
`- http://bit.ly/b8VNE2 : Rcpp at LondonR, oct 5th