2017-08-04 5 views
0

Je suis nouveau à OpenMP en C. Je l'ai utilisé pour mettre en parallèle mes boucles for dans une fonction, mais il s'avère considérablement ralenti mes boucles for en comparaison avec un seul thread Cas. Par exemple, la boucle for peut être effectuée autour de 10s pour chaque point (halo) mais cela prend quelques minutes avec OpenMP.C OpenMP parallèle pour les boucles font beaucoup plus lent que le simple thread

Dans cette fonction, j'essaie de calculer la densité de plusieurs coques pour chaque point (halo), en comptant les particules qui se trouvent à l'intérieur de la coque, puis les écarter en un réseau. Il y a 512^3 particules, et environ 200 points (halos) que je veux calculer. Je veux diviser les points (halos) pour différents threads pour le rendre plus rapide.

#define ArrayAccess2D_n2(a, n1, n2, i1, i2) (a)[ i2+n2*i1 ] 


void halo_shell_rho(float boxsize, float *halo_pos, float *halo_R, int halo_number,\ 
int halo_start, int halo_end, float *par_pos, long long par_number,\ 
int shell_bins, float rmax_fac, float *out_shell_den){ 

    float temp; 

    long long iter_sfs=0, iter_sfc=0, iter_ufs=0, iter_ufc=0; 
    int dim=3; 

    float par_posx, par_posy, par_posz, dist; 
    float halo_posx, halo_posy, halo_posz, halo_rad; 
    int i=0, ini_j=0, vol_j=0; 
    int a=0, b=0; 
    long long k=0; 

    #pragma omp parallel for private(i, ini_j, vol_j, a, b, k) 
    for(i=halo_start; i<=halo_end; i++){ 
      printf("halo %d\n", i); 
      float count[shell_bins]; 
      float volume[shell_bins]; 

      for(ini_j=0; ini_j<shell_bins; ini_j++){ 
        count[ini_j] = 0; 
        volume[ini_j] = 0; } 

      halo_posx = ArrayAccess2D_n2(halo_pos, dim, halo_number, 0, i); 
      halo_posy = ArrayAccess2D_n2(halo_pos, dim, halo_number, 1, i); 
      halo_posz = ArrayAccess2D_n2(halo_pos, dim, halo_number, 2, i); 
      halo_rad = halo_R[i]; 

      for(vol_j=0; vol_j<shell_bins; vol_j++){ 

        volume[vol_j] = shell_volume((vol_j+1)*halo_rad*rmax_fac/(shell_bins*1000), vol_j*halo_rad*rmax_fac/(shell_bins*1000)); } 

      for(k=0; k<par_number; k++){ 

        par_posx = ArrayAccess2D_n2(par_pos, par_number, dim, k, 0); 
        par_posy = ArrayAccess2D_n2(par_pos, par_number, dim, k, 1); 
        par_posz = ArrayAccess2D_n2(par_pos, par_number, dim, k, 2); 

        dist = pb_distance(boxsize*1000, halo_posx, halo_posy, halo_posz, par_posx, par_posy, par_posz); //1000 for boxsize in Mpc 

        if(dist <= 2*rmax_fac*halo_rad){ 

          for(a=0; a<shell_bins; a++){ 

            if((dist <= halo_rad*(a+1)*rmax_fac/shell_bins) && (dist >= halo_rad*a*rmax_fac/shell_bins)){ 

              count[a] += 1; } 
          } 
        } 
      } 

      for(b=0; b<shell_bins; b++){ 

      out_shell_den[(i-halo_start+0*(1+halo_end-halo_start))*shell_bins+b] = count[b]/volume[b]; 
      //out_shell_den has shape (2, halo_number, shell_bins), 0 for edge, 1 for density 
      out_shell_den[(i-halo_start+1*(1+halo_end-halo_start))*shell_bins+b] = (2*b+1)*rmax_fac/(shell_bins*2); 
      } 
    } 

}

Quelqu'un pourrait-il me aider avec ça? Je sais que c'est une question très fréquente qui est posée, mais je n'ai pas trouvé de solutions d'autres postes. Je cours sur un cluster avec 32 threads si cela aide.

Merci!

+0

Les threads ne vont-ils pas se battre sur 'halo_rad'? –

+0

Oui, vous devez regarder la portée des variables et corriger les conditions de course. – tim18

+0

@DavidSchwartz Ce qui signifie que halo_posx, halo_posy, halo_posz sont également combattus? –

Répondre

0

Merci pour @DavidSchwartz et @ tim18.

Les variables comme halo_rad et par_posx sont déclarées avant le parallèle, ce qui signifie qu'elles sont implicitement supposées être publiques. Donc, il ralentit parce que tous les threads se disputent le droit de les utiliser. Une façon de résoudre ceci est d'ajouter toutes les variables à private(). Mais je pense que le meilleur moyen est de simplement déclarer les variables à l'intérieur du parallèle comme ceci:

void halo_shell_rho(float boxsize, float *halo_pos, float *halo_R, int halo_number, int halo_start, int halo_end, float *par_pos, long long par_number, int shell_bins, float rmax_fac, float *out_shell_den){ 

    int dim=3; 
    int i=0, ini_j=0, vol_j=0, a=0, b=0; 
    long long k=0; 

    #pragma omp parallel for private(i, ini_j, vol_j, a, b, k) 
    for(i=halo_start; i<=halo_end; i++){ 
      printf("halo %d\n", i); 

      float halo_posx, halo_posy, halo_posz, halo_rad; 
      float count[shell_bins]; 
      float volume[shell_bins]; 

      for(ini_j=0; ini_j<shell_bins; ini_j++){ 
        count[ini_j] = 0; 
        volume[ini_j] = 0; } 

      halo_posx = ArrayAccess2D_n2(halo_pos, dim, halo_number, 0, i); 
      halo_posy = ArrayAccess2D_n2(halo_pos, dim, halo_number, 1, i); 
      halo_posz = ArrayAccess2D_n2(halo_pos, dim, halo_number, 2, i); 
      halo_rad = halo_R[i]; 

      for(vol_j=0; vol_j<shell_bins; vol_j++){ 

        volume[vol_j] = shell_volume((vol_j+1)*halo_rad*rmax_fac/(shell_bins*1000), vol_j*halo_rad*rmax_fac/(shell_bins*1000)); } 


      for(k=0; k<par_number; k++){ 
        float par_posx, par_posy, par_posz, dist; 

        par_posx = ArrayAccess2D_n2(par_pos, par_number, dim, k, 0); 
        par_posy = ArrayAccess2D_n2(par_pos, par_number, dim, k, 1); 
        par_posz = ArrayAccess2D_n2(par_pos, par_number, dim, k, 2); 

        dist = pb_distance(boxsize*1000, halo_posx, halo_posy, halo_posz, par_posx, par_posy, par_posz); //1000 for boxsize in Mpc 

        if(dist <= 2*rmax_fac*halo_rad){ 

          for(a=0; a<shell_bins; a++){ 

            if((dist <= halo_rad*(a+1)*rmax_fac/shell_bins) && (dist >= halo_rad*a*rmax_fac/shell_bins)){ 

              count[a] += 1; } 
          } 
        } 
      } 

      for(b=0; b<shell_bins; b++){ 

      out_shell_den[(i-halo_start+0*(1+halo_end-halo_start))*shell_bins+b] = count[b]/volume[b]; //out_shell_den has shape (2, halo_number, shell_bins), 0 for edge, 1 for density 
      out_shell_den[(i-halo_start+1*(1+halo_end-halo_start))*shell_bins+b] = (2*b+1)*rmax_fac/(shell_bins*2); 
      } 
    } 
}