/**************************pi.c**************************/
#include <stdio.h>
#include <stdlib.h>
#include <stddef.h>
#include <math.h>
#include <time.h>
#include "mpi.h"

#define PI ((double)(4.0*atan(1.0)))
#define TOTALDARTS 240000000UL

/* #define TOTALDARTS 240000000UL */


MPI_Status status;

void srand48();
long lrand48();
double drand48();
double mpdarts(long darts, int seedkey);

/* dartboard subroutine */

double mpdarts(long darts, int seedkey)
{
    double x, y, pi, rand;
    long hits, i;

    /* Seed the random number generator using the procnum to create separate
       series for each independent calculation */

    srand48 ((long) time(0) >> seedkey);
    srand48(lrand48());
    hits = 0;

    /* "throw darts at board" */

    for (i = 1; i <= darts; i++) {

      /* generate random x and y coordinates in the range -1.0 to 1.0 */

        x = (2.0 * drand48()) - 1.0;
        y = (2.0 * drand48()) - 1.0;

	/* if dart lands in the circle, increment hit counter */

        if (x*x + y*y <= 1.0)
            hits++;
    }

    /* calculate pi from the ratio of areas */

  pi = 4.0 * (double)hits/(double)darts;
  return(pi);
}

main(int argc, char **argv)
{
    double localpi, pisum, picomm, gavgpi;
    int i, j, mynum, nprocs, source, check;
    int type=2;
    double sttime, edtime;

    /* All instances call startup  and get their number*/
    MPI_Init(&argc, &argv);
    MPI_Comm_rank(MPI_COMM_WORLD, &mynum);
    MPI_Comm_size(MPI_COMM_WORLD, &nprocs);

    if (mynum==0) { 
      printf ("Starting to calculate pi\n");
      sttime=MPI_Wtime();
    }

    /* Do the computation. Each processor uses the dartboard to estimate pi */
    localpi = mpdarts(TOTALDARTS/nprocs, mynum);

    /* Host processor collects answers from other processors and averages them */

    if (mynum==0) {
      pisum = 0;
      printf ("Node 0 calculated pi = %lf\n", localpi);
      pisum = localpi;
      for (i=1; i<nprocs; i++) {
         source = i;
         check = MPI_Recv(&picomm, 1, MPI_DOUBLE, source, type, MPI_COMM_WORLD, &status);
         printf ("Node 0 got pi = %lf from Node %d\n", picomm,i);
         pisum=pisum+picomm;
      }
      gavgpi = pisum/nprocs;

      printf ("Average value for pi = %lf err %lf\n", gavgpi,gavgpi-PI);

      edtime=MPI_Wtime();
      printf("Total time %.1fs\n", edtime-sttime);
      fflush(stdout);
    }
    /* Other processors just send their value and stop.*/
    else {
      check = MPI_Send(&localpi, 1, MPI_DOUBLE, 0, type, MPI_COMM_WORLD);
      if (check != 0) {
         printf ("Node %d failed to send\n", mynum);
         exit(0);
      }
    }

    MPI_Finalize();
}