Параллелизация с openMP: общие и критические предложения

Question

Ниже приведена часть кода, распараллеленная через openMP. Массивы, ap[] и sc[], имеют отношение к добавлению, поэтому я решил сделать их разделяемыми, а затем поместить их в раздел статьи critical, так как статья reduction не принимает массивы. Но это дает другой результат, чем его последовательный аналог. В чем проблема?

Vector PN, Pf, Nf; // Vector is user-defined structure 
Vector NNp, PPp; 
Vector gradFu, gradFv, gradFw; 
float dynVis_eff, SGSf; 
float Xf_U, Xf_H; 
float mf_P, mf_N; 
float an_diff, an_conv_P, an_conv_N, an_trans; 
float sc_cd, sc_pres, sc_trans, sc_SGS, sc_conv_P, sc_conv_N; 
float ap_trans; 

#pragma omp parallel for 
for (int e=0; e<nElm; ++e) 
{ 
    ap[e] = 0.f; 
    sc[e] = 0.f; 
} 

#pragma omp parallel for shared(ap,sc) 
for (int f=0; f<nFaces; ++f) 
{ 
    PN = cntE[face_N[f]] - cntE[face_P[f]]; 
    Pf = cntF[f] - cntE[face_P[f]]; 
    Nf = cntF[f] - cntE[face_N[f]]; 
    PPp = Pf - (Pf|norm(PN))*norm(PN); 
    NNp = Nf - (Nf|norm(PN))*norm(PN); 

    mf_P = mf[f]; 
    mf_N = -mf[f]; 

    SGSf = (1.f-ifac[f]) * SGSvis[face_P[f]] 
     +  ifac[f] * SGSvis[face_N[f]]; 

    dynVis_eff = dynVis + SGSf; 

    an_diff = dynVis_eff * Ad[f]/mag(PN); 

    an_conv_P = -neg(mf_P); 
    an_conv_N = -neg(mf_N); 

    an_P[f] = an_diff + an_conv_P; 
    an_N[f] = an_diff + an_conv_N; 

    // cross-diffusion 
    sc_cd = an_diff * ((gradVel[face_N[f]]|NNp) - (gradVel[face_P[f]]|PPp)); 

    #pragma omp critical 
    {  
     ap[face_P[f]] += an_N[f]; 
     ap[face_N[f]] += an_P[f]; 

     sc[face_P[f]] += sc_cd + sc_conv_P; 
     sc[face_N[f]] += -sc_cd + sc_conv_N; 
    } 

Answer 1

Вы не указали, должны ли общие переменные в вашем параллельном предложении делиться или нет. Вы можете сделать это в целом с помощью предложения default. Если значение по умолчанию не указано, переменные все разделяются, что вызывает проблемы в вашем коде.

В вашем случае, я предполагаю, что вы должны пойти на

#pragma omp parallel for default(none), shared(ap,sc,face_N,face_P,cntF,cntE,mf,ifac,Ad,an_P,an_N,SGSvis,dynVis), private(PN,Pf,Nf,PPp,NNp,mf_P,mf_N,SGSf,dynVis_eff,an_diff,an_conv_P,an_conv_N,sc_cd) 

Я настоятельно рекомендую всегда использовать default(none) так, что компилятор жалуется каждый раз, когда вы не объявлять переменную явно и заставляет вас думать о том, это явно.