与 omp stucks 并行

Question

我对以下代码有疑问：

int *chosen_pts = new int[k];
std::pair<float, int> *dist2 = new std::pair<float, int>[x.n];
// initialize dist2
for (int i = 0; i < x.n; ++i) {
    dist2[i].first = std::numeric_limits<float>::max();
    dist2[i].second = i;
}

// choose the first point randomly
int ndx = 1;
chosen_pts[ndx - 1] = rand() % x.n;
double begin, end;
double elapsed_secs;
while (ndx < k) {
    float sum_distribution = 0.0;
    // look for the point that is furthest from any center
    begin = omp_get_wtime();
    #pragma omp parallel for reduction(+:sum_distribution)
    for (int i = 0; i < x.n; ++i) {

        int example = dist2[i].second;
        float d2 = 0.0, diff;
        for (int j = 0; j < x.d; ++j) {
            diff = x(example,j) - x(chosen_pts[ndx - 1],j);
            d2 += diff * diff;
        }
        if (d2 < dist2[i].first) {
            dist2[i].first = d2;
        }

        sum_distribution += dist2[i].first;

    }

    end = omp_get_wtime() - begin;

    std::cout << "center assigning -- " 
            << ndx << " of " << k << " = " 
            << (float)ndx / k * 100 
            << "% is done. Elasped time: "<< (float)end <<"\n";        

    /**/
    bool unique = true;

    do {
        // choose a random interval according to the new distribution
        float r = sum_distribution * (float)rand() / (float)RAND_MAX;
        float sum_cdf = dist2[0].first;
        int cdf_ndx = 0;
        while (sum_cdf < r) {
            sum_cdf += dist2[++cdf_ndx].first;
        }
        chosen_pts[ndx] = cdf_ndx;

        for (int i = 0; i < ndx; ++i) {
            unique = unique && (chosen_pts[ndx] != chosen_pts[i]);
        }
    } while (! unique);


    ++ndx;
}

如您所见，我使用 omp 使 for 循环并行化。它工作正常，我可以实现显着的加速。但是，如果我将 x.n 的值增加到 20000000 以上，该函数将在 8-10 次循环后停止工作：

它不产生任何输出 (std::cout)
只有一个核心工作
没有任何错误

如果我注释掉 do while 循环，它会按预期再次运行。所有核心都很忙，每次迭代后都有一个输出，我可以根据需要增加 k.n 超过 1 亿。

Answer 1

卡住的不是 OpenMP 并行，显然是在串行 do-while 循环中。

我看到的一个特殊问题是在访问 dist2 的内部 while 循环中没有数组边界检查。理论上，越界访问永远不应该发生；但在实践中它可能 - 请参阅下面的原因。所以首先我会重写 cdf_ndx 的计算以保证在检查所有元素时循环结束：

    float sum_cdf = 0;
    int cdf_ndx = 0;
    while (sum_cdf < r && cdf_ndx < x.n ) {
        sum_cdf += dist2[cdf_ndx].first;
        ++cdf_ndx;
    }

现在，sum_cdf 怎么可能达不到 r？这是由于浮点运算的特殊性以及 sum_distribution 是并行计算而 sum_cdf 是串行计算的事实。问题是一个元素对总和的贡献可能低于浮点数的精度；换句话说，当您将相差超过 ~8 个数量级的两个浮点值相加时，较小的值不会影响总和。

因此，在某个点之后有 20M 的浮点数时，可能会发生下一个要添加的值与累积的 sum_cdf 相比太小以至于添加这个值不会改变它！另一方面，sum_distribution 本质上是作为几个独立的部分和（每个线程一个）计算的，然后组合在一起。因此它更准确，并且可能比 sum_cdf 所能达到的更大。

一种解决方案是分部分计算 sum_cdf，具有两个嵌套循环。例如：

    float sum_cdf = 0;
    int cdf_ndx = 0;
    while (sum_cdf < r && cdf_ndx < x.n ) {
        float block_sum = 0;
        int block_end = min(cdf_ndx+10000, x.n); // 10000 is arbitrary selected block size
        for (int i=cdf_ndx; i<block_end; ++i ) {
            block_sum += dist2[i].first;
            if( sum_cdf+block_sum >=r ) {
                block_end = i; // adjust to correctly compute cdf_ndx
                break;
            }
        }
        sum_cdf += block_sum;
        cdf_ndx = block_end;
    }

在循环之后你需要检查 cdf_ndx < x.n，否则用新的随机间隔重复。

与 omp stucks 并行

parallel for with omp stucks

c++

openmp