#include <math.h>
#include <stdio.h>
#include <stdlib.h>
#include <omp.h>

int main(int argc, char *argv[])
{
  int i, n, id, num_threads;
  double PI25DT = 3.141592653589793238462643;
  double mypi, pi, h, sum, x, a;

  if (argc != 2) {
    printf("Usage: a.out N\n");
    exit(0);
  }
  n = atoi(argv[1]);

  h   = 1.0 / (double) n;
  sum = 0.0;

#pragma omp parallel private (i, id, x, num_threads) shared (h,n) reduction(+:sum)
    {
    id = omp_get_thread_num();
    num_threads= omp_get_num_threads();
        
  for (i = id*(n/num_threads)+1; i <= (id+1)*(n/num_threads); i ++) {
    x = h * ((double)i - 0.5);
    sum += 4.0 / (1.0 + x*x);
  }
}        

  mypi = h * sum;
  printf("pi is approximately %.16f, Error is %.16f\n", mypi, fabs(mypi - PI25DT));
  return 0;
}