/* ==========================================================================
   OpenMP program: compute Pi by integrating the function 1/SQRT(1-x^2)
   ========================================================================== */

#include <iostream>
#include <math.h>
#include <omp.h>
#include <stdlib.h>
#include <sys/time.h>

using namespace std;


double f(double a)
{
   return( 2.0 / sqrt(1 - a*a) );
}

/* =======================
   MAIN
   ======================= */

int main(int argc, char *argv[])
{
   int    N;
   double sum;
   double w, x;

   struct timeval start_time, stop_time;
   int elapsed;

   /* -----
      Check command line
      ----- */
   if ( argc != 2 )
   {
      cout << "Usage: " << argv[0] << " Num_Intervals" << endl;
      exit(1);
   }

   /* -----
      Get number of intervals and number of threads
      ----- */
   N = atoi(argv[1]);
   w = 1.0/(double) N;

   sum = 0.0;

   gettimeofday(&start_time, NULL);

#pragma omp parallel 
{
   int i, num_threads;
   double x;

   num_threads = omp_get_num_threads() ;

   for (i = omp_get_thread_num(); i < N; i = i + num_threads)
   {
      x = w*(i + 0.5);
      #pragma omp critical
      {
         sum = sum +  w*f(x);
      }
   }
}

   cout << "Computed Pi = " << sum << endl << endl;

   gettimeofday(&stop_time, NULL);
   elapsed = (stop_time.tv_sec*1000000 + stop_time.tv_usec) -
                (start_time.tv_sec*1000000 + start_time.tv_usec);
   cout << "Elapsed time = " << elapsed << " microseconds" << endl;

}