Parallel Computing and Random Numbers
When writing programs that work together on different machines. Programmers encounter numerous types of problems ranging from endianess, file storage format, executable formats and many others. One possible problem is generating random numbers.
Random number generators are implemented in many different ways. I will not aim to discuss the details of generating random numbers here. Just a few solutions to some common problems.
Simple C Program for generating random numbers:
#include <stdio.h>
#include <stdlib.h>
#include <time.h>
int main(int argc, char *argv[])
{
double rn;
int i;
time_t mytime;
srand (time (&mytime));
for (i=0;i<10;i++) {
rn = 1.0 * (rand() / (RAND_MAX + 1.0));
printf("Process 0, random number %d: %.14fn", i, rn);
}
return 0;
}
An obvious parallel program would look like this:
#include <stdio.h>
#include <stdlib.h>
#include <time.h>
#include <mpi.h>
int main(int argc, char *argv[])
{
time_t mytime;
double rn;
int i, myid, nprocs;
/* initialize MPI */
MPI_Init(&argc, &argv);
MPI_Comm_rank(MPI_COMM_WORLD, &myid);
MPI_Comm_size(MPI_COMM_WORLD, &nprocs);
srand(time(&mytime));
for (i=0;i<10;i++) {
rn = 1.0 * (rand() / (RAND_MAX + 1.0));
printf("Process %d, random number %d: %.14fn", myid, i+1, rn);
}
MPI_Finalize();
return 0;
}
What is wrong with the picture? Well, if you look at the output since the random number generators in each process is seeded with more or less the same seed then … they generate the same sequence of random numbers! Oh no! Then if all the numbers are taken together … they are not random!
So what should we do?
#include <stdio.h>
#include <stdlib.h>
#include <time.h>
#include <mpi.h>
int main(int argc, char *argv[])
{
time_t mytime;
double rn;
int i, myid, nprocs;
/* initialize MPI */
MPI_Init(&argc, &argv);
MPI_Comm_rank(MPI_COMM_WORLD, &myid);
MPI_Comm_size(MPI_COMM_WORLD, &nprocs);
srand(time(&mytime) * myid);
for (i=0;i<10;i++) {
rn = 1.0 * (rand() / (RAND_MAX + 1.0));
printf("Process %d, random number %d: %.14fn", myid, i+1, rn);
}
MPI_Finalize();
return 0;
}
The simplest solution is the code above. Notice, we just made the seed vary by using the rank of the process. The output of this solution seems random. But No! Each process generates a set of random numbers that are distributed in the same way. The distribution is not common for the entire communicator! This has dire implications when using these random numbers for monte carlo computations and others that need properly distributed random numbers.
Fortunately, a bunch of generous and smart folks have written a library that allows you to generate random numbers in parallel. The library is pretty powerful that it allows you to tune the parameters of your pseudo-random number generator. This library is called SPRNG.
SPRNG 1.0 provides the user the various SPRNG random number generators each in its own library. For most users this is acceptable, as one rarely uses more than one type of generator in a single program. However, if the user desires this added flexibility, SPRNG 2.0 provides it. In all other respects, SPRNG 1.0 and SPRNG 2.0 are identical. Both versions only uses the GNU Multi Precision (GMP) package for one of their generators. SPRNG 3.0 uses GMP for all generators. SPRNG 4.0 is a C++ version with the GMP package removed. It is not backwards compatible with any of the previous SPRNG versions, except for its default FORTRAN interface.
You can download a package with a pre-built SPRNG 3.0 library and some sample code. It comes with a Makefile for your convenience. You can download the code from here.
Please make sure you have GMP and its development libraries installed before building this code. I wrote a feature on GMP here. In a Fedora system, you can install GMP with the following commands:
yum install gmp gmp-devel
The sample code can be found here:
#include <stdio.h>
#include <mpi.h>
#include "sprng.h"
#define SIMPLE_SPRNG
#define SEED 985456376
int main(int argc, char *argv[])
{
int *stream;
double rn;
int i, myid, nprocs;
/* initialize MPI */
MPI_Init(&argc, &argv);
MPI_Comm_rank(MPI_COMM_WORLD, &myid);
MPI_Comm_size(MPI_COMM_WORLD, &nprocs);
/**
* initial PRNG streams.
* 1 - 48 bit LCRNG.
* 2 - 64 bit LCRNG.
* the steam id is the same as the rank and the
* number of streams is the same as the number of
* processes.
**/
stream = init_sprng(2, myid, nprocs, SEED, SPRNG_DEFAULT);
for (i=0;i<3;i++) {
rn = sprng(stream);
printf("Process %d, random number %d: %.14fn", myid, i+1, rn);
}
/* free random number stream */
free_sprng(stream);
MPI_Finalize();
return 0;
}
See it was not that hard. Now you can write your Monte Carlo simulator. Good luck!
PS. I left a little bug in the sample code above. It will compile and run but I left something that is not normally a best practice when using pseudo-random number generators.
