Numerical Integration : TechWeb : Boston University

The five basic examples in this group are all based on a single application: numerical integration. Numerical integration is chosen because it is trivially parallelizable and at the same time a problem that is very narrowly focused. The same application is employed in each example; different MPI library routines, however, are used. In addition, disadvantages with the use of specific MPI functionalities are pointed out in each example and are remedied in the next example. Hence, improvements in implementation proficiency can be seen from example to example.

If you have not deal with numerical integration before, please spend a few minutes getting familiar with the integration procedure illustrated below. If you are familiar with it, please spend the necessary time to be acquainted with the terminologies and variable names used before going to the parallel application examples.

Example 1. Numerical Integration

The example application is to integrate cos(x) from a to b numerically. There are various ways to perform numerical integrations of this type. Among them, the Mid-point rule is the least accurate but is chosen nevertheless for its simplicity. Essentially, the integrand, cos(x), is assumed to be constant within the upper and lower limit of integration and is taken to be the mid-point between the limits.

Normally, the integration range need only be divided into a series of smaller intervals so that the mid-point rule can be applied. Here, the integration range is first divided into a number of “partitions”, each of which is assigned to a processor. Each processor will then subdivide its own sub-range into smaller segments for mid-point rule integration. The final integral sum is obtained by adding the integral sums from all processors. For single processor, the number of partitions is set to unity and the sub-range is the full range, from a to b.

Example 1. Fortran Code

     Program Example1

c#######################################################################
c#                                                                     #
c# This is an MPI example on parallel integration                      #
c# Demonstrate here is the serial program.                             #
c#                                                                     #
c# Dr. Kadin Tseng                                                     #
c# Scientific Computing and Visualization                              #
c# Boston University                                                   #
c# 1998                                                                #
c#                                                                     #
c#######################################################################
      implicit none
      integer n, p, i, j
      real h, integral_sum, a, b, integral, pi

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

      integral_sum = 0.0      ! stores answer to the integral
      do i=0,p-1        ! sum of integrals over all processes
        ai = a + i*n*h  ! lower limit of integration for partition i
        integral_sum = integral_sum + integral(ai,h,n)
      enddo

      print *,'The integral_sum =', integral_sum

      stop
      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. C code

#include <mpi.h>
#include <math.h>
#include <stdio.h>
float integral(float ai, float h, int n);
void main(void)
{
/***********************************************************************
 *                                                                     *
 * This is a serial C-version of the MPI example on integration        *
 *                                                                     *
 * Dr. Kadin Tseng                                                     *
 * Scientific Computing and Visualization                              *
 * Boston University                                                   *
 * 1998                                                                *
 *                                                                     *
 ***********************************************************************/

      int n, p, i, j, ierr;
      float h, integral_sum, a, b, pi, my_int;

      pi = acos(-1.0);  /* = 3.14159... */
      a = 0.;           /* lower limit of integration */
      b = pi*1./2.;     /* upper limit of integration */
      p = 4;            /* number of processes (partitions) */
      n = 500;          /* number of increment within each process */
      h = (b-a)/n/p;    /* length of increment */

      integral_sum = 0.0;
/* sum of integrals over all processes */
      for (i=0; i<p; i++)  {
        ai = a + i*n*h; /* lower limit of integration for partition i */
        integral_sum += integral(ai,h,n);
      }

      printf("The integral sum =%fn",integral_sum);
}
float integral(float ai, float h, int n)
{
      int j;
      float aij, integ;

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