¿Cómo generar números aleatorios en paralelo?

Question

Quiero generar números pseudoaleatorios en paralelo usando OpenMP, algo como esto:

int i; 
#pragma omp parallel for 
for (i=0;i<100;i++) 
{ 
    printf("%d %d %d\n",i,omp_get_thread_num(),rand()); 
} 
return 0; 

que he probado en ventanas y tengo enorme aumento de velocidad, pero cada hilo generó exactamente los mismos números. Lo probé también en Linux y tuve una gran desaceleración, la versión paralela en el procesador 8core fue aproximadamente 10 veces más lenta que la secuencia, pero cada hilo generó números diferentes.

¿Hay alguna forma de tener aceleración y números diferentes?

Editar 27.11.2010
creo que he resuelto usando una idea de Jonathan Dursi puesto. Parece que el siguiente código funciona rápido tanto en Linux como en Windows. Los números también son pseudoaleatorios. ¿Qué piensa usted al respecto?

int seed[10]; 

int main(int argc, char **argv) 
{ 
int i,s; 
for (i=0;i<10;i++) 
    seed[i] = rand(); 

#pragma omp parallel private(s) 
{ 
    s = seed[omp_get_thread_num()]; 
    #pragma omp for 
    for (i=0;i<1000;i++) 
    { 
     printf("%d %d %d\n",i,omp_get_thread_num(),s); 
     s=(s*17931+7391); // those numbers should be choosen more carefully 
    } 
    seed[omp_get_thread_num()] = s; 
} 
return 0; 
} 

PD .: Todavía no he aceptado ninguna respuesta, porque necesito estar seguro de que esta idea es buena.

Answer 1

voy a publicar aquí lo que he publicado a Concurrent random number generation:

Creo que estás buscando rand_r(), que toma forma explícita el estado actual RNG como parámetro. Entonces cada hilo debe tener su propia copia de los datos iniciales (si quieres que cada hilo comience con la misma semilla o diferentes depende de lo que estés haciendo, aquí quieres que sean diferentes o obtendrás la misma fila una y otra vez). Hay una discusión sobre rand_r() y seguridad de hilos aquí: whether rand_r is real thread safe?.

Digamos que quería que cada hilo comenzara con su número de hilo (que probablemente no es el que desea, ya que daría los mismos resultados cada vez que ejecutara el mismo número de hilos, pero igual un ejemplo):

#pragma omp parallel default(none) 
{ 
    int i; 
    unsigned int myseed = omp_get_thread_num(); 
    #pragma omp for 
    for(i=0; i<100; i++) 
      printf("%d %d %d\n",i,omp_get_thread_num(),rand_r(&myseed)); 
} 

Editar: Sólo por intentarlo, comprueba para ver si lo anterior sería conseguir cualquier aumento de velocidad. código completo era

#define NRANDS 1000000 
int main(int argc, char **argv) { 

    struct timeval t; 
    int a[NRANDS]; 

    tick(&t); 
    #pragma omp parallel default(none) shared(a) 
    { 
     int i; 
     unsigned int myseed = omp_get_thread_num(); 
     #pragma omp for 
     for(i=0; i<NRANDS; i++) 
       a[i] = rand_r(&myseed); 
    } 
    double sum = 0.; 
    double time=tock(&t); 
    for (long int i=0; i<NRANDS; i++) { 
     sum += a[i]; 
    } 
    printf("Time = %lf, sum = %lf\n", time, sum); 

    return 0; 
} 

donde tic tac y son sólo envoltorios a gettimeofday() y tac() devuelve la diferencia en cuestión de segundos. La suma se imprime solo para asegurarse de que no se optimice nada, y para demostrar un pequeño punto; Obtendrás diferentes números con diferentes números de subprocesos porque cada subproceso obtiene su propio número de subproceso como una semilla; si ejecuta el mismo código una y otra vez con el mismo número de hilos obtendrá la misma suma, por el mismo motivo. De todos modos, el tiempo (que se ejecuta en una caja Nehalem de 8 núcleos con ningún otro usuario):

