/* This is the MPI program computing the value of PI */
#include "mpi.h"
#include <math.h>
#include <stdio.h>
#include <stdlib.h>

int main(int argc, char *argv[])
{
  int done = 0, n, myid, numprocs, i, rc;
  double PI25DT = 3.141592653589793238462643;
  double mypi, pi, h, sum, x, a, tmp;
  MPI_Status status;

  MPI_Init(&argc,&argv);
  MPI_Comm_size(MPI_COMM_WORLD,&numprocs);
  MPI_Comm_rank(MPI_COMM_WORLD,&myid);
  if (argc != 2) {
    if (myid == 0) {
      printf("Usage: a.out N\n");
    } 
    MPI_Finalize();
    exit(0);
  }
  n = atoi(argv[1]);

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

  if (myid != 0) {
	  MPI_Send(&mypi, 1, MPI_DOUBLE, 0, 100, MPI_COMM_WORLD);
  } else {
	  for (i = 1; i<numprocs; i++) {
		  MPI_Recv(&tmp, 1, MPI_DOUBLE, i, 100, MPI_COMM_WORLD, &status);
		  mypi += tmp;
	  }
	  pi = mypi;
  }
  
  if (myid == 0) 
    printf("pi is approximately %.16f, Error is %.16f\n", pi, fabs(pi - PI25DT));
  MPI_Finalize();
  return 0;
}