Multiprocessing by Message Passing MPI

Example 1.5 More on Integration with MPI Collective Communications

This example demonstrates parallel integration (see Example 1 for detail) via the following MPI library routines. Click on them for details:

      Program Example1_5
c#######################################################################
c#
c# This is an MPI example on parallel integration
c# It demonstrates the use of :
c#
c# * MPI_Init
c# * MPI_Comm_rank
c# * MPI_Comm_size
c# * MPI_Bcast
c# * MPI_Gather
c# * MPI_Scatter
c# * MPI_Finalize
c#
c# Dr. Kadin Tseng
c# Scientific Computing and Visualization
c# Boston University
c# 1998
c#
c#######################################################################
      implicit none
      integer n, p, i, j, ierr, master
      real h, result, a, b, integral, pi

      include "mpif.h"  !! This brings in pre-defined MPI constants, ...
      integer Iam
      real my_result, buf(50), tmp
      parameter (master=0)

c**Starts MPI processes ...

      call MPI_Init(ierr)                               !! starts MPI
      call MPI_Comm_rank(MPI_COMM_WORLD, Iam, ierr)    !! get current process id
      call MPI_Comm_size(MPI_COMM_WORLD, p, ierr)       !! get number of processes

      pi = acos(-1.0)   !!  = 3.14159...
      a = 0.0           !! lower limit of integration
      b = pi*1./2.      !! upper limit of integration
      n = 500
      h = (b-a)/n/p     !! length of increment

      my_result = integral(a,Iam,h,n)

      write(*,"('Process ',i2,' has the partial result of',f10.6)")
     &          Iam,my_result

      call MPI_Gather(my_result, 1, MPI_REAL, buf, 1, MPI_REAL, 0,
     &                  MPI_COMM_WORLD, ierr)

      call MPI_Scatter(buf, 1, MPI_REAL, tmp, 1, MPI_REAL, 0,
     &                  MPI_COMM_WORLD, ierr)
      print *,'Result sent back from buf =', tmp

      if(Iam .eq. master) then
        result = 0.0
        do i=1,p
          result = result + buf(i)
        enddo
        print *,'The result =',result
      endif

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

      stop
      end
      real function integral(a,i,h,n)
      implicit none
      integer n, i, j
      real h, h2, aij, a
      real fct, x

      fct(x) = cos(x)                !! kernel of the integral
      integral = 0.0                 !! initialize integral
      h2 = h/2.
      do j=0,n-1                     !! sum over all "j" integrals
        aij = a + (i*n +j)*h         !! lower limit of "j" integral
        integral = integral + fct(aij+h2)*h
      enddo

      return
      end