线程总结变量在OpenMP中给出了错误的答案

Question

为了练习并行化do循环,我在Fortran中执行以下操作

$\integral{0}{1} \frac{4}{1+x^{2}} = \pi$

以下是我实现的代码:

      program mpintegrate

      integer i,nmax,nthreads,OMP_GET_NUM_THREADS
      real xn,dx,value
      real X(100000)

      nthreads = 4
      nmax = 100000
      xn = 0.0
      dx = 1.0/nmax
      value = 0.0

      do i=1,nmax
         X(i) = xn
         xn = xn + dx
      enddo

      call OMP_SET_NUM_THREADS(nthreads)

!$OMP Parallel 

!$OMP Do Schedule(Static) Private(i,X)

      do i=1,nmax
         value = value + dx*(4.0/(1+X(i)*X(i)))
      enddo

!$OMP End DO NoWait

!$OMP End Parallel

      print *, value

      end

编译程序没问题

gfortran -fopenmp -o mpintegrate mpintegrate.f

问题是当我执行程序时.当我按原样运行程序时,我得到的值(1,4).然而,当我使用omp do循环取消注释print语句时,最终值大约是应该是什么,pi.

为什么答案value不正确？

Answer 1

这里的一个问题是X需要不是私有的(需要在并行线上指定,而不是do线); 每个人都需要看到它,并且每个线程都有单独的副本没有意义.更糟糕的是,您在此处访问私有副本所获得的结果是未定义的,因为一旦您进入私有区域,该私有变量尚未初始化.您可以使用firstprivate而不是private,使用并行区域之前的内容为您初始化它,但这里最简单/最好的只是shared.

最终没有多大意义no wait,因为end parallel必须等待每个人都完成.

然而,话虽如此,你仍然有一个非常重要(和经典)的正确性问题.如果你在循环中更明确一点,那么这里发生的事情会更清楚(因为问题不依赖于所选的时间表而放弃了清晰度):

!$OMP Parallel do Private(i) Default(none) Shared(value,X,dx,nmax)

      do i=1,nmax
         value = value + dx*(4.0/(1+X(i)*X(i)))
      enddo

!$OMP End Parallel Do
print *, value

重复运行这会给出不同的值:

$ ./foo
   1.6643878
$ ./foo
   1.5004054
$ ./foo
   1.2746993

问题是所有线程都写入相同的共享变量value.这是错误的-每个人都在一次写,结果是乱码,作为一个线程可以计算出它自己的贡献,准备将其添加到value,只是因为它是约而至,另一个线程可以做它的写入value,然后把它迅速重挫.对同一共享变量的并发写入是一种典型的竞争条件,这是一种标准的错误系列,特别是在与OpenMP一样的共享内存编程中.

除了错误之外,它还很慢.争用相同的几个字节的内存的许多线程 - 内存足够接近以便落入同一缓存行 - 由于内存系统中的争用而非常慢.即使它们不是完全相同的变量(在这种情况下也是如此),这种内存争用 - 在它们恰好是相邻变量的情况下的假共享 - 可以显着减慢速度.取出显式线程数设置,并使用环境变量:

$ export OMP_NUM_THREADS=1
$ time ./foo
   3.1407621

real    0m0.003s
user    0m0.001s
sys 0m0.001s

$ export OMP_NUM_THREADS=2
$ time ./foo
   3.1224852

real    0m0.007s
user    0m0.012s
sys 0m0.000s

$ export OMP_NUM_THREADS=8
$ time ./foo
   1.1651508

real    0m0.008s
user    0m0.042s
sys 0m0.000s

因此,使用更多线程运行的速度几乎要慢3倍(而且越来越多).

那么我们可以做些什么来解决这个问题？我们可以做的一件事就是确保每个人的补充都不会被atomic指责:

!$OMP Parallel do Schedule(Static) Private(i) Default(none) Shared(X,dx, value, nmax)
      do i=1,nmax
!$OMP atomic
         value = value + dx*(4.0/(1+X(i)*X(i)))
      enddo
!$OMP end parallel do

这解决了正确性问题:

$ export OMP_NUM_THREADS=8
$ ./foo
   3.1407621

但对速度问题没有任何作用:

$ export OMP_NUM_THREADS=1
$ time ./foo
   3.1407621

real    0m0.004s
user    0m0.001s
sys 0m0.002s

$ export OMP_NUM_THREADS=2
$ time ./foo
   3.1407738

real    0m0.014s
user    0m0.023s
sys 0m0.001s

(注意,对于不同的线程数,你会得到略有不同的答案.这是由于最终总和的计算顺序与串行情况不同.对于单精度实数,由于操作顺序不同,差异显示在第7位很难避免,这里我们做了10万次操作.)

那么我们还能做些什么呢？一种方法是让每个人跟踪他们自己的部分金额,然后在我们完成时将它们全部加起来:

!... 

      integer, parameter :: nthreads = 4
      integer, parameter :: space=8
      integer :: threadno
      real, dimension(nthreads*space) :: partials
!...

      partials=0

!...

!$OMP Parallel Private(value,i,threadno) Default(none) Shared(X,dx, partials)
      value = 0
      threadno = omp_get_thread_num()
!$OMP DO
      do i=1,nmax
         value = value + dx*(4.0/(1+X(i)*X(i)))
      enddo
!$OMP END DO
      partials((threadno+1)*space) = value
!$OMP end parallel

      value = sum(partials)
      print *, value
      end

这是有效的 - 我们得到了正确的答案,如果你玩线程的数量,你会发现它非常活泼 - 我们将部分和数组中的条目隔开以避免错误共享(并且它是错误的,这时间,因为每个人都在写入阵列中的不同条目 - 没有覆盖.

尽管如此,这是一项愚蠢的工作,只是为了在线程之间获得正确的总和!有一种更简单的方法可以做到这一点--OpenMP有一个自动执行此操作的减少构造(并且比上面的手工版本更有效:)

!$OMP Parallel do reduction(+:value) Private(i) Default(none) Shared(X,dx)
      do i=1,nmax
         value = value + dx*(4.0/(1+X(i)*X(i)))
      enddo
!$OMP end parallel do
      print *, value 

现在程序运行正常,速度快,代码相当简单.在更现代的Fortran中,最终代码看起来像这样:

      program mpintegrate
      use omp_lib
      integer, parameter :: nmax = 100000
      real :: xn,dx,value
      real :: X(nmax)
      integer :: i

      integer, parameter :: nthreads = 4

      xn = 0.0
      dx = 1.0/nmax
      value = 0.0
      partials=0

      do i=1,nmax
         X(i) = xn
         xn = xn + dx
      enddo

      call omp_set_num_threads(nthreads)

!$OMP Parallel do reduction(+:value) Private(i) Default(none) Shared(X,dx)
      do i=1,nmax
         value = value + dx*(4.0/(1+X(i)*X(i)))
      enddo
!$OMP end parallel do
      print *, value

      end