title: Sending distributed chunks of a 2D array to the root process in MPI
body: <p>I have a 2D array which is distributed across a MPI process grid (3 x 2 processes in this example). The values of the array are generated within the process which that chunk of the array is distributed to, and I want to gather all of those chunks together at the root process to display them.</p>
  
  <p>So far, I have the code below. This generates a cartesian communicator, finds out the co-ordinates of the MPI process and works out how much of the array it should get based on that (as the array need not be a multiple of the cartesian grid size). I then create a new MPI derived datatype which will send the whole of that processes subarray as one item (that is, the stride, blocklength and count are different for each process, as each process has different sized arrays). However, when I come to gather the data together with MPI_Gather, I get a segmentation fault.</p>
  
  <p>I think this is because I shouldn't be using the same datatype for sending and receiving in the MPI_Gather call. The data type is fine for sending the data, as it has the right count, stride and blocklength, but when it gets to the other end it'll need a very different derived datatype. I'm not sure how to calculate the parameters for this datatype - does anyone have any ideas?</p>
  
  <p>Also, if I'm approaching this from completely the wrong angle then please let me know!</p>
  
  <pre><code>#include&lt;stdio.h&gt;
  #include&lt;array_alloc.h&gt;
  #include&lt;math.h&gt;
  #include&lt;mpi.h&gt;
  
  int main(int argc, char ** argv)
  {
      int size, rank;
      int dim_size[2];
      int periods[2];
      int A = 2;
      int B = 3;
      MPI_Comm cart_comm;
      MPI_Datatype block_type;
      int coords[2];
  
      float **array;
      float **whole_array;
  
      int n = 10;
      int rows_per_core;
      int cols_per_core;
      int i, j;
  
      int x_start, x_finish;
      int y_start, y_finish;
  
      /* Initialise MPI */
      MPI_Init(&amp;argc, &amp;argv);
  
      /* Get the rank for this process, and the number of processes */
      MPI_Comm_size(MPI_COMM_WORLD, &amp;size);
      MPI_Comm_rank(MPI_COMM_WORLD, &amp;rank);
  
      if (rank == 0)
      {
          /* If we're the master process */
          whole_array = alloc_2d_float(n, n);
  
          /* Initialise whole array to silly values */
          for (i = 0; i &lt; n; i++)
          {
              for (j = 0; j &lt; n; j++)
              {
                  whole_array[i][j] = 9999.99;
              }
          }
  
          for (j = 0; j &lt; n; j ++)
          {
              for (i = 0; i &lt; n; i++)
              {
                  printf("%f ", whole_array[j][i]);
              }
              printf("\n");
          }
      }
  
      /* Create the cartesian communicator */
      dim_size[0] = B;
      dim_size[1] = A;
      periods[0] = 1;
      periods[1] = 1;
  
      MPI_Cart_create(MPI_COMM_WORLD, 2, dim_size, periods, 1, &amp;cart_comm);
  
      /* Get our co-ordinates within that communicator */
      MPI_Cart_coords(cart_comm, rank, 2, coords);
  
      rows_per_core = ceil(n / (float) A);
      cols_per_core = ceil(n / (float) B);
  
      if (coords[0] == (B - 1))
      {
          /* We're at the far end of a row */
          cols_per_core = n - (cols_per_core * (B - 1));
      }
      if (coords[1] == (A - 1))
      {
          /* We're at the bottom of a col */
          rows_per_core = n - (rows_per_core * (A - 1));
      }
  
      printf("X: %d, Y: %d, RpC: %d, CpC: %d\n", coords[0], coords[1], rows_per_core, cols_per_core);
  
      MPI_Type_vector(rows_per_core, cols_per_core, cols_per_core + 1, MPI_FLOAT, &amp;block_type);
      MPI_Type_commit(&amp;block_type);
  
      array = alloc_2d_float(rows_per_core, cols_per_core);
  
      if (array == NULL)
      {
          printf("Problem with array allocation.\nExiting\n");
          return 1;
      }
  
      for (j = 0; j &lt; rows_per_core; j++)
      {
          for (i = 0; i &lt; cols_per_core; i++)
          {
              array[j][i] = (float) (i + 1);
          }
      }
  
      MPI_Barrier(MPI_COMM_WORLD);
  
      MPI_Gather(array, 1, block_type, whole_array, 1, block_type, 0, MPI_COMM_WORLD);
  
      /*
      if (rank == 0)
      {
          for (j = 0; j &lt; n; j ++)
          {
              for (i = 0; i &lt; n; i++)
              {
                  printf("%f ", whole_array[j][i]);
              }
              printf("\n");
          }
      }
      */
      /* Close down the MPI environment */
      MPI_Finalize();
  }
  </code></pre>
  
  <p>The 2D array allocation routine I have used above is implemented as:</p>
  
  <pre><code>float **alloc_2d_float( int ndim1, int ndim2 ) {
  
    float **array2 = malloc( ndim1 * sizeof( float * ) );
  
    int i;
  
    if( array2 != NULL ){
  
      array2[0] = malloc( ndim1 * ndim2 * sizeof( float ) );
  
      if( array2[ 0 ] != NULL ) {
  
        for( i = 1; i &lt; ndim1; i++ )
      array2[i] = array2[0] + i * ndim2;
  
      }
  
      else {
        free( array2 );
        array2 = NULL;
      }
  
    }
  
    return array2;
  
  }
  </code></pre>