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

#define index3D(a,b,c,LD1,LD2)  a + (b*LD1) + (c*LD1*LD2)

int MINVAL(int s, int*a)
{
   int v, e;

   v = a[0];

   for ( e = 0; e < s; e++ )
   {
       if ( v > a[e] ) v = a[e];
   
   }
   
   return(v);
}

int MAXVAL(int s, int*a)
{
   int v, e;

   v = a[0];

   for ( e = 0; e < s; e++ )
   {
       if ( v < a[e] ) v = a[e];
   
   }
   
   return(v);
}
   
   
void Terminate(int p, int err)
{
#include "mpi.h"
   int ierr, rc;
   
   fprintf(stdout,"Process %d  terminates\n",p);
   if ( err != 0 ) {
   	   ierr = MPI_Abort(MPI_COMM_WORLD,err);
   } else {
   	   ierr = MPI_Finalize();
   }
   
   return;
}
   
void IntData2Mpgm(int s1, int s2, int* idata, char* name, int Mag) 
{
   /*
      Simple subroutine to dump integer data in a PGM format
   */
   int i, j, r, m;
   FILE* ouni;
   int  vp, vs;
   float  rmin, rmax;
   char  fname[80];
   /*
      Write on unit 700 with  PGM format
   */
   strcpy(fname,name);
   strcat(fname,".pgm\0");
   ouni = fopen(fname,"w");
   if ( ! ouni ) {
      fprintf(stderr,"!!!! Error write access to file %s\n",fname);
   }
   /*  Magic code */
   fprintf(ouni,"P2\n");
   /*  Dimensions */
   fprintf(ouni,"%d %d\n",s1*Mag, s2*Mag);
   /*  Maximum value */
   fprintf(ouni,"255\n");
   /*  Values from 0 to 255 */
   rmin = MINVAL(s1*s2,idata); rmax = MAXVAL(s1*s2,idata);
   vs = 0;
   for ( i = 0; i < s1; i++ ) {
   	   for ( r = 0; r < Mag; r++ ) {
   	   	   for ( j = 0; j < s2; j++ ) {
   	   	   	   vp = (int) ( (idata[i+(j*s1)] - rmin) * 255.0 / (rmax - rmin) );
   	   	   	   for ( m = 0; m < Mag; m++) {
   	   	   	   	   vs++;
   	   	   	   	   fprintf(ouni,"%4d ", vp);
   	   	   	   }
   	   	   	   if ( vs >= 10 ) {
   	   	   	   	   fprintf(ouni," \n");
   	   	   	   	   vs = 0;
   	   	   	   }
   	   	   }
   	   }
      fprintf(ouni," ");
      vs = 0;
   }
   fclose(ouni);
   return;
}


int ProcessBand(int length, int* start, int* end, int intervals, int n)
{
    int delta;

    delta = ( length + (intervals - 1) ) / intervals;
    *start = delta * n; *end = delta * ( n + 1) - 1;
    if ( *end > length-1 ) *end = length-1;
   
   return(0);
}

