多线程并行快速求解Pi的方法

于朋

摘 要 本文首先分别应用不同的求Pi方法分析讲解了WinAPI、OpenMP、MPI三大并行算法中常遇到的一些问题，根据问题提出了相应的解决方案。另外实现了对BBP算法的初步并行化，大大缩减了该算法的运行时间。文章最后利用三种方法求解了蒙特卡洛算法求Pi，并对比分析了三种算法的优缺点。

【关键词】WinAPI OpenMP MPI

1 问题背景与提出

随着技术的发展，单片机越来越难满足人们对于大量数据处理的需求，因此，人们越来越依赖于利用并行计算技术来解决程序规模庞大，运算时间长及数据量大的课题。本文即对当下比较常用的几种并行技术WinAPI、OpenMP、MPI三种并行模式进行讨论课研究，并以计算π值为例，将并行模式与串行模式进行对比，研究并行计算的机理、优缺点及一些常见问题。

2 模型建立

2.1 蒙特卡羅思想，蒲风投针实验

（1）取白纸一张，在上面画许多间距为d的等距平行线；

（2）取一根长为l（l

（3）直线与针相交概率p的近似值可用m/n得到，进而可得到圆周率的近似值为2nl/md。

2.2 级数方法

2.3 并行BBP算法

当今世界在进行计算机性能测试时往往选择在一定时间内计算Pi值得位数，而最常用的算法之一就是BBP算法。该算法不需要多精度浮点算术运算的支持，在支持IEEE浮点运算的通用计算机上即可进行计算而且该算法只需要非常少的内存。BBP算法也是目前算法中非常适用于并行计算的算法。

公式为：

3 算法描述与实现

3.1 API实现

Windows实现多线程并行计算时会产生数据竞争，数据竞争（Data racing）导致计算结果不准确。产生错误的原因在于两个线程在同时访问同一内存区域时，且一个线程在进行写操作。为此采用同步技术，利用临界区解决此问题。由于本次试验中，数据计算量较小，所以基本不会产生数据竞争的错误，但是一旦运行计算量较大或者在共享量上运算时间较长的程序时，就会产生数据竞争，在本例中，我们以cout输出程序为引子，观察使用临界区和不使用临界区的区别：

代码如下：

}

pi= count*4.0 / numSteps；

//EnterCriticalSection（&cs）；

cout << piSum<< " " << threadID << endl；

piSum+= pi；

// LeaveCriticalSection（&cs）；

我们将输出程序置于piSum值之前，由于cout运行时，在屏幕上显示需要消耗较大时间，所以就会有机会产生数据重叠或者两者读到了一样的piSum值。（如图1所示）。

所以，当我们使用临界区后结果为：（如图2所示）。

另外，临界区的设定也直接影响计算速度，所以尽量缩小临界区有利于提高运行速度。

线程取值过大的影响：

3.2 OpenMP实现级数方法

OpenMP是一种面向共享内存以及分布共享内存的多处理器多线程并行编程语言，是一种能够被应用于显式指导多线程、共享内存并行的应用程序编程接口（如图3所示）。OpenMP具有良好可移植性，支持多种编程语言。parallel for 循环并行化利用编译指导语句parallel for 并行原理是采用工作分配的执行方式，将循环所需要的工作量按一定方式分配到各个线程，所有线程执行工作总和是原串行完成的工作量。parallel for 根据CUP的大小（线程数多少）按工作量平均分配，对于线程数为8个，这里蒙特卡洛方法工作量取numSteps=10，000，000步，则每个线程要计算（1/8）*numSteps步的工作量。要特别注意的是在这里必须避免数据竞争（Data racing），采用的方法是对竞争变量私有化，使用private（t，fun，...）。

核心代码：

double Pi（double x）

{

omp_set_num_threads（numThreads）；

int sign = 1；

#pragma omp parallel for reduction（+：c[omp_get_thread_num（）]）

for （int i = 2； i < N； i++）

{

c[omp_get_thread_num（）] += 4 * x；

sign = sign * （-1）；

x = sign /double （2 * i - 1）；

}

这个算法中，我们发现运行后，串行时间和并行时间是基本一致的：

这个问题在于x，sign的值为共享量，当不同的线程在调用时存在等待的时间，所以如果将这些值私有化，则可以解决这些问题。由于使用private时，各个私有变量必须在for中有定义，即在进入并行区时，for循环中的所有私有变量是需要重新定义的，因此我们必须使得x，sign等变量在定义时全部使用已知量来定义

#pragma omp parallel for private（x，y）reduction（+：pi）

for （int i = 1； i < N； i += 2）

{

y = 1 / double（2 * i - 1）；

pi = pi + 4 * y；

x = -1 / double（2 * （1 + i） - 1）；

pi = pi + 4 * x；

}

因此程序完美的解决了并行问题。（如表1所示）。

3.3 BBP算法的实现

MPI是专门为集群和多节点并行计算语言，运行效率高，实现方式多样，可以进行主从式、并列式以及流水线式等方式的实现。在利用计算机多核CPU，模拟多个节点的实现方式。改造成主从式程序，利用0号节点作为主节点收发数据，也参与计算，而其他节点只进行计算，为负载均衡选择详尽的计算量运行到不同的节点上。

根据BBP算法（如图4所示）我们将整个运算分成4个部分，分别用一个线程进行运算，其运算总时间约等于线程中运行时间最长的那个（计算结果如表2所示）。

3.4 API、OpenMP、MPI对比实验

我们以蒙塔卡罗算法为例，分别采用API、OpenMP、MPI三种并行方式对值进行计算，來对比其运算效率和运算精度（计算结果如表3所示）。

4 结论

（1）MPI模拟多节点计算的速度是最快的，而且计算结果也是最接近穿行结果的，所以MPI适合计算大规模的较为精确的求解问题。

（2）WinAPI实现用临界区的效果最差，而且WinAPI的计算速度很大程度也和临界区的设定有关，尽量缩小临界区有利于提高并行速度。

（3）WinAPI实现用线程号为每个线程分配不同的任务，分配过程需要人为干预，对于函数复杂的程序来说，实现繁琐，不利于使用。

（4）数值计算时要根据实际情况选择和改造算法。比如在计算值就比e更适合并行。而且每种并行方法都有它的特使要求，比如在使用parallel for private时，由于变量进入for循环后属于重新定义，所以不能出现自己为自己赋值的情况，需要按照程序一步步展开来写，相对繁琐。

参考文献

[1]多核系列教材编写组.多核程序设计[M].北京：清华大学出版社，2007.

[2]朱建伟，刘荣.多线程并行快速求解e值的六种方法[J].现代计算机，2013.

[3]张翔圆.周率Pi的BBP多核并行算法实现[J].普洱学院学报，2013.

作者单位

中国石油大学（华东）理学院 山东省青岛市 266500