思路:应用积分的思想求解4/(x2+1)在0-1区间积分的值。将0-1积分区间分成n份,使n尽可能大。这时,被划分的每个小区域就可近似为一个小矩形,其宽为1/n,长可取小区域内的任意高度。求出每个小矩形的面积相加后即得到π值。
在程序中,共n=4个线程,第i个线程负责计算[i/n,(i+1)/n]部分,各线程维护自己的部分和,最后求和得到总的部分和。
OpenMP 是非常方便的
#include <stdio.h>
#include <omp.h>
#define int64_t long long
const int64_t nThreads = 4;
const int64_t num_steps = 1e9;
double step;
double answer_pi, answer = 0.0;
double startTime, runTime;
double sum[nThreads];
int main()
{
step = 1.0 / (double)num_steps;
omp_set_num_threads(nThreads);
answer = 0.0;
startTime = omp_get_wtime();
#pragma omp parallel
{
int64_t id = omp_get_thread_num();
int64_t numthreads = omp_get_num_threads();
double x;
double sum = 0;
for (int64_t i = id * num_steps / nThreads; i < (id + 1) * num_steps / nThreads; i++)
{
x = (i + 0.5) * step;
sum += 4.0 / (1.0 + x * x);
}
#pragma omp critical
answer += sum;
}
answer_pi = step * answer;
runTime = omp_get_wtime() - startTime;
printf("pi=%.14lf
time=%.10lf sec
%lld threads
", answer_pi, runTime, nThreads);
}