c
c compute pi by integrating f(x) = 4/(1 + x**2)
c
c each process:
c - receives the # of intervals used in the approximation
c - calculates the areas of it's rectangles
c - synchronizes for a global summation
c process 0 prints the result and the time it took
c
program main
include 'mpif.h'
double precision PIX
parameter (PIX = 4*atan(1.0))
double precision mypi, pi, h, sum, x, f, a
integer n, myid, numprocs, i, ierr
c
c function to integrate
c
f(a) = 4.d0 / (1.d0 + a*a)
call MPI_INIT(ierr)
call MPI_COMM_RANK(MPI_COMM_WORLD, myid, ierr)
call MPI_COMM_SIZE(MPI_COMM_WORLD, numprocs, ierr)
print *, "Process ", myid, " of ", numprocs, " is alive"
if (myid .eq. 0) then
print *,"Enter the number of intervals: (0 to quit)"
read(5,*) n
endif
call MPI_BCAST(n, 1, MPI_INTEGER, 0, MPI_COMM_WORLD, ierr)
c check for n > 0
IF (N.GT.0) THEN
c
c calculate the interval size
c
h = 1.0d0 / n
sum = 0.0d0
do 20 i = myid + 1, n, numprocs
x = h * (dble(i) - 0.5d0)
sum = sum + f(x)
20 continue
mypi = h * sum
c
c collect all the partial sums
c
call MPI_REDUCE(mypi, pi, 1, MPI_DOUBLE_PRECISION,
+ MPI_SUM, 0, MPI_COMM_WORLD, ierr)
c
c process 0 prints the result
c
if (myid .eq. 0) then
write(6, 97) pi, abs(pi - PIX)
97 format(' pi is approximately: ', F18.16,
+ ' Error is: ', F18.16)
endif
ENDIF
call MPI_FINALIZE(ierr)
stop
end