/* mpi.collective.example.c SJ */

#include <math.h>
#include <stdio.h>
#include <stdlib.h>
#include <unistd.h>
#include <mpi.h>

/* these variables are used by prefix sum also, hence here */
static int P, PID;

int prefix_sum_28(int value);

int main(int argc, char **argv)
{
  int i, num, psum, ksum, *arr, *input;

  MPI_Init(&argc, &argv);
  MPI_Comm_rank(MPI_COMM_WORLD, &PID);
  MPI_Comm_size(MPI_COMM_WORLD, &P);

  arr = (int*)malloc(P*sizeof(int));

  /* 0 creates input, sends */
  if (PID == 0) {
	  srand(getpid());
	  input = (int*)malloc(P*sizeof(int));
	  printf("Input is:");
	  for (i = 0; i < P; i++) {
		  input[i] = rand()%10;
		  printf(" %d", input[i]);
	  }
	  printf("\n");
  }

  /* 0 sends pieces of data to everyone */
  MPI_Scatter(input, 1, MPI_INT, &num, 1, MPI_INT, 0, MPI_COMM_WORLD);
	  
  printf("%d: num = %d\n", PID, num);

  /* count something */
  psum = prefix_sum_28(num);

  /* collect the prefix sum to the first process */
  MPI_Gather(&psum, 1, MPI_INT, arr, 1, MPI_INT, 0, MPI_COMM_WORLD);
  if (PID == 0) {
	  printf("Gather:");
	  for (i = 0; i < P; i++)
		  printf(" %d", arr[i]);
	  printf("\n");
  }

  /* make sum of all values to process 0 */
  MPI_Reduce(&num, &ksum, 1, MPI_INT, MPI_SUM, 0, MPI_COMM_WORLD);
  if (PID == 0) {
	  printf("Reduce: %d\n", ksum);
  }

  /* same again, but for all processes */
  MPI_Allreduce(&num, &ksum, 1, MPI_INT, MPI_SUM, MPI_COMM_WORLD);
  printf("%d: Allreduce: %d\n", PID, ksum);


  /* prefix sum over all processes */
  MPI_Scan(&num, &psum, 1, MPI_INT, MPI_SUM, MPI_COMM_WORLD);
  MPI_Gather(&psum, 1, MPI_INT, arr, 1, MPI_INT, 0, MPI_COMM_WORLD);
  if (PID == 0) {
	  printf("Scan&Gather:");
	  for (i = 0; i < P; i++)
		  printf(" %d", arr[i]);
	  printf("\n");
  }


  free(arr);

  MPI_Finalize();

  exit(0);

}
















/* prefix sum */
int prefix_sum_28(int value)
{
  int sum = value, i, mess;
  MPI_Status stat;

  /* repeat logP times */
  /* i is like 2^(1-i) at original algorithm */
  i = 1;
  while (i < P) {

	  /* for j = i to N-1 pardo, PID = j */

	  /* if someone needs this value, send */
	  if (PID + i < P)
		  MPI_Send(&sum, 1, MPI_INT, PID + i, 0, MPI_COMM_WORLD);
			   
	  /* if this processor has still work to do, receive */
	  if (PID >= i) {
		  MPI_Recv(&mess, 1, MPI_INT, PID - i, 0, MPI_COMM_WORLD, &stat);
		  sum += mess;
	  }

	  i *= 2;
  }

  return sum;
}