我怎样才能避免这个循环中的“i”依赖？福特兰语言

Question

我在代码中使用 OpenMP，但为了做到这一点，我必须解决此依赖性：

do q=1,pppp
    i=0
    
    DO j=1, pppp
        do c1=1,3
            vect(c1)=xx(q,c1)-xx(j,c1)
        end do
        dist=sqrt(vect(1)**2+vect(2)**2+vect(3)**2)
        if(dist<0.0001)then
            i=i+1
            if(i>10)i=10
            caravec(q,i)=j
        endif
    ENDDO
ENDDO  

我试图避免使用该ordered条款，因为它很昂贵，但我不知道如何消除依赖性。我怎样才能做到这一点？
感谢所有的帮助

Answer 1

免责声明：我从未使用 Fortan 编程。所以我的代码很可能有错误。

但是使用逻辑并假设您所说的是正确的，即唯一的依赖项是“i”变量。您可以尝试以下操作：

do q=1,pppp        
    DO j=1, pppp
        do c1=1,3
            vect(c1)=xx(q,c1)-xx(j,c1)
        end do
        dist(q,j)=sqrt(vect(1)**2+vect(2)**2+vect(3)**2)
    ENDDO
ENDDO  

因此，您可以使用 OpenMP 并行化外循环。这样，您就可以并行运行循环中计算要求最高的部分，而无需依赖“i”。请注意，我将距离结果保存到向量中。

然后你可以依次进行剩余的循环计算，如下所示：

do q=1,pppp
    i=0    
    DO j=1, pppp
        if(dist(q,j) < 0.0001)then
            i=i+1
            if(i>10)i=10
            caravec(q,i)=j
        endif
    ENDDO
ENDDO

当前方法的缺点是，除了需要重构代码之外，它还比原始方法使用更多的内存。另一方面，它消除了线程之间同步的需要。

正如@Ian Bush在评论中提到的，您可以通过以下方式进一步改进您的顺序代码：

测试 n 距离平方而不是距离将避免调用 sqrt

Answer 2

您可以采取一些 Fortran 措施来改进@dreamcrash的优秀答案：

• Fortran 是列优先的，因此您希望迭代内部循环中最左边的索引，而不是最右边的索引。所以你可能应该转置xx，所以它是一个3*pppp索引为的数组，xx(c1,q)而不是pppp*3索引为的数组xx(q,c1)。

• 您可能希望使用整个数组操作而不是单元素操作，因为它们更有可能被矢量化。

• 您可以存储结果dist<0.0001而不是dist.

• 您可以将 替换if (i>10)...为min语句，该语句将是无分支的，因此可能运行得更快。

所以代码看起来像这样：

do q=1,pppp        
  do j=1,pppp
    vect = xx(:,q)-xx(:,j)
    within_tolerance(j, q) = dot_product(vect, vect) < 1.0e-8
  enddo
enddo

do q=1,pppp
  i=0    
  do j=1,pppp
    if (within_tolerance(j, q)) then
      i = min(i+1, 10)
      caravec(q,i)=j
    endif
  enddo
enddo

您可能还想转置caravec，但这取决于它在其他地方的使用方式。

XY 问题的解决方案

如果这是代码的瓶颈，您可能需要研究基于体素的方法来查找向量集中的最近邻。快速谷歌会带来例如这个。

q体素方法消除了对和的双循环的需要j，这可以让您更快地进行此类比较。

另一方面，体素方法非常复杂，我认为它们并不适合所有情况（例如，如果您有非常密集的点区域和非常稀疏的点区域，那么我认为体素方法会很困难）。