[R-sig-hpc] parallel random numbers: set.seed(i), rsprng, rlecuyer, [one solution]

Ross Boylan ross at biostat.ucsf.edu
Mon Jan 17 22:52:08 CET 2011

On Tue, 2010-12-14 at 13:51 -0800, Ross Boylan wrote:
> The question is how to generate reproducible streams of parallel
> random
> numbers in a way that is insensitive to the number of nodes used.  If
> I
> can run 10 jobs one time, and 60 the next, I want to get the same
> random
> numbers.
> Scenarios
> A) Each job generates its own random numbers.
> B) A specialized subset of jobs generates random numbers; the subset
> expands as the total number of jobs expands.
> C) A fixed subset, e.g., 5 jobs, are responsible for generating the
> random number.
> C) seems the most likely to be achievable. rsprng, which we use,
> initializes streams with a call that includes the stream number and
> the
> total number of streams.  I don't know if the stream number alone
> determines the sequence, or if the stream number, and possibly
> messaging
> between processes, comes into play.  Even if rsprng came with some
> guarantees, other parallel generators (e.g., rlecuyer) might not.
I have some good news to report with rsprng, based on sprng v2.
Scenario A) is quite possible.  That is, if you have 1,000 simulations
you can operate using 1,000 streams.  You just have to remember to kill
the old stream if the same process is handling multiple simulations.

The stream number one to pass to the initialization routine is the job
or simulation number, rather than the rank of the process on which it is

According to the sprng docs, MPI is used only in in transmitting a seed
to all nodes.  My code set the seed though other means, and so there was
never any messaging.  This made it possible to start the streams at
different times.

I also did a small experiment to see if stream 1 was the same regardless
of whether the total number of streams was 500 or 1000.  It appeared to
be.  However, the documentation provides no guarantees that this is the

Here are some excerpts from my code, which runs under snow:
# master and slave setup
setup0 <- function(){
  # set these in the global environment for use by slaveDo
  nTotalSimulations <<- 1000
  seed <<- 1066

# setup to run on master
# cl is cluster
setupMaster <- function(cl, nsim, fname){
  clusterEvalQ(cl, source("ascertain.R")) # file with this code
  clusterEvalQ(cl, setup0())
  system.time(r <- clusterApplyLB(cl, seq(nsim), "slaveDo"))
  invisible(save(r, file=fname))

# sample slave job
# uses global nTotalSimulation and seed
slaveDo <- function(k){
  init.sprng(nTotalSimulations, k-1, seed) 
  # do and return simulation results

cl <- getMPIcluster()
if (mpi.comm.rank(0) == 0) {
  setupMaster(cl, 1000, "r1.RData")

More information about the R-sig-hpc mailing list