多线程计算PI碰到的问题

例子如下，用于计算PI的值。gIterations是计算PI的迭代次数，gThreadCount是线程的个数。方法是这样子的，把PI分成gThreadCount个段，分别让一个线程来执行PI的求值操作。求得PI值有两种方法，一种是直接把各个线程每一步所求得的值加到gSum上去，另一种是把各个线程所求得的值加到一个与之对应的全局变量中去。对每个线程i，输出Thread number:I aaaaaaaa，表示线程开始执行，输出Thread number:I bbbbbbb则表示线程执行完毕。有些地方还可以优化的，不过这里只是为了演示多线程的问题，所以就不予关注了。恩。

代码如下。当只有一个thread的时候，结果是OK的（gSum==sum==3.14159*，用等号有点问题，但是结果差异在十万分之一以内）。当有三个threads的时候，问题就开始出现了！gSum计算出来只有2.*！怎么会这样子呢？各位有兴趣的话，可以运行下面的代码试试看。接着看下面的分析。

#include <windows.h>

#include <stdio.h>

#include <time.h>

 

const int gIterations = 100000000;

const int gThreadCount = 3;

double gSum = 0.0;

double gPart[gThreadCount];

 

DWORD WINAPI threadFunction(LPVOID pArg)

{

    int threadNum = (int)pArg;//starts from 0

    printf("Thread number:%d: aaaaaaaaaaaa/n", threadNum);

    for ( int i=threadNum; i<gIterations; i+=gThreadCount )

    {

        double dx = (i + 0.5f) / gIterations;

        gSum += 4.0f / (1.0f + dx*dx);//cause problems here!

        gPart[threadNum] += 4.0f / (1.0f + dx*dx);

    }

 

    printf("part%d value:%.6f/n", threadNum, gPart[threadNum]/gIterations);

    printf("Thread number:%d: bbbbbbbbbbbb/n", threadNum);

    return 0;

}

 

int main()

{

    memset(gPart, 0.0, sizeof(gPart)/sizeof(double));//init to 0

 

    printf("Computing value of Pi: /n");

    clock_t start = clock();

 

    HANDLE threadHandles[gThreadCount];

    for ( int i=0; i<gThreadCount; i++ )

    {

        threadHandles[i] = CreateThread( NULL,           // Security attributes

                                         0,              // Stack size

                                         threadFunction, // Thread function

                                         (LPVOID)i, // Data for thread func()

                                         0,              // Thread start mode

                                         NULL);          // Returned thread ID

    }

 

    WaitForMultipleObjects(gThreadCount, threadHandles, TRUE, INFINITE);

 

    clock_t finish = clock();

    printf("Executing time:%d/n", finish-start);

 

    printf("global: %f/n", gSum / gIterations);

 

    double sum = 0.0;

    for(int i=0; i<gThreadCount; i++)

        sum += gPart[i];

    printf("parts: %f/n", sum / gIterations);

 

    return 0;

}

 

输出信息：

Computing value of Pi:

Thread number:1: aaaaaaaaaaaa

Thread number:0: aaaaaaaaaaaa

Thread number:2: aaaaaaaaaaaa

part1 value:1.047198

Thread number:1: bbbbbbbbbbbb

part0 value:1.047198

Thread number:0: bbbbbbbbbbbb

part2 value:1.047198

Thread number:2: bbbbbbbbbbbb

Executing time:19109

global: 2.711738

parts: 3.141593

Press any key to continue

以上是输出信息通过gSum求出来的值在2.7左右，事实上有的时候还会更低。WHY？问题出现在哪里呢？通过各个线程计算出来的值是对的，说明问题不是出现在这里，也就是说问题是出现在线程切换的时候使得gSum少加了一些值！什么时候切换会导致这个问题呢？问题出现在下面这一句里面：

        gSum += 4.0f / (1.0f + dx*dx);//cause problems here!

这一行等价于：

                   gSum = gSum + value;

这一行代码相当于两行代码：

         temp = gSum + value;

         gSum = temp;

如果有两个线程的话：

线程A:

1、             temp = gSum + value;

2、             gSum = temp;

线程B:

3、             temp = gSum + value;

4、             gSum = temp;

由于线程切换的任意性，这几条指令的执行顺序有以下几种可能：

1 2 3 4，1 3 2 4，1 3 4 2，3 1 2 4，3 1 4 2，3 4 1 2

其中1 3 2 4顺序就是会出错的，很显然按照1 3 2 4顺序的时候1中的value就没有被加进来了。这就是问题所在！同样1 3 4 2，3 1 2 4，3 1 4 2都是有问题。

那如何解决这个问题呢？要把1和2捆绑在一起作为一个单位操作，即所谓原子操作，要么不执行，要么就全都执行了。

正确的代码如下。给gSum+=操作放到一个critical section中，保证此时不会被线程切换干扰。关于critical section的详细信息请参考MSDN。Good luck & have fun.

#include <windows.h>

#include <stdio.h>

 

const int gIterations = 100000;

const int gThreadCount = 4;

double gSum = 0.0;

CRITICAL_SECTION gCS;

 

DWORD WINAPI threadFunction(LPVOID pArg)

{

     double partialSum = 0.0;

 

     for ( int i=(int)pArg+1; i<gIterations; i+=gThreadCount )

     {

         double dx = (i - 0.5f) / gIterations;

         partialSum += 4.0f / (1.0f + dx*dx);

     }

 

     EnterCriticalSection(&gCS);

     gSum += partialSum;

     LeaveCriticalSection(&gCS);

 

     return 0;

}

 

int main()

{

     printf("Computing value of Pi: /n");

 

     InitializeCriticalSection(&gCS);

     HANDLE threadHandles[gThreadCount];

     for ( int i=0; i<gThreadCount; ++i )

     {

         threadHandles[i] = CreateThread( NULL,           // Security attributes

                                          0,              // Stack size

                                          threadFunction, // Thread function

                                          (LPVOID)i,      // Data for thread func()

                                          0,              // Thread start mode

                                          NULL);          // Returned thread ID

     }

     WaitForMultipleObjects(gThreadCount, threadHandles, TRUE, INFINITE);

     DeleteCriticalSection(&gCS);

 

     printf("%f/n", gSum / gIterations);

 

     return 0;

}