title: LU Decomposition MPI
body: <p>This is a MPI code for LU Decomposition. </p>
  
  <p>I have used the following strategy - </p>
  
  <p>There is a master(rank 0) and others are slaves. The master sends rows to each slave.
  Since each slave might receive more than row, I store all the received rows in a 
  buffer and then perform LU Decomposition on it. After doing that I send back the 
  buffer to the master. The master does not do any computation. It just sends and receives.</p>
  
  <pre><code>for(i=0; i&lt;n; i++)
      map[i] = i%(numProcs-1) + 1;
  
  for(i=0; i&lt;n-1; i++)
  {
      if(rank == 0)
      {
          status = pivot(LU,i,n);
  
          for(j=0; j&lt;n; j++)
              row1[j] = LU[n*i+j];
      }
  
      MPI_Bcast(&amp;status, 1, MPI_INT, 0, MPI_COMM_WORLD);
  
      if(status == -1)
          return -1;
  
      MPI_Bcast(row1, n, MPI_DOUBLE, 0, MPI_COMM_WORLD);
  
      int tag1 = 1, tag2 = 2, tag3 = 3, tag4 = 4;
  
      if(rank == 0)
      {
          int pno, start, index, l, rowsReceived = 0;
          MPI_Request req;
          MPI_Status stat;
  
          for(j=i+1; j&lt;n; j++)
              MPI_Isend(&amp;LU[n*j], n, MPI_DOUBLE, map[j], map[j], MPI_COMM_WORLD, &amp;req);
  
          if(i&gt;=n-(numProcs-1))
              cnt++;
  
          for(j=0; j&lt;numProcs-1-cnt; j++)
          {
              MPI_Recv(&amp;pno, 1, MPI_INT, MPI_ANY_SOURCE, tag2, MPI_COMM_WORLD, &amp;stat);
              //printf("1. Recv from %d and j : %d and i : %d\n",pno,j,i);
              MPI_Recv(&amp;rowsReceived, 1, MPI_INT, pno, tag3, MPI_COMM_WORLD, &amp;stat);
              MPI_Recv(rowFinal, n*rowsReceived, MPI_DOUBLE, pno, tag4, MPI_COMM_WORLD, &amp;stat);
  
              /* Will not go more than numProcs anyways */
              for(k=i+1; k&lt;n; k++)
              {
                  if(map[k] == pno)
                  {
                      start = k;
                      break; 
                  }
              }
  
              for(k=0; k&lt;rowsReceived; k++)
              {
                  index = start + k*(numProcs-1);
  
                  for(l=0; l&lt;n; l++)
                      LU[n*index+l] = rowFinal[n*k+l];
              }
          }
      }
  
      else
      {
          int rowsReceived = 0;
          MPI_Status stat, stats[3];
          MPI_Request reqs[3];
  
          for(j=i+1; j&lt;n; j++)
              if(map[j] == rank)
                  rowsReceived += 1;
  
  
          for(j=0; j&lt;rowsReceived; j++)
          {
              MPI_Recv(&amp;rowFinal[n*j], n, MPI_DOUBLE, 0, rank, MPI_COMM_WORLD, &amp;stat);
          }
  
          for(j=0; j&lt;rowsReceived; j++)
          {
              double factor = rowFinal[n*j+i]/row1[i];
  
              for(k=i+1; k&lt;n; k++)
                  rowFinal[n*j+k] -= (row1[k]*factor);
  
              rowFinal[n*j+i] = factor;
          }
  
          if(rowsReceived != 0)
          {
              //printf("Data sent from %d iteration : %d\n",rank,i);
              MPI_Isend(&amp;rank, 1, MPI_INT, 0, tag2, MPI_COMM_WORLD, &amp;reqs[0]);
              MPI_Isend(&amp;rowsReceived, 1, MPI_INT, 0, tag3, MPI_COMM_WORLD, &amp;reqs[1]);
              MPI_Isend(rowFinal, n*rowsReceived, MPI_DOUBLE, 0, tag4, MPI_COMM_WORLD, &amp;reqs[2]);
          }
          //MPI_Waitall(3,reqs,stats);
      }
  }
  </code></pre>
  
  <p>The problem that I am facing is that sometimes the program hangs. My guess is
  that the sends and receives are not being matched but I am not being able to 
  figure out where the problem lies.</p>
  
  <p>I ran test cases on matrices of size 1000x1000, 2000x2000, 3000x3000, 5000x5000 
  and 7000x7000. <strong>Presently the code hangs for 7000x7000</strong>. Could someone please help
  me out?</p>
  
  <p><strong>Things to note :-</strong></p>
  
  <ol>
  <li><p><strong>map</strong> implements the mapping scheme, which row goes to which slave.</p></li>
  <li><p><strong>rowsReceived</strong> tells each slave the no of rows it will receive. I dont need to 
  calculate that each and every time, but I will fix it later.</p></li>
  <li><p><strong>row1</strong> is the buffer in which the active row will be stored.</p></li>
  <li><p><strong>rowFinal</strong> is the buffer of the rows being received and being modified.</p></li>
  <li><p><strong>cnt</strong> is not important and can be ignored. For that the check for
  rowReceived!=0 needs to be removed.</p></li>
  </ol>