在并行环境中调用子程序

Question

我认为我的问题与 here 描述的问题相关甚至相同。但是我不明白到底发生了什么。

我将 openMP 与 gfortran 编译器一起使用，我有以下任务要做：我在二维表面上有一个密度分布 F(X, Y)，x 坐标 X 和 y -坐标 Y。矩阵 F 的大小为 Nx x Ny.

我现在有一组坐标 Xp(i) 和 Yp(i)，我需要将密度 F 插值到这些点上。这个问题是为了并行化。

!$OMP PARALLEL DO DEFAULT(SHARED) PRIVATE(i)
    do i=1, Nmax

        ! Some stuff to be done here

        Fint(i) = interp2d(Xp(i), Yp(i), X, Y, F, Nx, Ny)

        ! Some other stuff to be done here

    end do
!$OMP END PARALLEL DO

除了 i 之外的所有内容都是共享的。函数 interp2d 正在做一些简单的线性插值。

单线程运行良好，但多线程运行失败。我将问题追溯到取自 Numerical Recipes 的 hunt-子例程，它被 interp2d 调用。 hunt-子例程基本上计算索引 ix 使得 X(ix) <= Xp(i) < X(ix+1)。这是获取插值起点所必需的。

对于多线程，它时不时会发生，一个线程从 hunt 获得正确的索引 ix，而调用 hunt 的线程下一个获得完全相同的索引，尽管 Xp(i) 还差得远。

我可以使用 CRITICAL 环境来防止这种情况发生：

!$OMP PARALLEL DO DEFAULT(SHARED) PRIVATE(i)
    do i=1, Nmax

        ! Some stuff to be done here

  !$OMP CRITICAL
        Fint(i) = interp2d(Xp(i), Yp(i), X, Y, F, Nx, Ny)
  !$OMP END CRITICAL

        ! Some other stuff to be done here

    end do
!$OMP END PARALLEL DO

但是这样会降低效率。例如，如果我使用三个线程，则 CRITICAL 环境的平均负载为 1.5。如果没有，我的平均负载为 2.75，但结果错误，甚至有时 SIGSEGV 运行时错误。

这里到底发生了什么？在我看来，所有线程都在调用相同的 hunt-子例程，如果它们同时调用，就会发生冲突。这有意义吗？

我该如何防止这种情况发生？

Answer 1

在 Fortran 90+ 中组合变量声明和初始化具有赋予变量 SAVE 属性的副作用。

integer :: i = 0

大致相当于：

integer, save :: i

if (first_invocation) then
  i = 0
end if

SAVE' 变量在例程的多次调用之间保留它们的值，因此通常作为静态变量实现。根据管理 OpenMP 中隐式数据共享类的规则，除非在 threadprivate 指令中列出，否则此类变量将被共享。

OpenMP 要求兼容的编译器应应用上述语义，即使底层语言是 Fortran 77。

在并行环境中调用子程序

Calling subroutine in parallel environment

fortran

openmp

gfortran

fortran90