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

void main() {

  int i, niter;
  double a[1000000], error, sum;

/* OpenMP not needed here since it will be evaluated at compile time */
  for(i = 0; i <  999999; i+=2) {a[i]   = 0.0; a[i+1] = 1.0;}

  niter = 0;

  do { niter++;

#omp parallel for shared(a) private(i)
       for (i = 1; i < 1000000; i+=2) a[i] = (a[i] + a[i-1]) / 2.0;

#omp parallel for shared(a) private(i)
       for (i = 0; i <  999999; i+=2) a[i] = (a[i] + a[i+1]) / 2.0;

       error = 0.0;
#omp parallel for private(i) shared (a) reduction(+:error)
       for (i = 0; i < 1000000; i++) error += fabs(a[i] - a[i+1]);

  } while (error >= 1.0);

  sum = 0.0;
#omp parallel for private(i) shared(a) reduction(+:sum)
  for (i = 0; i < 1000000; i++) sum += a[i];

  printf("Average value = %lf\n", sum / 1000000);
  printf("Number of iterations = %d\n", niter);

}