2

我有一个程序，可以将粒子分配到单元格中的云网格中。简单地循环粒子总数 (Ntot) 并填充 256^3 网格（即每个粒子分布在 8 个单元上）。

% gfortran -fopenmp cic.f90 -o ./cic

编译得很好。但是当我运行它（./cic）时，我遇到了分段错误。我的循环是一个经典的 omp do 问题。当我不在 openmp 中编译它时，该程序可以工作。

!$omp parallel do
 do i = 1,Ntot
   if (x1(i).gt.0.and.y1(i).gt.0.and.z1(i).gt.0) then
     dense(int(x1(i)),int(y1(i)),int(z1(i))) = dense(int(x1(i)),int(y1(i)),int(z1(i))) &
     + dx1(i) * dy1(i) * dz1(i) * mpart
   end if

   if (x2(i).le.Ng.and.y1(i).gt.0.and.z1(i).gt.0) then
     dense(int(x2(i)),int(y1(i)),int(z1(i))) = dense(int(x2(i)),int(y1(i)),int(z1(i))) &
     + dx2(i) * dy1(i) * dz1(i) * mpart
   end if

   if (x1(i).gt.0.and.y2(i).le.Ng.and.z1(i).gt.0) then
     dense(int(x1(i)),int(y2(i)),int(z1(i))) = dense(int(x1(i)),int(y2(i)),int(z1(i))) &
     + dx1(i) * dy2(i) * dz1(i) * mpart
   end if

   if (x2(i).le.Ng.and.y2(i).le.Ng.and.z1(i).gt.0) then
     dense(int(x2(i)),int(y2(i)),int(z1(i))) = dense(int(x2(i)),int(y2(i)),int(z1(i))) &
     + dx2(i) * dy2(i) * dz1(i) * mpart
   end if

   if (x1(i).gt.0.and.y1(i).gt.0.and.z2(i).le.Ng) then
     dense(int(x1(i)),int(y1(i)),int(z2(i))) = dense(int(x1(i)),int(y1(i)),int(z2(i))) &
     + dx1(i) * dy1(i) * dz2(i) * mpart
   end if

   if (x2(i).le.Ng.and.y1(i).gt.0.and.z2(i).le.Ng) then
     dense(int(x2(i)),int(y1(i)),int(z2(i))) = dense(int(x2(i)),int(y1(i)),int(z2(i))) &
     + dx2(i) * dy1(i) * dz2(i) * mpart
   end if

   if (x1(i).gt.0.and.y2(i).le.Ng.and.z2(i).le.Ng) then
     dense(int(x1(i)),int(y2(i)),int(z2(i))) = dense(int(x1(i)),int(y2(i)),int(z2(i))) &
     + dx1(i) * dy2(i) * dz2(i) * mpart
   end if

   if (x2(i).le.Ng.and.y2(i).le.Ng.and.z2(i).le.Ng) then
     dense(int(x2(i)),int(y2(i)),int(z2(i))) = dense(int(x2(i)),int(y2(i)),int(z2(i))) &
     +  dx2(i) * dy2(i) * dz2(i) * mpart
   end if
  end do
!$omp end parallel do

迭代之间没有依赖关系。想法？

4

2 回答 2

2

这个问题以及您的另一个问题中的问题来自启用 OpenMP 时禁用自动堆数组这一事实。这意味着如果没有-fopenmp，大数组将自动放置在静态存储（称为.bss段）中，而小数组则分配在堆栈上。当您打开 OpenMP 支持时，不会使用自动静态分配，并且您的dense数组会在例程的堆栈上分配。OS X 上的默认堆栈限制非常严格，因此存在分段错误。

您在这里有几个选择。第一个选项是通过给它属性来进行dense静态分配。SAVE另一种选择是通过创建它ALLOCATABLE然后使用ALLOCATE语句在堆上显式分配它，例如：

REAL, DIMENSION(:,:,:), ALLOCATABLE :: dense

ALLOCATE(dense(256,256,256))

! Computations, computations, computations

DEALLOCATE(dense)

SAVE较新的 Fortran 版本支持在超出范围时自动释放不带属性的数组。

请注意，您的 OpenMP 指令很好，不需要额外的数据共享子句。您不需要i在PRIVATE子句中声明，因为循环计数器具有预先确定的私有数据共享类。您不需要将其他变量放在 SHARED子句中，因为它们是隐式共享的。然而，您执行的操作应该与（或仅在较旧的 OpenMP 实现上）dense同步，或者您应该使用. 原子更新被转换为锁定的添加，并且不会导致太多的减速，与在循环中使用条件的巨大减速相比：ATOMIC UPDATEATOMICREDUCTION(+:dense)

INTEGER :: xi, yi, zi

!$OMP PARALLEL DO PRIVATE(xi,yi,zi)
...
if (x1(i).gt.0.and.y1(i).gt.0.and.z1(i).gt.0) then
  xi = int(x1(i))
  yi = int(y1(i))
  zi = int(z1(i))
  !$OMP ATOMIC UPDATE
  dense(xi,yi,zi) = dense(xi,yi,zi) &
                  + dx1(i) * dy1(i) * dz1(i) * mpart
end if
...

复制代码并针对其他情况进行适当的更改。如果您的编译器抱怨构造UPDATE中的子句ATOMIC，只需将其删除即可。

REDUCTION(+:dense)将dense在每个线程中创建一个副本，这将消耗大量内存，并且最终应用的减少会随着dense. 对于小型数组，它会比原子更新更好。

于 2012-12-14T18:09:21.437 回答
-1

有关如何使变量共享和私有的描述，请参见https://computing.llnl.gov/tutorials/openMP/#Clauses 。

看起来您的所有变量都应该共享，除了i必须是私有的循环变量。这将建议使用以下行：

!$omp parallel do default(shared) private(i)

这应该可以解决您的分段错误（假设我得到了正确的所有变量

但是，存在不同线程会尝试同时覆盖相同部分的风险dense，从而导致总数不正确。为了防止这种情况，您需要将每个分配包装到dense类似!$omp atomic 或!$omp critical部分的内容中。

但是，您可能会发现这样的关键部分会导致线程花费大部分时间等待，因此您可能看不到与纯串行代码相比有任何改进。

原则上你可以通过声明关键字来解决这个问题dense，reduction但不幸的是它不能用于数组。

于 2012-12-14T02:33:47.890 回答