multithreading - 单/多线程 (OpenMP) 模式下计算精度的差异

Question

有人可以解释/理解单/多线程模式下计算结果的不同吗？

这是一个大约的例子。圆周率的计算:

#include <iomanip>
#include <cmath>
#include <ppl.h>

const int itera(1000000000);

int main()
{
    printf("PI calculation \nconst int itera = 1000000000\n\n");

    clock_t start, stop;

    //Single thread
    start = clock();
    double summ_single(0);
    for (int n = 1; n < itera; n++)
    {
        summ_single += 6.0 / (static_cast<double>(n)* static_cast<double>(n));
    };
    stop = clock();
    printf("Time single thread             %f\n", (double)(stop - start) / 1000.0);


    //Multithread with OMP
    //Activate OMP in Project settings, C++, Language
    start = clock();
    double summ_omp(0);
#pragma omp parallel for reduction(+:summ_omp)
    for (int n = 1; n < itera; n++)
    {
        summ_omp += 6.0 / (static_cast<double>(n)* static_cast<double>(n));
    };
    stop = clock();
    printf("Time OMP parallel              %f\n", (double)(stop - start) / 1000.0);


    //Multithread with Concurrency::parallel_for
    start = clock();
    Concurrency::combinable<double> piParts;
    Concurrency::parallel_for(1, itera, [&piParts](int n)
    {
        piParts.local() += 6.0 / (static_cast<double>(n)* static_cast<double>(n)); 
    }); 

    double summ_Conparall(0);
    piParts.combine_each([&summ_Conparall](double locali)
    {
        summ_Conparall += locali;
    });
    stop = clock();
    printf("Time Concurrency::parallel_for %f\n", (double)(stop - start) / 1000.0);

    printf("\n");
    printf("pi single = %15.12f\n", std::sqrt(summ_single));
    printf("pi omp    = %15.12f\n", std::sqrt(summ_omp));
    printf("pi comb   = %15.12f\n", std::sqrt(summ_Conparall));
    printf("\n");

    system("PAUSE");

}

结果:

PI calculation VS2010 Win32
Time single thread 5.330000
Time OMP parallel 1.029000
Time Concurrency:arallel_for 11.103000

pi single = 3.141592643651
pi omp = 3.141592648425
pi comb = 3.141592651497


PI calculation VS2013 Win32
Time single thread 5.200000
Time OMP parallel 1.291000
Time Concurrency:arallel_for 7.413000

pi single = 3.141592643651
pi omp = 3.141592648425
pi comb = 3.141592647841


PI calculation VS2010 x64
Time single thread 5.190000
Time OMP parallel 1.036000
Time Concurrency::parallel_for 7.120000

pi single = 3.141592643651
pi omp = 3.141592648425
pi comb = 3.141592649319


PI calculation VS2013 x64
Time single thread 5.230000
Time OMP parallel 1.029000
Time Concurrency::parallel_for 5.326000

pi single = 3.141592643651
pi omp = 3.141592648425
pi comb = 3.141592648489

测试是在 AMD 和 Intel CPU、Win 7 x64 上进行的。

单核和多核PI计算不同的原因是什么？
为什么 Concurrency::parallel_for 的计算结果在不同版本(编译器、32/64 位平台)上不是恒定的？

附言
Visual Studio Express 不支持 OpenMP。

Answer 1

由于舍入错误，浮点加法是非关联运算，因此运算顺序很重要。让您的并行程序给出与串行版本不同的结果是正常的。理解和处理它是编写(可移植)并行代码艺术的一部分。这在 32 位和 64 位版本中更加严重，因为在 32 位模式下，VS 编译器使用 x87 指令，而 x87 FPU 以 80 位的内部精度执行所有操作。在 64 位模式下使用 SSE 数学。

在串行情况下，一个线程计算 s1+s2+...+sN，其中 N 是展开中的项数。

在 OpenMP 情况下，有 n 个部分和，其中 n 是 OpenMP 线程的数量。哪些项进入每个部分和取决于迭代在线程之间的分布方式。许多 OpenMP 实现的默认设置是静态调度，这意味着线程 0(主线程)计算 ps0 = s1 + s2 + ... + sN/n；线程 1 计算 ps1 = sN/n+1 + sN/n+2 + ... + s2N/n；等等。最后，减少以某种方式结合了这些部分总和。
parallel_for情况与 OpenMP 非常相似。不同之处在于默认情况下迭代以动态方式分布 - 请参阅 auto_partitioner 的文档。 ，因此每个部分和包含或多或少随机选择的项。这不仅给出了略有不同的结果，而且每次执行也给出了略有不同的结果，即来自两个连续 parallel_for 的结果。具有相同线程数的 's 可能会有所不同。如果您用 simple_partitioner 的实例替换分区器并将块大小设置为 itera / number-of-threads ，如果以相同的方式执行缩减，您应该得到与 OpenMP 情况相同的结果。

您可以使用 Kahan summation并使用 Kahan 求和实现您自己的归约。然后并行代码应该产生与串行代码相同(更相似)的结果。