OpenMP e GSL RNG - Problema di prestazioni - implementazione di 4 thread 10 volte più lenta rispetto a quella sequenziale pura (CPU quadcore)

Question

Sto cercando di trasformare un mio progetto C da sequenziale in programmazione parallela. Sebbene la maggior parte del codice sia stata ridisegnata da zero per questo scopo, la generazione di numeri casuali è ancora al centro. Pertanto, le cattive prestazioni del generatore di numeri casuali (RNG) influiscono molto negativamente sulle prestazioni generali del programma.

Ho scritto alcune linee di codice (vedi sotto) per mostrare il problema che sto affrontando senza molta verbosità.

Il problema è il seguente: ogni volta che aumenta il numero di thread nt, la performance peggiora notevolmente. A questa workstation (kernel Linux 2.6.33.4; gcc 4.4.4; CPU quadcore intel) il ciclo parallelo richiede circa 10 volte più tempo per terminare con nt = 4 che con nt = 1, indipendentemente dal numero di iterazioni n.

Questa situazione sembra essere descritta here ma il focus è principalmente in fortran, una lingua di cui so molto poco, quindi apprezzerei molto l'aiuto.

Ho cercato di seguire la loro idea di creare diversi RNG (con un seme diverso) a cui accedere da ciascun thread, ma le prestazioni sono ancora pessime. In realtà, questo diverso punto di seeding per ogni thread mi infastidisce anche perché non riesco a vedere come sia possibile garantire la qualità dei numeri generati alla fine (mancanza di correlazioni, ecc.).

Ho già pensato di eliminare completamente il GSL e di implementare un algoritmo generatore casuale (come Mersenne-Twister) ma sospetto che mi piacerebbe imbattersi nello stesso problema più avanti.

Grazie mille in anticipo per le vostre risposte e consigli. Per favore, chiedi qualcosa di importante che potrei aver dimenticato di menzionare.

EDIT: correzioni suggerite da Implemented lucas1024 (pragma per-loop di dichiarazione) e JonathanDursi (semina, di impostazione "a" come una variabile privata). Le prestazioni sono ancora molto lente nella modalità multithread.

MODIFICA 2: Soluzione implementata suggerita da Jonathan Dursi (vedere i commenti).

#include <stdio.h> 
    #include <stdlib.h> 
    #include <time.h> 
    #include <gsl/gsl_rng.h> 
    #include <omp.h> 

    double d_t (struct timespec t1, struct timespec t2){ 

     return (t2.tv_sec-t1.tv_sec)+(double)(t2.tv_nsec-t1.tv_nsec)/1000000000.0; 
    } 

    int main (int argc, char *argv[]){ 

     double a, b; 

     int i,j,k; 

     int n=atoi(argv[1]), seed=atoi(argv[2]), nt=atoi(argv[3]); 

     printf("\nn\t= %d", n); 
     printf("\nseed\t= %d", seed); 
     printf("\nnt\t= %d", nt); 

     struct timespec t1, t2, t3, t4; 

     clock_gettime(CLOCK_PROCESS_CPUTIME_ID, &t1); 

     //initialize gsl random number generator 
     const gsl_rng_type *rng_t; 
     gsl_rng **rng; 
     gsl_rng_env_setup(); 
     rng_t = gsl_rng_default; 

     rng = (gsl_rng **) malloc(nt * sizeof(gsl_rng *)); 

      #pragma omp parallel for num_threads(nt) 
     for(i=0;i<nt;i++){ 
      rng[i] = gsl_rng_alloc (rng_t); 
      gsl_rng_set(rng[i],seed*i); 
     } 

     clock_gettime(CLOCK_PROCESS_CPUTIME_ID, &t2); 

     for (i=0;i<n;i++){ 
      a = gsl_rng_uniform(rng[0]); 
     } 

     clock_gettime(CLOCK_PROCESS_CPUTIME_ID, &t3); 

     omp_set_num_threads(nt); 
     #pragma omp parallel private(j,a) 
     { 
      j = omp_get_thread_num(); 
      #pragma omp for 
      for(i=0;i<n;i++){ 
       a = gsl_rng_uniform(rng[j]); 
      } 
     } 

     clock_gettime(CLOCK_PROCESS_CPUTIME_ID, &t4); 

     printf("\n\ninitializing:\t\tt1 = %f seconds", d_t(t1,t2)); 
     printf("\nsequencial for loop:\tt2 = %f seconds", d_t(t2,t3)); 
     printf("\nparalel for loop:\tt3 = %f seconds (%f * t2)", d_t(t3,t4), (double)d_t(t3,t4)/(double)d_t(t2,t3)); 
     printf("\nnumber of threads:\tnt = %d\n", nt); 

     //free random number generator 
     for (i=0;i<nt;i++) 
      gsl_rng_free(rng[i]); 
     free(rng); 

     return 0; 

    } 

Answer 1

Il problema è nella seconda riga #pragma omp. Il primo #pragma omp genera 4 thread. Dopodiché dovresti semplicemente dire #pragma omp for - non #pragma omp parallel for.

Con il codice corrente, in base alle impostazioni di nidificazione omp, si creano 4 x 4 thread che eseguono lo stesso lavoro e accedono agli stessi dati.