LAM/MPI logo

LAM/MPI General User's Mailing List Archives

  |   Home   |   Download   |   Documentation   |   FAQ   |   all just in this list

From: Stephan Mertens (Stephan.Mertens_at_[hidden])
Date: 2005-05-12 05:37:06


On May 11, 2005, at 16:00, Jeff Squyres wrote:

> I assume you're using some flavor of BSEND, right?

Yes.

> Second: even if you do use BSEND, it should trigger at least some
> progress every time you invoke BSEND (or most any other communication
> MPI function).

No, MPI_Bsend does never trigger any progress on pending messages if each
individual message is larger than the short protocol! This leads to the
bizarre situation that a sender that is *exclusively* using MPI_Bsend will
never get a message through before he calls MPI_Finalize. Before that
he will die from buffer exhaustion, of course. Run the program
below and you will see...

/* Bsend.c - demonstrate the futility of buffered sends */

#include <stdio.h>
#include <getopt.h>
#include <stdlib.h>
#include <time.h>
#include <unistd.h>
#include <mpi.h>

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

  int myrank, nprocs, i, k;
  char *buf, *mpibuf, c;
  MPI_Request req;
  MPI_Status status;
  time_t T;
  long bufsize;

  MPI_Init(&argc, &argv);

  MPI_Comm_size(MPI_COMM_WORLD, &nprocs);
  MPI_Comm_rank(MPI_COMM_WORLD, &myrank);

  bufsize = 100000; /* override with Bsend -b <bufsize> */
  while ((c = getopt (argc, argv, "b:")) != -1)
  {
    switch (c)
    {
    case 'b':
      bufsize = atol (optarg);
      break;
    }
  }

  buf = (char*)malloc (bufsize*sizeof(char));

  MPI_Barrier (MPI_COMM_WORLD);

  if (myrank == 0) {
    printf ("Buffer size is %ld\n", bufsize*sizeof(char));
    mpibuf = (char*)malloc (3*bufsize*sizeof(char));
    /* let MPI buffer at least 3 messages */
    MPI_Buffer_attach (mpibuf, 4*bufsize*sizeof(char));
    for (k = 0; k < 10; k++) {
      time(&T);
      printf("sending message %d at %s", k+1, ctime(&T));
      MPI_Bsend (buf, bufsize, MPI_CHAR, nprocs-1, 99, MPI_COMM_WORLD);
      sleep (5);
    }
    MPI_Buffer_detach (mpibuf, &i);
    free (mpibuf);
  }
  else if (myrank == nprocs-1) {
    for (k = 0; k < 10; k++) {
      MPI_Recv (buf, bufsize, MPI_CHAR, 0, 99, MPI_COMM_WORLD, &status);
      time(&T);
      printf ("****** message %d received at %s", k+1, ctime(&T));
    }
  }

  free (buf);

  MPI_Finalize();
  return EXIT_SUCCESS;
}

-- 
Stephan Mertens @ http://www.uni-magdeburg.de/mertens
Supercomputing in Magdeburg @ http://tina.nat.uni-magdeburg.de