# include <cstdlib>
# include <iostream>
# include <iomanip>
# include <cmath>
# include <omp.h>


#define N 1000
using namespace std;

int main ( int argc, char *argv[] );
double test01 ( int n, double x[], double y[] );
double test02 ( int n, double x[], double y[] );

int main ( int argc, char *argv[] )
{
  double factor;
  int i;
  int n;
  double wtime;
  double *x;
  double xdoty;
  double *y;

  cout << "\n";
  cout << "DOT_PRODUCT\n";
  cout << "  C++/OpenMP version\n";
  cout << "\n";
  cout << "  A program which computes a vector dot product.\n";

  cout << "\n";
  cout << "  Number of processors available = " << omp_get_num_procs ( ) << "\n";
  cout << "  Number of threads =              " << omp_get_max_threads ( ) << "\n";
  n = N;

  while ( n < 1000000000 )
  {

    x = new double[n];         //starting from 1000 the code is comparing serial vs parallel runs till 1000000000
    y = new double[n];

    for ( i = 0; i < n; i++ )
    {
      x[i] = 1;
    }

    for ( i = 0; i < n; i++ )
    {
      y[i] = 1;
    }

    cout << "\n";
//
//  Test #1
//
    wtime = omp_get_wtime ( );

    xdoty = test01 ( n, x, y ); //sequential testing

    wtime = omp_get_wtime ( ) - wtime;

    cout << "  Sequential"
         << "  " << setw(8) << n
         << "  " << setw(14) << xdoty
         << "  " << setw(15) << wtime << "\n";
//
//  Test #2
//
    wtime = omp_get_wtime ( );

    xdoty = test02 ( n, x, y ); //parallel testing

    wtime = omp_get_wtime ( ) - wtime;

    cout << "  Parallel  "
         << "  " << setw(8) << n
         << "  " << setw(14) << xdoty
         << "  " << setw(15) << wtime << "\n";

    delete [] x;
    delete [] y;
     n = n * 10;
  }
  cout << "\n";
  cout << "DOT_PRODUCT\n";
  cout << "  Normal end of execution.\n";

  return 0;
}

double test01 ( int n, double x[], double y[] )
{
  int i;
  double xdoty;

  xdoty = 0.0;

  for ( i = 0; i < n; i++ )
  {
    xdoty = xdoty + x[i] * y[i];
  }

  return xdoty;
}

double test02 ( int n, double x[], double y[] )
{
  int i;
  double xdoty;

  xdoty = 0.0;

# pragma omp parallel \
  shared ( n, x, y ) \
  private ( i )

# pragma omp for reduction ( + : xdoty )

  for ( i = 0; i < n; i++ )
  {
    xdoty = xdoty + x[i] * y[i];
  }

  return xdoty;
}