$ export OMP_NUM_THREADS=1 
$ ./rand 
Time = 0.008639, sum = 1074808568711883.000000 

$ export OMP_NUM_THREADS=2 
$ ./rand 
Time = 0.006274, sum = 1074093295878604.000000 

$ export OMP_NUM_THREADS=4 
$ ./rand 
Time = 0.005335, sum = 1073422298606608.000000 

$ export OMP_NUM_THREADS=8 
$ ./rand 
Time = 0.004163, sum = 1073971133482410.000000 

Así aumento de velocidad, si no es genial; como señala @ruslik, este no es realmente un proceso de cálculo intensivo, y otros problemas como el ancho de banda de la memoria comienzan a desempeñar un papel. Por lo tanto, solo un tono de más de 2 veces la aceleración en 8 núcleos.

Answer 2

Obtenga cada subproceso para establecer una semilla diferente en función de su Id. De subproceso, p. srand(omp_get_thread_num() * 1000);

Answer 3

Los números aleatorios se pueden generar muy rápido, por lo que generalmente la memoria sería el cuello de botella. Al dividir esta tarea entre varios hilos, crea una comunicación adicional y gastos indirectos de sincronización (y la sincronización de cachés de diferentes núcleos no es barata).

Sería mejor utilizar un solo hilo con una mejor función random().

Answer 4

Parece que rand tiene un estado global compartido entre todos los hilos en Linux y un estado de almacenamiento local de hilos para él en Windows. El estado compartido en Linux está causando su ralentización debido a la sincronización necesaria.

No creo que exista una forma portátil en la biblioteca C para usar el paralelo RNG en múltiples hilos, por lo que necesita otro. Puede usar un Mersenne Twister. Como Marcog dijo, necesitas inicializar la semilla para cada hilo de manera diferente.

Answer 5

No puede usar la función C rand() desde varios hilos; esto da como resultado un comportamiento indefinido. Algunas implementaciones pueden darle bloqueo (lo que lo hará lento); otros pueden permitir que los hilos toquen el estado de cada uno, posiblemente bloqueando su programa o simplemente dando números aleatorios "malos".

Para resolver el problema, escriba su propia implementación de PRNG o utilice una existente que permita al llamante almacenar y pasar el estado a la función de repetidor de PRNG.

Answer 6

en Linux/Unix puede utilizar

long jrand48(unsigned short xsubi[3]); 

donde xsubi [3] codifica el estado del generador de números aleatorios, como esto:

#include<stdio.h> 
#include<stdlib.h> 
#include <algorithm> 
int main() { 
    unsigned short *xsub; 
#pragma omp parallel private(xsub) 
    { 
    xsub = new unsigned short[3]; 
    xsub[0]=xsub[1]=xsub[2]= 3+omp_get_thread_num(); 
    int j; 
#pragma omp for 
    for(j=0;j<10;j++) 
     printf("%d [%d] %ld\n", j, omp_get_thread_num(), jrand48(xsub)); 
    } 
} 

compilar con

g++-mp-4.4 -Wall -Wextra -O2 -march=native -fopenmp -D_GLIBCXX_PARALLEL jrand.cc -o jrand 

(reemplace g ++ - mp-4.4 con lo que necesite para llamar a la versión 4.4 o 4.3 de g ++) y obtiene

$ ./jrand 
0 [0] 1344229389 
1 [0] 1845350537 
2 [0] 229759373 
3 [0] 1219688060 
4 [0] -553792943 
5 [1] 360650087 
6 [1] -404254894 
7 [1] 1678400333 
8 [1] 1373359290 
9 [1] 171280263 

es decir, 10 números pseudoaleatorios diferentes sin ningún tipo de bloqueo de mute o condiciones de carrera.