如何使用OpenMP正确地并行化这个数组？

在我尝试用openmp并行化代码后，数组中的元素是错误的，因为元素的顺序并不是很重要。还是用c++标准向量而不是数组来并行化更方便，能不能给我个简单的建议？

#include <stdio.h>
#include <math.h>

int main()
{
    int n = 100;
    int a[n*(n+1)/2]={0};
    int count=0;

    #pragma omp parallel for reduction(+:a,count)
    for (int i = 1; i <= n; i++) {
        for (int j = i + 1; j <= n; j++) {
            double k = sqrt(i * i + j * j);
            if (fabs(round(k) - k) < 1e-10) {
                a[count++] = i;
                a[count++] = j;
                a[count++] = (int) k;
            }
        }
    }
    
    for(int i=0;i<count;i++)
        printf("%d %s",a[i],(i+1)%3?"":", ");
    printf("
count: %d", count);
    return 0;
}

原始输出：

3 4 5、5 12 13、6 8 10、7 24 25、8 15 17、9 12 15、9 40 41、10 24 26、11 60 61、12 16 20、12 35 37、13 84 85、14 48 50、15 20 25、15 36 39、16 30 34、16 63 65、18 24 30、18 80 82、20 21 29、20 48 52、20 99 101、21 28 35、21 72 75、24 32 40、24 45 51、24 70、25 60 65、27 36 45、2832 60 68，33 44 55，33 56 65，35 84 91，36 48 60，36 77 85，39 52 65，39 80 89，40 42 58，40 75 85，40 96 104，42 56 70，45 60 75，48 55 73，48 64 80，48 90 102，51 68 85，54 72 90，56 90 106，57 76 95，60 63 87，60 80 100，60 91 109，63 84 105，65 72 97，66 88 110，69 92 115，72 96 120，75 100 125，80 84 116 计数：189

使用openmp(GCC文件c-fopenmp)后：

411 538 679,344 609 711,354 533 649,218 387 449,225 475 534,182 283339，81 161 182，74 190 204，77 138 159，79 176 195，18 24 30，18 80 82，0 0 0，0 0 0，0 0 0，0 0，0 0，0 0，0 0，0 0，0 0，0 0，0 0，0 0，0 0，0 0，0 0，0 0，0 0，0 0，0 0，0 00 0 0、0 0 0、0 0 0、0 0 0、0 0 0、0 0 0、0 0 0、0 0 0、0 0 0、0 0 0、0 0 0、0 0 0、0 0 0、0 0 0、0 0 0、0 0 0、0 0 0、0 0 0、0 0 0、0 0、0 0、0 0、0 0、0 0、0 0、0 0、0 0、0 0、0 0、0 0、0 0、0 0、0 0、0 0、0 0、0 0、0 0 计数：189

---

解决方案

作为使用critical节的替代方案，此解决方案使用原子，因此可以更快。

以下代码可能会由于内存消耗而冻结您的计算机。小心！

#include <cstdio>
#include <cmath>

#include <vector>

int main() {
    int const n = 100;
    // without a better (smaller) upper_bound this is extremely
    // wasteful in terms of memory for big n 
    long const upper_bound = 3L * static_cast<long>(n) *
                             (static_cast<long>(n) - 1L) / 2l; 
    std::vector<int> a(upper_bound, 0);
    int count = 0;

    #pragma omp parallel for schedule(dynamic) shared(a, count)
    for (int i = 1; i <= n; ++i) {
        for (int j = i + 1; j <= n; ++j) {
            double const k = std::sqrt(static_cast<double>(i * i + j * j));

            if (std::fabs(std::round(k) - k) < 1e-10) {
                int my_pos;
                #pragma omp atomic capture
                my_pos = count++;

                a[3 * my_pos] = i;
                a[3 * my_pos + 1] = j;
                a[3 * my_pos + 2] = static_cast<int>(std::round(k));
            }
        }
    }
    count *= 3;
    
    for(int i = 0; i < count; ++i) {
        std::printf("%d %s", a[i], (i + 1) % 3 ? "" : ", ");
    }
    printf("
count: %d", count);

    return 0;
}

编辑： 我的答案最初是对使用critical节次优方式删除的答案的反应。在下面，我将介绍另一个解决方案，它将critical部分与使用std::vector::emplace_back()相结合，以避免需要类似于Toby Speight的解决方案的upper_bound。通常，使用Toby Speight的解决方案中的reduce子句应该比使用critical节和atomics更可取，因为减少的线程数越多，伸缩性越好。在这种特定情况下(相对较少的计算将写入a)，并且没有大量内核可在其上运行，下面的代码可能仍然更可取。

#include <cstdio>
#include <cmath>

#include <tuple>
#include <vector>

int main() {
    int const n = 100;

    std::vector<std::tuple<int, int, int>> a{};
    
    // optional, might reduce number of reallocations
    a.reserve(2 * n); // 2 * n is an arbitrary choice

    #pragma omp parallel for schedule(dynamic) shared(a)
    for (int i = 1; i <= n; ++i) {
        for (int j = i + 1; j <= n; ++j) {
            double const k = std::sqrt(static_cast<double>(i * i + j * j));

            if (std::fabs(std::round(k) - k) < 1e-10) {
                #pragma omp critical
                a.emplace_back(i, j, static_cast<int>(std::round(k)));
            }
        }
    }
    long const count = 3L * static_cast<long>(a.size());
    
    for(unsigned long i = 0UL; i < a.size(); ++i) {
        std::printf("%d %d %d
",
                    std::get<0>(a[i]), std::get<1>(a[i]), std::get<2>(a[i]));
    }
    printf("
count: %ld", count);

    return 0;
}