Scientific Computing & Visualization | Documentation | Tutorials | Multiprocessing by Message Passing MPI

This Site

Multiprocessing by Message Passing MPI

Example 1.2 Numerical Integration with MPI Blocking Send/Receive

The major drawback of example 1.1 is the potential communication deadlock stemming from sending and receiving a message on processor 0 concurrently. Since my_int is computed on processor 0, there is no need to perform message passing for that. The elimination of message passing on processor 0 leads to a code that is no longer prone to deadlock.

For detail explanations of functions used in this example, please read MPI: The Complete Reference.

Example 1.2 Fortran code

Example 1.1 Fortran Code

      Program Example1_2
c##################################################################################
c#                                                                                #
c# This is an MPI example on parallel integration to demonstrate the use of:      #
c#                                                                                #
c# * MPI_Init, MPI_Comm_rank, MPI_Comm_size, MPI_Finalize                         #
c# * MPI_Recv, MPI_Send                                                           #
c#                                                                                #
c# Dr. Kadin Tseng                                                                #
c# Scientific Computing and Visualization                                         #
c# Boston University                                                              #
c# 1998                                                                           #
c#                                                                                #
c##################################################################################
      implicit none
      integer n, p, i, j, proc, ierr, master, myid, tag, comm
      real h, a, b, integral, pi, ai, my_int, integral_sum
      include "mpif.h"  ! brings in pre-defined MPI constants, ...
      integer status(MPI_STATUS_SIZE)  ! size defined in mpif.h
      data master/0/    ! processor 0 collects integral sums from other processors

      comm = MPI_COMM_WORLD      
      call MPI_Init(ierr)                        ! starts MPI
      call MPI_Comm_rank(comm, myid, ierr)      ! get current proc ID
      call MPI_Comm_size(comm, p, ierr)         ! get number of procs

      pi = acos(-1.0)   !  = 3.14159...
      a = 0.0           ! lower limit of integration
      b = pi/2.         ! upper limit of integration
      n = 500           ! number of increment within each process
      tag = 123         ! set the tag to identify this particular job
      h = (b-a)/n/p     ! length of increment

      ai = a + myid*n*h ! lower limit of integration for partition myid
      my_int = integral(ai, h, n) 
      write(*,"('Process ',i2,' has the partial sum of',f10.6)")
     &  myid,my_int

      if(myid .eq. master) then      ! do following only on master ...
        integral_sum = my_int           !starts with value on master 
        do proc=1,p-1   ! loop on procs to collect local sums, exclude master
          call MPI_Recv( 
     &       my_int, 1, MPI_REAL, 
     &       proc,        ! message source
     &       tag,         ! message tag
     &       comm, status, ierr)    ! status reports source, tag
          integral_sum = integral_sum + my_int   ! sum my_int from processors
        enddo
        print *,'The integral =',integral_sum
      else
        call MPI_Send(  
     &     my_int, 1, MPI_REAL,   ! buffer, size, datatype
     &     master,     ! where to send message
     &     tag,         ! message tag
     &     comm, ierr)
      endif

      call MPI_Finalize(ierr)                          ! MPI finish up ...

      end
      real function integral(ai, h, n)
      implicit none
      integer n, j
      real h, ai, aij

      integral = 0.0                ! initialize integral
      do j=0,n-1                    ! sum integrals
        aij = ai +(j+0.5)*h         ! abscissa mid-point
        integral = integral + cos(aij)*h
      enddo

      return
      end

Example 1.2 C code

Example 1.1 C code

#include <mpi.h>
#include <math.h>
#include <stdio.h>
/* Prototype */
float integral(float ai, float h, int n);
void main(int argc, char* argv[])
{
/*###############################################################################
 #                                                                              #
 # This is an MPI example on parallel integration to demonstrate the use of:    #
 #                                                                              #
 # * MPI_Init, MPI_Comm_rank, MPI_Comm_size, MPI_Finalize                       #
 # * MPI_Recv, MPI_Send                                                         #
 #                                                                              #
 # Dr. Kadin Tseng                                                              #
 # Scientific Computing and Visualization                                       #
 # Boston University                                                            #
 # 1998                                                                         #
 #                                                                              #
 ###############################################################################*/
      int n, p, myid, tag, proc, ierr;
      float h, integral_sum, a, b, ai, pi, my_int;
      int master = 0;  /* processor performing total sum */
      MPI_Comm comm;
      MPI_Status status;

      comm = MPI_COMM_WORLD;      
      ierr = MPI_Init(&argc,&argv);           /* starts MPI */
      MPI_Comm_rank(comm, &myid);            /* get current process id */
      MPI_Comm_size(comm, &p);               /* get number of processes */

      pi = acos(-1.0);  /* = 3.14159... */
      a = 0.;           /* lower limit of integration */
      b = pi*1./2.;     /* upper limit of integration */
      n = 500;          /* number of increment within each process */
      tag = 123;        /* set the tag to identify this particular job */
      h = (b-a)/n/p;    /* length of increment */

      ai = a + myid*n*h;  /* lower limit of integration for partition myid */
      my_int = integral(ai, h, n);   /* 0<=myid<=p-1 */

      printf("Process %d has the partial integral of %f\n", myid,my_int);

      MPI_Send( 
              &my_int, 1, MPI_FLOAT, 
              master,        /* message destination */
              tag,           /* message tag */
              comm);

      if(myid == master) {  /* Receives serialized */
        integral_sum = 0.0;
        for (proc=0;proc<p;proc++) {
          MPI_Recv( 
                &my_int, 1, MPI_FLOAT, 
                proc,        /* message source */
                tag,         /* message tag */
                comm, &status);     /* status reports source, tag */
          integral_sum += my_int;
        }
        printf("The Integral =%f\n",integral_sum); /* sum of my_int */
      }
      MPI_Finalize();                        /* let MPI finish up ... */
}
float integral(float ai, float h, int n)
{
      int j;
      float aij, integ;

      integ = 0.0;                 /* initialize */
      for (j=0;j<n;j++) {          /* sum integrals */
        aij = ai + (j+0.5)*h;      /* mid-point */
        integ += cos(aij)*h;
      }
      return integ;
}

Discussion

In the next example, a non-blocking point-to-point send routine is introduced to replace the blocking send. This enables any work that follows the send to proceed while the send is pending or in progress.

OIT | CCS | June 26, 2009