main( int argc, char *argv[] )   /* Program LIFE */
{
//     MaxLoop iterations of an NxN life board

      int N, MaxLoop, L1, L2;
      int *step;
      int loop, i, j, cr, nx, r, c, nb, ls;
      char name[16];
      double  wall1, wall2, wall;
//    
//     MPI variables
	  int my_rank, numprocs, ierr, tag;
	  MPI_Status status;
	  char server[MPI_MAX_PROCESSOR_NAME];
//  
	  int p, my_start, my_end, my_out_col0, my_out_col1, my_cols, loc_out_col0, loc_out_col1;
	  int my_inn_col0, my_inn_col1, loc_inn_col0, loc_inn_col1, istart, iend;
//
//     MPI initialization
//
	  ierr = MPI_Init(&argc,&argv);
	  ierr = MPI_Comm_rank(MPI_COMM_WORLD,&my_rank);
      ierr = MPI_Comm_size(MPI_COMM_WORLD,&numprocs);
	  ierr = MPI_Get_processor_name(server,&ls);
	  printf(" Activated process %d of %d on server %s\n", my_rank,numprocs,server);
//
//     Initialization
//
	  if ( my_rank == 0 ) {
	  	  printf(" N, MaxLoop ?\n");
                  N = 200; MaxLoop = 32;
	//  	  fscanf(stdin,"%d %d",&N,&MaxLoop);
	  	  printf(" N, MaxLoop = %d %d",N, MaxLoop);
	  	  // Allocate full step board
	  	  step = (int*)malloc(sizeof((int)1)*N*N*2);
                  L1 = N; L2 = N;
	  	  if ( step == NULL ) {
	  	  	  fprintf(stderr," Error allocating step(N,N,2)\n");
	  	  	  ierr = -1;
	  	  	  Terminate(my_rank,ierr);
	  	  	  return(-1);
	  	  }
//         Initialize step
          for ( i = 0; i < N*N*2; i++ ) step[i] = 0;
//         Ten active cells in a line

          for ( i = (N/2)-5; i < (N/2)+5; i++ ) {
          	  step[index3D(N/2,i,0,L1,L2)] = 1;
          }
//
//         Print starting config to file Life0000.pgm
//
	      loop = 0;
	      sprintf(name,"Life%4.4d",loop);
          IntData2Mpgm(N,N,&step[index3D(0,0,0,L1,L2)],name,4);
      }
      //  Send initial data to processes
      ierr = MPI_Bcast(&N,1,MPI_INT,0,MPI_COMM_WORLD);
      ierr = MPI_Bcast(&MaxLoop,1,MPI_INT,0,MPI_COMM_WORLD);
      // print*," Process ",my_rank,": N, MaxLoop = ",N, MaxLoop
      if ( my_rank == 0 ) {
         ProcessBand(N,&my_start,&my_end,numprocs,my_rank);    //  Global columns
         my_out_col0 = (my_start-1 >= 0) ? (my_start-1) : 0;   // First allocated column in global numeration
         my_out_col1 = (my_end+1 < N) ? (my_end+1) : N-1;   // Last allocated column in global numeration
         my_inn_col0 = my_out_col0 + 1;
         my_inn_col1 = my_out_col1 - 1;
         my_cols = (my_out_col1-my_out_col0+1); // Total allocated columns
         for ( p = 1; p < numprocs; p++ ) {
         //  Send process data
            ProcessBand(N,&istart,&iend,numprocs,p);    //  Global columns
            tag = 100;   // istart
            ierr = MPI_Send(&istart,1,MPI_INT,p,tag,MPI_COMM_WORLD);    
            if ( ierr != 0 ) {
	  	  	   fprintf(stderr," Error MPI_Send(istart,  from 0 to %d\n",p);
	  	  	   ierr = -1;
	  	  	   Terminate(my_rank,ierr);
	  	  	   return(-1);
	  	  	}
            // print*," Process ",my_rank,": sent istart"
//
            tag = 110;   // iend
            ierr = MPI_Send(&iend,1,MPI_INT, p,tag,MPI_COMM_WORLD);  
            if ( ierr != 0 ) {
	  	  	   fprintf(stderr," Error MPI_Send(iend,  from 0 to %d\n",p);
	  	  	   ierr = -1;
	  	  	   Terminate(my_rank,ierr);
	  	  	   return(-1);
	    }
            // print*," Process ",my_rank,": sent iend"

	  	    tag = 120;   // step board
            ierr = MPI_Send(&step[index3D(0,istart,0,L1,L2)],(iend-istart+1)*N,MPI_INT, 
            	p,tag,MPI_COMM_WORLD);      
            if ( ierr != 0 ) {
	  	  	   fprintf(stderr," Error MPI_Send(step(1,istart,1) from 0 to %d\n",p);
	  	  	   ierr = -1;
	  	  	   Terminate(my_rank,ierr);
	  	  	   return(-1);
	  	  	}
            // print*," Process ",my_rank,": sent step(1,",istart,",1),",(iend-istart+1),","
         }
      } else {   // my_rank > 0   
      	  tag = 100;   // istart
          ierr = MPI_Recv(&my_start,1, MPI_INT, 0, tag, MPI_COMM_WORLD, &status);   
          if ( ierr != 0 ) {
	  	  	   fprintf(stderr," Error MPI_Recv(istart,1, from 0 in process %d\n",
	  	  	   	   my_rank);
	  	  	   ierr = -1;
	  	  	   Terminate(my_rank,ierr);
	  	  	   return(-1);
	  	  }
          // print*," Process ",my_rank,": recv my_start"

	  	  tag = 110;   // iend
          ierr = MPI_Recv(&my_end,1, MPI_INT, 0, tag, MPI_COMM_WORLD, &status);     	  
          if ( ierr != 0 ) {
	  	  	   fprintf(stderr," Error MPI_Recv(iend,1, from 0 in process %d",my_rank);
	  	  	   ierr = -1;
	  	  	   Terminate(my_rank,ierr);
	  	  	   return(-1);
	  }
          // print*," Process ",my_rank,": recv my_end"

	  my_out_col0 = (my_start-1 >= 0) ? (my_start-1) : 0;   // First allocated column in global numeration
          my_out_col1 = (my_end+1 < N ) ? (my_end+1) : N-1;     // Last allocated column in global numeration
          my_inn_col0 = my_out_col0 + 1;
          my_inn_col1 = my_out_col1 - 1;
          my_cols = (my_out_col1-my_out_col0+1); // Total allocated columns
      	  step = (int*)malloc(sizeof((int)1)*N*my_cols*2);
          L1 = N; L2 = my_cols;
	  	  if ( step == NULL ) {
	  	  	  fprintf(stderr," Error allocating step(N,my_cols,2) in process %d",
	  	  	  	  my_rank);
	  	  	   ierr = -1;
	  	  	   Terminate(my_rank,ierr);
	  	  	   return(-1);
	  	  }
	  	  for ( i = 0; i < N*my_cols*2; i++ ) step[i] = 0;   // board initialization
	  	  tag = 120;   // step board
          ierr = MPI_Recv(&step[index3D(0,(my_start-my_out_col0),0,L1,L2)],(my_end-my_start+1)*N, 
          	  MPI_INT, 0, tag, MPI_COMM_WORLD, &status);
          if ( ierr != 0 ) {
	  	  	   fprintf(stderr," Error MPI_Recv( step(1,(my_start-my_out_col0+1),1) from 0 in process %d",my_rank);
	  	  	   ierr = -1;
	  	  	   Terminate(my_rank,ierr);
	  	  	   return(-1);
	  	    }
          // print*," Process ",my_rank,": recv step(1,",(my_start-my_out_col0+1),",1),",(my_end-my_start+1),","
      }
      loc_out_col0 = 0;   // First allocated column in local numeration
      loc_inn_col0 = 1;   // First column to be computed in local numeration
      loc_out_col1 = loc_out_col0 + my_out_col1 - my_out_col0;
      loc_inn_col1 = loc_inn_col0 + my_inn_col1 - my_inn_col0;   // Last column to be computed in local numeration
//
//     Perform MaxLoop updates
//  
      for (loop = 1; loop <= MaxLoop; loop++) {
         cr = (loop-1)%2;
         nx = loop%2;
         // Send border values
         if ( my_end < N-1 ) {
 	  	    tag = 310;   // Send my rightmost column to left border of receiver
            ierr = MPI_Send(&step[index3D(0,loc_inn_col1,cr,L1,L2)],N, MPI_INT, my_rank+1, tag, 
            	MPI_COMM_WORLD);
            if ( ierr != 0 ) {
	  	  	   fprintf(stderr," Error MPI_Send( step(1,loc_inn_col1+1,1),N, in process %d\n",my_rank);
	  	  	   ierr = -1;
	  	  	   Terminate(my_rank,ierr);
	  	  	   return(-1);
	  	  	}
            // print*," Process ",my_rank,": sent step(1,",loc_inn_col1,",cr),N, "
         }  
         if ( my_start > 0 ) {   
 	  	    tag = 320;   // Send my leftmost column to right border of receiver
            ierr = MPI_Send(&step[index3D(0,loc_inn_col0,cr,L1,L2)],N, MPI_INT, my_rank-1, tag, 
            	MPI_COMM_WORLD);
            if ( ierr != 0 ) {
	  	  	   fprintf(stderr," Error MPI_Send( step(1,loc_inn_col0-1,1),N, in process %d\n",my_rank);
	  	  	   ierr = -1;
	  	  	   Terminate(my_rank,ierr);
	  	  	   return(-1);
	  	  	}
            // print*," Process ",my_rank,": sent step(1,",loc_inn_col0,",cr),N, "
         }
         // Acquire border values
         // print*," Process ",my_rank,": my_start, my_end = ",my_start, my_end
         if ( my_start > 0 ) {   
 	  	    tag = 310;   // Acquire left border
            ierr = MPI_Recv(&step[index3D(0,loc_out_col0,cr,L1,L2)],N, MPI_INT, my_rank-1, tag, 
            	MPI_COMM_WORLD, &status);
            if ( ierr != 0 ) {
	  	  	   fprintf(stderr," Error MPI_Recv( step(1,loc_out_col0,cr),N, in process %d\n",my_rank);
	  	  	   ierr = -1;
	  	  	   Terminate(my_rank,ierr);
	  	  	   return(-1);
	  	  	}
            // print*," Process ",my_rank,": received step(1,",loc_out_col0,",cr),N, "
         }  
         if ( my_end < N-1 ) {   
 	  	    tag = 320;   // Acquire right border
            ierr = MPI_Recv(&step[index3D(0,loc_out_col1,cr,L1,L2)],N, MPI_INT, my_rank+1, tag, 
            	MPI_COMM_WORLD, &status);
            if ( ierr != 0 ) {
	  	  	   fprintf(stderr," Error MPI_Recv( step(1,loc_out_col1,1),N, in process %d\n",my_rank);
	  	  	   ierr = -1;
	  	  	   Terminate(my_rank,ierr);
	  	  	   return(-1);
	  	  	}
            // print*," Process ",my_rank,": received step(1,",loc_out_col1,",cr),N, "
         }
         for (r = 1; r < N-1; r++) {
            for ( c = loc_inn_col0; c <= loc_inn_col1; c++ ) {   // local columns
               step[index3D(r,c,nx,L1,L2)] = step[index3D(r,c,cr,L1,L2)];   // Inizialize cell
               //  Per each cell
               nb = 0;
               for ( i = r-1; i <= r+1; i++ ) {
               	   for ( j = c-1; j <= c+1; j++ ) {
                      // Count active neighbours
                  	  if ( ! (i == r && j == c) ) {
                  	  	  nb = nb + step[index3D(i,j,cr,L1,L2)];
                  	  }
                   }
               }
               if  (nb < 2 || nb > 3)
                  step[index3D(r,c,nx,L1,L2)] = 0;
               else if ( nb == 3 )
                  step[index3D(r,c,nx,L1,L2)] = 1;
            }
         }
//
//        Send data to root process
//        
	 tag = 200;
	 if ( my_rank == 0 ) {
	 	 for ( p = 1; p <= numprocs-1; p++ ) {
               ProcessBand(N,&istart,&iend,numprocs,p);
               ierr = MPI_Recv(&step[index3D(0,istart,nx,L1,L2)],(iend-istart+1)*N,MPI_INT, p,tag,
               	   MPI_COMM_WORLD,&status);      
               if ( ierr != 0 ) {
	  	           fprintf(stderr," Error MPI_Recv(step(1,istart,nx) from %d\n",p);
	  	  	       ierr = -1;
	  	  	       Terminate(my_rank,ierr);
	  	  	       return(-1);
	           }
            }
	 } else {   // Send data to root 
            ierr = MPI_Send(&step[index3D(0,(my_start-my_out_col0),nx,L1,L2)],(my_end-my_start+1)*N, MPI_INT, 0,tag,MPI_COMM_WORLD);   
            if ( ierr != 0 ) {
     	        fprintf(stderr," Error MPI_Send(step(1,(my_start-my_out_col0+1),1),(my_end-my_start+1) from process %d\n",my_rank);
	  	  	   ierr = -1;
	  	  	   Terminate(my_rank,ierr);
	  	  	   return(-1);
	  	  	}
	 }	    
     if ( my_rank == 0 ) {
//
//           Write out computed state of board
//
       	    sprintf(name,"Life%4.4d\0",loop);
            IntData2Mpgm(N,N,&step[index3D(0,0,nx,L1,L2)],name,4);
         }
      }   // End loop over steps

          ierr = MPI_Barrier(MPI_COMM_WORLD);
	  Terminate(my_rank,ierr);
	  return(0);

}  /* end program LIFE */