/* Example MPI program in C for use with MPE */
/* 2008-01-15*/

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

# define nmax 1000000

void main (int argc, char *argv[]) {

  float x[nmax], y[nmax], sum;
  int npart, nrest, index, iproc;
  int error, rank, size, i;
  MPI_Status status;

  error=MPI_Init(&argc, &argv);                 /* Initialize MPI       */
  error=MPI_Comm_rank(MPI_COMM_WORLD, &rank);   /* Get my rank          */
  error=MPI_Comm_size(MPI_COMM_WORLD, &size);   /* Get the total 
                                                number of processors */
  npart = nmax/size;
  nrest = nmax - npart*size;
  if (nrest != 0){
    printf(" size = %d should be a factor of nmax = %d\n",size,nmax);
    MPI_Finalize();
    exit(0);
  }

  if(rank == 0){
  	for(i=0; i<nmax; i++){
		x[i] = i+1;
	}
	if(size > 1){
		for(iproc=1;iproc<size;iproc++){
			index = iproc*npart;
			error=MPI_Send(&x[index],npart,MPI_FLOAT,iproc,1,MPI_COMM_WORLD);
		}
	}
	for(i=0; i<npart; i++){
		y[i] = powf(x[i],2.3);
	}
	if(size > 1){
		for(iproc=1;iproc<size;iproc++){
			index = iproc*npart;
			error=MPI_Recv(&y[index],npart,MPI_FLOAT,iproc,2,MPI_COMM_WORLD,&status);
	printf(" %d %d %d\n",status.MPI_SOURCE,status.MPI_TAG,status.MPI_ERROR);
		}
	}
	sum = 0;
	for(i=0; i<nmax; i++){
		if(y[i] == 0)
			printf(" y[i] ==0 for i= %d\n",i);
		sum = sum + y[i];
	}
	printf(" sum = %f\n",sum);
  }
  else {
	error=MPI_Recv(x,npart,MPI_FLOAT,0,1,MPI_COMM_WORLD,&status);

	for(i=0; i<npart; i++){
                y[i] = powf(x[i],2.3);
        }
	error=MPI_Send(y,npart,MPI_FLOAT,0,2,MPI_COMM_WORLD);
  }
  error = MPI_Finalize();
}