MPI程序实例：数值积分（π值计算）

news2025/1/13 13:51:26

一、算法介绍

二、代码实现

一、算法介绍

本程序实例通过下面的公司计算定积分 $\pi = \int_{0}^{1} 4/(1+x^{2})dx$ 的近似值：

$h\sum_{i=0}^{n-1}f(x_{i})$ (1)

其中 $n>0$ 为积分区间数， $h=1/n$ 为积分步长， $x_{i}=(i+0.5)h\;\;\;(i=0,1,\cdots,n-1)$ 为积分区间的中点，被积函数 $f(x)=4/(1+x^{2})$ 。

假设总进程数为p（程序中的numprocs变量），各进程分别负责计算公式（1）中的一部分计算区间，然后再调用MPI_Reduce将各进程的结果加起来。代码中计算区间采用循环分配的方式，即将计算公式写成：

$\sum_{k=0}^{p-1} h \sum _{0\leqslant j<n,j\;mod \;p=k}f(x_{j})$

每个进程独立计算上式中的一个内层求和，然后再将这项结果加起来。

二、代码实现

#include <stdio.h>
#include "mpi.h"

double f(double a) {return (4.0/(1.0 + a*a));}

int main(int argc, char *argv[])
{
    int n, myid, numprocs, i, namelen;
    double PI25DT = 3.141592653589793238462643;
    double mypi, pi, h, sum, x;
    double startwtime, endwtime;
    char processor_name[MPI_MAX_PROCESSOR_NAME];

    MPI_Init(&argc, &argv);
    MPI_Comm_size(MPI_COMM_WORLD, &numproces);
    MPI_Comm_rank(MPI_COMM_WORLD, &myid);
    MPI_Get_processor_name(processor_name, &namelen)

    printf("Process %d on %s.\n", myid, processor_name);

    if(myid == 0)
    {
        n = 10000;
        startwtime = MPI_Wtime();
    }

    MPI_Bcast(&n, 1, MPI_INT, 0, MPI_COMM_WORLD);

    h = 1.0/(double) n;
    sum = 0.0;
    for(i = myid; i < n; i += numprocs)
    {
        x = h * ((double) i + 0.5);
        sum += f(x);
    }

    mypi = h * sum;

    MPI_Reduce(&mypi, &pi, 1, MPI_DOUBLE, MPI_SUM, 0, MPI_COMM_WORLD);

    if(myid == 0)
    {
        endwtime = MPI_Wtime();
        printf("PI is approximately %.16f, error is %.16f. \n", pi, pi - PI25DT);
        printf("Wall clock time = %f.\n", endwtime - startwtime);
    }
    
    MPI_Finalize();

    return 0;
}

编译并运行上述代码：