c - 从随机分布的粒子到规则网格的通信的最佳并行化

我正在并行化我用于模拟地球内部 2D 和 3D 变形的细胞内粒子代码。代码的几个例程很容易使用 OpenMP 并行化并且扩展性很好。但是，我在处理从粒子到网格单元的插值的代码的关键部分遇到了问题。粒子在每次迭代中都在移动(根据速度场)。许多计算在规则的非变形网格上执行最有效。因此，每次迭代都涉及从“随机”分布的粒子到网格单元的通信。

这个问题可以用以下简化的一维代码来说明:

//EXPLANATION OF VARIABLES (all previously allocated and initialized, 1D arrays)
//double *markerval; // Size Nm. Particle values. Are to be interpolated to the grid
//double *grid; // Size Ng=Nm/100 Grid values.
//uint *markerpos; // Size Nm. Position of particles relative to grid (each particle
// knows what grid cell it belongs to) possible values are 0,1,...Ng-1

//#pragma omp parallel for schedule(static) private(e)
for (e=0; e<Nm; e++) {
    //#pragma omp atomic
    grid[markerpos[e]]+=markerval[e];
}

在最坏的情况下，粒子位置是随机的，但通常，粒子在内存中彼此相邻，在空间中也彼此相邻，因此也在网格内存中。

我如何有效地并行化这个过程？多个粒子映射到同一个网格单元，因此如果上述循环直接并行化，则出现竞争条件和缓存交换的可能性非零。使更新原子化可以防止竞争条件，但会使代码比顺序情况慢得多。

我还尝试为每个线程制作网格值的私有(private)副本，然后随后将它们相加。但是，这可能需要在要使用的代码中使用太多内存，并且对于此示例，它在线程数量上的扩展性不佳(原因我不确定)。

第三种选择可能是从网格映射到粒子，然后遍历网格索引而不是粒子索引。然而，我担心这会非常复杂并且需要对代码进行重大更改，我不确定它会有多大帮助，因为它还需要使用排序例程，这在计算上也会很昂贵。

有没有人对此或类似问题有任何经验？

最佳答案

一种选择是在线程上手动映射迭代:

#pragma omp parallel shared(Nm,Ng,markerval,markerpos,grid)
{
  int nthreads = omp_get_num_threads();
  int rank     = omp_get_thread_num();
  int factor   = Ng/nthreads;

  for (int e = 0; e < Nm; e++) {
    int pos = markerpos[e];
    if ( (pos/factor)%nthreads == rank )
      grid[pos]+=markerval[e];
  }
}

几点说明:

for 循环的迭代不在线程之间共享。相反 每个线程 执行所有迭代。

for 循环内的条件决定 哪个线程 将更新 pos 数组的位置 grid。这个位置 只属于一个线程 ，因此不需要 atomic 保护。

公式 (pos/factor)%nthreads 只是一个简单的启发式方法。任何返回 pos 范围内值的 0,...,nthreads-1 函数实际上都可以替换为该表达式，而不会影响最终结果的有效性(因此，如果您有更好的尝试，请随意更改它)。请注意，此功能选择不当可能会导致负载平衡问题。

关于c - 从随机分布的粒子到规则网格的通信的最佳并行化，我们在Stack Overflow上找到一个类似的问题：https://stackoverflow.com/questions/13277505/