/*
 * Reduces function T from a 2 D grid by summing over y direction into a 1 D grid
 * Following the standards of PAMR
 */
void Reduce(double *T, double *T_r) {
   int rc;
   int source, dest, tag;
   MPI_Status stat;
   int loc_shape[2];
   double loc_x;
   int j, shift;
   double loc_T[2000];
   int T_r_filled[2000];
   int loc_phys_bdy[4];
   int loc_ghost_width[4];
   int left_s, right_s;
   int ltrace = 0;
   int i;
   source = my_rank;
   dest = 0;
   tag = my_rank + 100;
   /******************** MPI SEND ************************************/
   rc = MPI_Send(&shape[0],2,MPI_INT,dest,tag,MPI_COMM_WORLD);
   rc = MPI_Send(&x[0],1,MPI_DOUBLE,dest,tag,MPI_COMM_WORLD);
   rc = MPI_Send(&T[0],shape[0],MPI_DOUBLE,dest,tag,MPI_COMM_WORLD);
   rc = MPI_Send(&phys_bdy[0],4,MPI_INT,dest,tag,MPI_COMM_WORLD);
   rc = MPI_Send(&ghost_width[0],4,MPI_INT,dest,tag,MPI_COMM_WORLD);
   /*******************END OF MPI SEND *********************************/
   if (my_rank == 0) {
            my_zero(T_r,Nx_r);
       if (ltrace) {
            my_zero_i(T_r_filled,2000);
       }
       for (i=0; i<numprocs; i++) {
          source = i;
          tag = i + 100;
             /************* MPI RECEIVE **********************************/
                  rc = MPI_Recv(&loc_shape[0], 2, MPI_INT,source, tag,
                  MPI_COMM_WORLD,&stat);
               rc = MPI_Recv(&loc_x, 1, MPI_DOUBLE,source, tag,
                    MPI_COMM_WORLD,&stat);
                    rc = MPI_Recv(&loc_T[0],loc_shape[0],
                  MPI_DOUBLE,source,tag, MPI_COMM_WORLD,&stat);
             rc = MPI_Recv(&loc_phys_bdy[0],4,MPI_INT,source,tag,
                MPI_COMM_WORLD,&stat);
             rc = MPI_Recv(&loc_ghost_width[0],4,MPI_INT,source,tag,
                MPI_COMM_WORLD,&stat);
           /***************** END OF MPI RECEIVE *************************/
            shift = round( (loc_x - base_bbox[0] ) / dx );
          
          if (loc_phys_bdy[0] == 0) {
         left_s = loc_ghost_width[0];
          } else {
             left_s = 0;
          }
          if (loc_phys_bdy[1] == 0) {
             //note that CPU's have one grid point in common
             //even without any ghost cells, that's why
             //I have +1 here:
             right_s = loc_ghost_width[1] + 1;
               } else {
            right_s = 0;
          }
          if (ltrace ) {
            printf("\n====================\n");
            printf("from CPU: %d\n",source);
            printf("shift = %d \n",shift);
            printf("loc_shape = %d x %d \n", loc_shape[0], loc_shape[1]);
            printf("left_s = %d, right_s = %d\n",left_s,right_s);
            printf("====================\n");
          }
          for(j=left_s;j<loc_shape[0]-right_s;j++){
               T_r[j+shift] += loc_T[j];
             if (ltrace) {
                 T_r_filled[j+shift]++;
             }
          }
       } //End of binning from all procs
       if (ltrace) {
      for ( j=0; j<Nx_r;j++ ){ printf("%d   ", T_r_filled[j]); }
      printf("\n");
       }
      } //End of rank=0 job
}

/*
 * Distributes a 1 D data to 2 D grids according to the position of each chunck of
 * 2D grid on each CPU. Following the standards of PAMR
 */

void Distribute(double *f_r, double *f) {
   int rc;
   int source, dest, tag;
   MPI_Status stat;
   int loc_shape[2];
   double loc_x;
   int shift;
   int ltrace = 0;
   int i;
   int k;
   source = my_rank;
   dest = 0;
   tag = my_rank + 100;
   /******************** MPI SEND ************************************/
   rc = MPI_Send(&shape[0],2,MPI_INT,dest,tag,MPI_COMM_WORLD);
   rc = MPI_Send(&x[0],1,MPI_DOUBLE,dest,tag,MPI_COMM_WORLD);
        /*******************END OF MPI SEND *********************************/
        if (my_rank == 0) {
          for (i=0; i<numprocs; i++) {
          source = i;
          tag = i + 100;
               /************* MPI RECEIVE **********************************/
                 rc = MPI_Recv(&loc_shape[0], 2, MPI_INT,source, tag,
                               MPI_COMM_WORLD,&stat);
               rc = MPI_Recv(&loc_x, 1, MPI_DOUBLE,source, tag,
                               MPI_COMM_WORLD,&stat);
            /***************** END OF MPI RECEIVE *************************/
            shift = round( (loc_x - base_bbox[0] ) / dx );
          
               if (ltrace ) {
            printf("\n====================\n");
            printf("from CPU: %d\n",source);
            printf("shift = %d \n",shift);
            printf("loc_shape = %d x %d \n", loc_shape[0], loc_shape[1]);
            printf("====================\n");
          }
          dest = i;
          tag = i + 200;
          rc = MPI_Send(&f_r[shift],loc_shape[0],MPI_DOUBLE,dest,tag,MPI_COMM_WORLD);
       } //End of distributiong from rank=0 to all Procs
      } //End of rank=0 job
      source = 0;
      tag = my_rank + 200;
      rc = MPI_Recv(&f[0],shape[0],MPI_DOUBLE,source,tag,MPI_COMM_WORLD,&stat);
      //Making the dummy copy:
      if (make_dummy_copy == 1) {
         for (k=shape[0];k<shape[0]*shape[1];k++) {
             f[k] = f[k % shape[0]];
         }
      }
      return;
}