[Rcpp-devel] Segfault error during simulation in Rcpp
Douglas Bates
bates at stat.wisc.edu
Mon May 6 15:34:54 CEST 2013
The segfaults seem to be related to PutRNGstate and I don't see that you
have declared an instance of the class that causes the RNGstate to be
accessed and restored (I have forgotten the name of the class but it should
be fairly easy to find in the examples). When you use random number
generators from within Rcpp you either need to call getRNGstate/putRNGstate
yourself or to simply declare an instance of this class.
On Mon, May 6, 2013 at 6:15 AM, Matteo Fasiolo <matteo.fasiolo at gmail.com>wrote:
> Hi All,
>
> I am trying to simulate from the following model:
>
> Y_t = Pois( phi * N_t );
> N_t = r * N_{t-1} * exp( -N_{t-1}^theta + e_t )
> e_t ~ Norm(0, sigma)
>
> so I have written a Rcpp function with prototype:
>
> genRickerCpp(int days, int nSimul, int nBurn, NumericMatrix params)
>
> where:
>
> days = length of each simulated path
> nSimul = number of simulated paths
> nBurn = number of simulations I'm going to discard before storing
> each path
> params = matrix of input parameters that can be either 1 by 4 (in which
> case all the paths are simulated using the same parameters)
> or nSimul by 4 (in which case each path uses a different
> vector of parameters)
>
> This is my code:
>
> #include <Rcpp.h>
>
> using namespace Rcpp;
>
> // [[Rcpp::export]]
> NumericMatrix genRickerCpp(int days, int nSimul, int nBurn, NumericMatrix
> params)
> {
> int nParams = params.ncol();
> int totDays = nBurn + days;
> bool multiParams = false;
>
> if(nParams != 4) stop("Wrong number of parameters");
> if(params.nrow() > 1) { multiParams = true; }
> if(multiParams == true && params.nrow() != nSimul)
> stop("Number of parameters vectors is different from the number of
> simulations");
>
> double r = exp(params(0, 0));
> double theta = exp(params(0, 1));
> double sigma = exp(params(0, 2));
> double phi = exp(params(0, 3));
>
> NumericVector procNoise( rnorm( totDays * nSimul ) );
> NumericVector initState( runif( nSimul ) );
> NumericMatrix output( nSimul, days );
>
> NumericVector::iterator noiseIter = procNoise.begin();
> NumericVector::iterator initIter = initState.begin();
>
> double currState;
>
> for(int iRow = 0; iRow < nSimul; iRow++, initIter++)
> {
>
> if( multiParams == true )
> {
> r = exp(params(iRow, 0));
> theta = exp(params(iRow, 1));
> sigma = exp(params(iRow, 2));
> phi = exp(params(iRow, 3));
> }
>
> currState = *initIter;
>
> for(int iCol = 1; iCol <= nBurn; iCol++, noiseIter++){
> currState = r * currState * exp( - pow( currState, theta ) +
> *noiseIter * sigma );
> }
>
> output(iRow, 0) = rpois(1, phi * currState)[0];
>
> for(int iCol = 1; iCol < days; iCol++, noiseIter++){
> currState = r * currState * exp( - pow( currState, theta ) +
> *noiseIter * sigma );
> output(iRow, iCol) = rpois(1, phi * currState)[0];
> }
>
> }
>
> return output;
>
> }
>
>
> the function seems to work well, I tried to compare the output the an
> equivalent R function
> and I get the same results. The problem is that if I run it a lot of times:
>
> library(Rcpp)
> sourceCpp("~/Desktop/genRickerCpp.cpp")
>
> for(ii in 1:10^6){
> data <- genRickerCpp(days = 1, nSimul = 1, nBurn = 1,
> params = matrix(log(c(r = exp(3.8), theta = 1, sigma
> = 0.3, phi = 10)), 1, 4))
> data <- as.numeric(data)
> }
>
> occasionally R crashes with error:
>
> *** caught segfault ***
> address 0x28, cause 'memory not mapped'
>
> the strange thing is that in most cases I can call it 10^6 times without
> any error.
> I tried to go through the code in gdb, but I didn't see anything wrong. I
> also ran the previous
> R code in valgrind, and there I get the following errors while the code is
> running:
>
> ==3031== Invalid read of size 1
> ==3031== at 0x4EF7BF1: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4EF7CE5: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4EF8894: Rf_duplicate (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4EA79E7: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F25B08: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F2779F: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F259AE: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F2791F: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F259AE: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F2958C: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F259AE: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F61FA2: Rf_ReplIteration (in /usr/lib/R/lib/libR.so)
> ==3031== Address 0xebd20c8 is 1,688 bytes inside a block of size 1,968
> free'd
> ==3031== at 0x4C2A82E: free (in
> /usr/lib/valgrind/vgpreload_memcheck-amd64-linux.so)
> ==3031== by 0x4F644AC: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F67965: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F6933E: Rf_allocVector (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4E7468C: PutRNGstate (in /usr/lib/R/lib/libR.so)
> ==3031== by 0xDD2017D: sourceCpp_74073_genRickerCpp (random.h:71)
> ==3031== by 0x4EEBEC5: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F25BCC: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F28E0C: Rf_applyClosure (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F2588F: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F2779F: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F259AE: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031==
> ==3031== Invalid read of size 8
> ==3031== at 0x4EF7F89: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4EF7CE5: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4EF8894: Rf_duplicate (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4EA79E7: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F25B08: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F2779F: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F259AE: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F2791F: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F259AE: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F2958C: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F259AE: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F61FA2: Rf_ReplIteration (in /usr/lib/R/lib/libR.so)
> ==3031== Address 0xebd20e8 is 1,720 bytes inside a block of size 1,968
> free'd
> ==3031== at 0x4C2A82E: free (in
> /usr/lib/valgrind/vgpreload_memcheck-amd64-linux.so)
> ==3031== by 0x4F644AC: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F67965: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F6933E: Rf_allocVector (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4E7468C: PutRNGstate (in /usr/lib/R/lib/libR.so)
> ==3031== by 0xDD2017D: sourceCpp_74073_genRickerCpp (random.h:71)
> ==3031== by 0x4EEBEC5: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F25BCC: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F28E0C: Rf_applyClosure (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F2588F: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F2779F: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F259AE: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031==
> ==3031== Invalid read of size 8
> ==3031== at 0x4EF7FA8: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4EF7CE5: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4EF8894: Rf_duplicate (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4EA79E7: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F25B08: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F2779F: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F259AE: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F2791F: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F259AE: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F2958C: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F259AE: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F61FA2: Rf_ReplIteration (in /usr/lib/R/lib/libR.so)
> ==3031== Address 0xebd20f8 is 1,736 bytes inside a block of size 1,968
> free'd
> ==3031== at 0x4C2A82E: free (in
> /usr/lib/valgrind/vgpreload_memcheck-amd64-linux.so)
> ==3031== by 0x4F644AC: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F67965: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F6933E: Rf_allocVector (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4E7468C: PutRNGstate (in /usr/lib/R/lib/libR.so)
> ==3031== by 0xDD2017D: sourceCpp_74073_genRickerCpp (random.h:71)
> ==3031== by 0x4EEBEC5: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F25BCC: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F28E0C: Rf_applyClosure (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F2588F: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F2779F: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F259AE: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031==
> ==3031== Invalid read of size 8
> ==3031== at 0x4EF7FC6: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4EF7CE5: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4EF8894: Rf_duplicate (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4EA79E7: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F25B08: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F2779F: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F259AE: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F2791F: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F259AE: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F2958C: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F259AE: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F61FA2: Rf_ReplIteration (in /usr/lib/R/lib/libR.so)
> ==3031== Address 0xebd20d0 is 1,696 bytes inside a block of size 1,968
> free'd
> ==3031== at 0x4C2A82E: free (in
> /usr/lib/valgrind/vgpreload_memcheck-amd64-linux.so)
> ==3031== by 0x4F644AC: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F67965: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F6933E: Rf_allocVector (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4E7468C: PutRNGstate (in /usr/lib/R/lib/libR.so)
> ==3031== by 0xDD2017D: sourceCpp_74073_genRickerCpp (random.h:71)
> ==3031== by 0x4EEBEC5: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F25BCC: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F28E0C: Rf_applyClosure (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F2588F: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F2779F: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F259AE: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031==
> ==3031== Invalid read of size 8
> ==3031== at 0x4EF7F7C: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4EF7CE5: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4EF8894: Rf_duplicate (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4EA79E7: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F25B08: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F2779F: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F259AE: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F2791F: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F259AE: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F2958C: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F259AE: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F61FA2: Rf_ReplIteration (in /usr/lib/R/lib/libR.so)
> ==3031== Address 0xebd20f0 is 1,728 bytes inside a block of size 1,968
> free'd
> ==3031== at 0x4C2A82E: free (in
> /usr/lib/valgrind/vgpreload_memcheck-amd64-linux.so)
> ==3031== by 0x4F644AC: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F67965: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F6933E: Rf_allocVector (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4E7468C: PutRNGstate (in /usr/lib/R/lib/libR.so)
> ==3031== by 0xDD2017D: sourceCpp_74073_genRickerCpp (random.h:71)
> ==3031== by 0x4EEBEC5: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F25BCC: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F28E0C: Rf_applyClosure (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F2588F: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F2779F: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F259AE: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031==
> ==3031== Invalid read of size 1
> ==3031== at 0x4EF7E34: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4EF7CE5: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4EF8894: Rf_duplicate (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4EA79E7: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F25B08: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F2779F: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F259AE: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F2791F: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F259AE: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F2958C: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F259AE: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F61FA2: Rf_ReplIteration (in /usr/lib/R/lib/libR.so)
> ==3031== Address 0xebd20c8 is 1,688 bytes inside a block of size 1,968
> free'd
> ==3031== at 0x4C2A82E: free (in
> /usr/lib/valgrind/vgpreload_memcheck-amd64-linux.so)
> ==3031== by 0x4F644AC: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F67965: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F6933E: Rf_allocVector (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4E7468C: PutRNGstate (in /usr/lib/R/lib/libR.so)
> ==3031== by 0xDD2017D: sourceCpp_74073_genRickerCpp (random.h:71)
> ==3031== by 0x4EEBEC5: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F25BCC: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F28E0C: Rf_applyClosure (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F2588F: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F2779F: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F259AE: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031==
> ==3031== Invalid read of size 1
> ==3031== at 0x4EF7C27: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4EF7CE5: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4EF8894: Rf_duplicate (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4EA79E7: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F25B08: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F2779F: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F259AE: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F2791F: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F259AE: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F2958C: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F259AE: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F61FA2: Rf_ReplIteration (in /usr/lib/R/lib/libR.so)
> ==3031== Address 0xebd20c9 is 1,689 bytes inside a block of size 1,968
> free'd
> ==3031== at 0x4C2A82E: free (in
> /usr/lib/valgrind/vgpreload_memcheck-amd64-linux.so)
> ==3031== by 0x4F644AC: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F67965: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F6933E: Rf_allocVector (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4E7468C: PutRNGstate (in /usr/lib/R/lib/libR.so)
> ==3031== by 0xDD2017D: sourceCpp_74073_genRickerCpp (random.h:71)
> ==3031== by 0x4EEBEC5: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F25BCC: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F28E0C: Rf_applyClosure (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F2588F: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F2779F: ??? (in /usr/lib/R/lib/libR.so)
> ==3031== by 0x4F259AE: Rf_eval (in /usr/lib/R/lib/libR.so)
> ==3031==
>
>
> The only thing that I understand is that probably there is an invalid read
> somewhere,
> but I went through the code several times and I even re-wrote everything
> from
> scratch and I can't find anything wrong. Given that I'm a C++/Rcpp
> beginner
> I guess that probably I'm doing a systematic error in my code, maybe you
> can point it
> out to me.
> Finally I'm running my code on Ubuntu 12.04 (precise) 64-bit and R version
> 2.15.3.
> Thanks a lot.
>
> Matteo
>
> _______________________________________________
> Rcpp-devel mailing list
> Rcpp-devel at lists.r-forge.r-project.org
> https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.r-forge.r-project.org/pipermail/rcpp-devel/attachments/20130506/b06760e1/attachment-0001.html>
More information about the Rcpp-devel
mailing list