[Rcpp-devel] R’s random number generator not initialized correctly when called at first

Paul Menzel paulepanter at users.sourceforge.net
Mon Dec 5 00:16:33 CET 2011


Dear Rcpp folks,


I am running an up to date Debian Sid/unstable system and I am
simulating how long a random walks stays negative and I am using Rcpp
and inline.

        library(inline)
        inc <- "
        #include <stdio.h>
        #include <stdlib.h>
        "
        
        rwRcpp <-cxxfunction(signature(ns="numeric", reps="numeric"),
        includes=inc, plugin="Rcpp", body='
          double n = Rcpp::as<double>(ns);
          double rep = Rcpp::as<double>(reps);
        
          double x = 0;
        
          Rcpp::NumericVector s(n + 2);
        
          for (int i = 0; i < rep; i++) {
            double len = 0;
            double csum = 0;
        
            do {
              x = rnorm(1, 0, 1)[0];
              Rprintf("x = %f\\n", x);
        
              csum = csum + x;
            } while ((csum < 0) && (len++ < n));
        
            s[len]++;
          }
        
          return s;
        ')

Now the interesting thing is that the random numbers do not seem to be
random when the function is called at the very beginning.

        $ R
        R version 2.14.0 (2011-10-31)
        Copyright (C) 2011 The R Foundation for Statistical Computing
        ISBN 3-900051-07-0
        Platform: i486-pc-linux-gnu (32-bit)
        
        > source("rw.r")
        > rwRcpp(4,10)
        x = -8.773321
        […]
        x = -8.773321
        [1]  0  0  0  0  0 10
        > rwRcpp(2,5)
        x = -8.773321
        x = -8.773321
        x = -8.773321
        x = -8.773321
        x = -8.773321
        x = -8.773321
        x = -8.773321
        x = -8.773321
        x = -8.773321
        x = -8.773321
        x = -8.773321
        x = -8.773321
        x = -8.773321
        x = -8.773321
        x = -8.773321
        [1] 0 0 0 5

Setting the seed value and it starts working.

        > set.seed(1)
        > rwRcpp(2,5)
        x = -0.626454
        x = 0.183643
        x = -0.835629
        x = 1.595281
        x = 0.329508
        x = -0.820468
        x = 0.487429
        x = 0.738325
        x = 0.575781
        [1] 3 0 1 1

Opening a now R session, I tested also the following.

        > source("rw.r")
        > rwRcpp(2,5)
        x = -8.773321
        x = -8.773321
        x = -8.773321
        x = -8.773321
        x = -8.773321
        x = -8.773321
        x = -8.773321
        x = -8.773321
        x = -8.773321
        x = -8.773321
        x = -8.773321
        x = -8.773321
        x = -8.773321
        x = -8.773321
        x = -8.773321
        [1] 0 0 0 5
        > rnorm(2)
        [1]  0.9520938 -2.0510458
        > rwRcpp(2,5)
        x = -0.438619
        x = 0.606704
        x = -0.658996
        x = -0.082486
        x = -2.124667
        x = -0.127444
        x = 0.447259
        x = -0.363092
        x = -0.822970
        x = 0.387696
        x = 1.454236
        [1] 1 2 0 2

So it looks like, Rcpp is not initializing(?) the random number
generator in R? Could that be the cause?


Thanks,

Paul
-------------- next part --------------
library(inline)
inc <- "
#include <stdio.h>
#include <stdlib.h>
"

rwRcpp <-cxxfunction(signature(ns="numeric", reps="numeric"), includes=inc, plugin="Rcpp", body='
  double n = Rcpp::as<double>(ns);
  double rep = Rcpp::as<double>(reps);

  double x = 0;

  Rcpp::NumericVector s(n + 2);

  for (int i = 0; i < rep; i++) {
    double len = 0;
    double csum = 0;

    do {
      x = rnorm(1, 0, 1)[0];
      Rprintf("x = %f\\n", x);

      csum = csum + x;
    } while ((csum < 0) && (len++ < n));

    s[len]++;
  }

  return s;
')
-------------- next part --------------
A non-text attachment was scrubbed...
Name: not available
Type: application/pgp-signature
Size: 198 bytes
Desc: This is a digitally signed message part
URL: <http://lists.r-forge.r-project.org/pipermail/rcpp-devel/attachments/20111205/2b93d442/attachment.pgp>


More information about the Rcpp-devel mailing list