mar18-pi3.c

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

#define CIRCLE(step, i) sqrt(4.0 - (double) (i) * (step) * (double) (i) * (step))

int main(int argc, char** argv) {
    int rc;
    rc = MPI_Init(&argc, &argv);
    L17_CHECK_ERROR(MPI_Init, rc);

    int size;
    rc = MPI_Comm_size(MPI_COMM_WORLD, &size);
    L17_CHECK_ERROR(MPI_Comm_size, rc);

    int rank;
    rc = MPI_Comm_rank(MPI_COMM_WORLD, &rank);
    L17_CHECK_ERROR(MPI_Comm_rank, rc);

    int N;

    if (rank == 0) {
        FILE* datp = fopen("mar18-pi3.dat", "r");
        if (datp == NULL) {
            printf("Error. Can't open the data file.\n");
            MPI_Abort(MPI_COMM_WORLD, 1);
        }
        rc = fscanf(datp, "%d", &N);
        if (rc != 1) {
            printf("Error. Can't read the N from the data file.\n");
            MPI_Abort(MPI_COMM_WORLD, 2);
        }
        fclose(datp);
    }

    rc = MPI_Bcast(&N, 1, MPI_INT, 0, MPI_COMM_WORLD);
    L17_CHECK_ERROR(MPI_Bcast, rc);

    int ibeg, iend;
    ibeg = (N / size) * rank;
    iend = (rank == size - 1) ? (N - 1) : ((N / size) * (rank + 1) - 1);
    double step = 2.0 / (double) N;
    double f_old = CIRCLE(step, ibeg);
    double f_new;
    double S = 0;
    int i;
    for (i = ibeg; i <= iend; ++i) {
        f_new = CIRCLE(step, i + 1);
        S += (f_old + f_new) * step / 2.0;
        f_old = f_new;
    }

    rc = MPI_Reduce((rank == 0) ? MPI_IN_PLACE : &S,
                    &S, 1, MPI_DOUBLE, MPI_SUM, 0, MPI_COMM_WORLD);
    L17_CHECK_ERROR(MPI_Reduce, rc);

    if (rank == 0)
        printf("Pi = %f\n", S);

    MPI_Finalize();

    return 0;
}