OpenMP 循环每次都会给出不同的结果

Question

#include <math.h>
#include <stdlib.h>
#include <stdio.h>
#include <omp.h>
#define KE 10
#define NSTEPS 100
main()
{
    float ex[KE],hy[KE];
    float bn[NSTEPS];
    int n,k,kc,ke;
    float ddx,dt,T;
    float t0,spread,pulse;
    omp_set_num_threads(4);
    kc = KE/2;
    t0 = 40.0;
    spread = 12;
    T = 0;
    FILE * temp = fopen("Exserial.txt", "w");
    /* Initialize */
    for ( k=0; k <= KE; k++ )
    {
        ex[k] = 0.0;
        hy[k] = 0.0;
    }
    double strt=omp_get_wtime();
#pragma omp parallel shared(ex,hy,bn) private(n,k,T,pulse,t0,spread)
    {
        for ( n=1; n<=NSTEPS ; n++) {
            T = T + 1;
            pulse = exp(-.5*(pow( (t0-T)/spread,2.0) ));
            ex[kc] = pulse;
    /* Main FDTD Loop */
    /* Calculate the Ex field */
#pragma omp critical
            {
#pragma omp for 
                for ( k=1; k < KE; k++ )
                { 
#pragma omp atomic update
                    ex[k] = ex[k] + .5*( hy[k-1] - hy[k] ) ;
                    bn[n]=ex[100];
                }
    /* Calculate the Hy field */
#pragma omp for
                for ( k=0; k < KE-1; k++ )
                {
#pragma omp atomic update
                    hy[k] = hy[k] + .5*( ex[k] - ex[k+1] ) ;
                }
            }
        }
    }
    double end=omp_get_wtime();
    printf( "------------------%f\n",end-strt);
    /* End of the Main FDTD Loop */
    /* At the end of the calculation, print out the Ex and Hy fields */
    for ( k=1; k <=KE; k++ )
    { 
        fprintf(temp, "%f\n",ex[k]);
        printf( "%f\n",ex[k]);
    }
}

并行输出（使用线程数=4）

 ------------------0.013850
    0.030408
    -0.130364
    0.107690
    0.061082
    0.023526
    0.066983
    -0.030821
    -0.024117
    0.037548
    0.000000

串行输出

-0.013813
-0.028405
-0.044722
-0.064258
-0.087550
-0.072226
-0.062273
-0.055638
-0.052620
0.000000

我不知道我在哪里犯了错误，但我在并行程序中得到了不同的结果。它每次都给出不同的输出。 还有一个问题，即当我从顶部删除 #pragma omp parallel 时，我得到了正确的值，但没有时间优势。我必须添加#pragma omp parallel 所以需要一些解决方案。

Answer 1

正如@1201ProgramAlarm 所建议的，您遇到了数据竞争。它涉及到并行区域内和临界区外的这个语句：

            ex[kc] = pulse;

每个线程在每次迭代时都会执行该值，并且每个线程也会读回该值。读取在关键部分内并不能避免数据竞争——写入也需要在关键部分（不一定是同一个），或者受某种同步保护，以避免数据竞赛。

但这里的问题要大得多，以至于数据竞赛几乎无关紧要。 main 问题是计算在循环迭代之间具有很强的数据依赖性。也就是说，当串行运行时，外循环的每次迭代都会根据前一次迭代中计算的值更新ex 和hy 的所有元素。无论您同步多少或您提出的关键原因有多广泛，如果您允许外部循环迭代以不同的顺序执行，您就不能期望相同的结果（甚至是一致的结果），就像并行化所做的那样。

这里的一线希望是并行化并没有让您获得太多好处，因为几乎所有外部循环执行的工作都在关键部分。仅删除 omp parallel 和 omp critical 指令将大大有助于改进此代码。对于如此小的KE 和如此简单的计算，离开内部的omp fors 可能也没有多大意义，但如果KE 更大那么你可以考虑离开内部循环的并行化。