1
votes

*The master task first initializes an array and then distributes an equal portion that array to the other tasks. After the other tasks receive their portion of the array, they perform an addition operation to each array element.They also maintain a sum for their portion of the array. The master task does likewise with its portion of the array. As each of the non-master tasks finish, they send their updated portion of the array to the master. An MPI collective communication call is used to collect the sums maintained by each task. Finally, the master task displays selected parts of the final array and the global sum of all array elements. *

#include "mpi.h"
#include <stdio.h>
#include <stdlib.h>
#define ARRAYSIZE 16000000
#define MASTER

float data[ARRAYSIZE];

int main (int argc, char *argv[])
{
int numtasks, taskid, rc, dest, source, offset, i, j, tag1,
    tag2, chunksize;
float mysum, sum;
float update(int myoffset, int chunk, int myid);
MPI_Status status;

/******Initializations******/
MPI_Init(&argc, &argv);
MPI_Comm_size(MPI_COMM_WORLD, &numtasks);
if (numtasks % 4 != 0) {
    printf("Quitting. Number of MPI tasks must be divisible by 4. \n");
    MPI_Abort(MPI_COMM_WORLD, rc);
    exit(0);
    }
MPI_Comm_rank(MPI_COMM_WORLD, &taskid);
chunksize = ARRAYSIZE / numtasks;
tag2 = 1;
tag1 = 2;

/******Master task only ******/
if (taskid == MASTER){
   /* Initialize the array */
    sum = 0;
    for (i = 0; i < ARRAYSIZE; i++){
         data[i] = i * 1.0;
         sum = sum + data[i];
      }
printf("Initialized array sum = %e\n", sum);

/* Send each task its portion of the array - mater keeps 1st part */
offset = chunksize;
for (dest = 1; dest < numtasks; dest++){
     MPI_Send(&offset, 1, MPI_INT, dest, tag1, MPI_COMM_WORLD);
     MPI_Send(&data[offset], chunksize, MPI_FLOAT, dest, tag2, MPI_COMM_WRLD);
     printf("Sent %d elements to task %d offset = %d\n, chunksize, dest, offset);
     offset = offset + chunksize;
  }

/* Master does its part of the work */
offset = 0;
mysum = update(offset, chunksize, taskid);

/* Get final sum */
MPI_Reduce(&mysum, &sum, 1, MPI_FLOAT, MPI_SUM, MASTER, MPI_COMM_WORLD);
printf("***Final sum = %e ***\n", sum);
} /* end of master section */

/******Non-master tasks only ******/
if (taskid > MASTER){

   /* Receive my portion of array from the master task */
   source = MASTER;
   MPI_Recv(&offset, 1, MPI_INT, source, tag1, MPI_COMM_WORLD, &status);
   MPI_Recv(&data[offset], chunksize, MPI_FLOAT, source, tag2, MPI_COMM_WORLD, &status);
   mysum = update(offset, chunksize, taskid);
   MPI_Reduce(&mysum, &sum, 1, MPI_FLOAT, MPI_SUM, MASTER, MPI_COMM_WORLD);
   }  /* end of non-master */

MPI_Finalize();
} /* end of main */

float update(int myoffset, int chunk, int myid){
int i;
float mysum;
/* Perform addition to each of my array elements and keep my sum */
mysum = 0;
for (i = myoffset; i < myoffset + chunk; i++){
     mysum = mysum + data[i];
  }
printf("Task %d mysum = %e\n", myid, mysum);
return mysum;
}


/******The result of this program is: ******/
MPI task 0 has started...
MPI task 1 has started...
MPI task 2 has started...
MPI task 3 has started...
Initialized array sum = 1.335708e+14
Sent 4000000 elements to task 1 offset= 4000000
Sent 4000000 elements to task 2 offset= 8000000
Task 1 mysum = 2.442024e+13
Sent 4000000 elements to task 3 offset= 12000000
Task 2 mysum = 3.991501e+13
Task 3 mysum = 5.809336e+13
Task 0 mysum = 7.994294e+12
Sample results: 
   0.000000e+00  1.000000e+00  2.000000e+00  3.000000e+00  4.000000e+00
   4.000000e+06  4.000001e+06  4.000002e+06  4.000003e+06  4.000004e+06
   8.000000e+06  8.000001e+06  8.000002e+06  8.000003e+06  8.000004e+06
   1.200000e+07  1.200000e+07  1.200000e+07  1.200000e+07  1.200000e+07
*** Final sum= 1.304229e+14 ***

*So my question is why these two sum don't hold the same value**

1

1 Answers

0
votes

You are storing the result in a 32-bit floating-point number (i.e. a float) which simply isn't enough to maintain all the accuracy you need. What you are seeing is a classic example of how rounding errors accumulate differently depending on what order you add numbers together.

If you just replace all your floats by doubles then it is OK:

mpiexec -n 4 ./arraysum
Initialized array sum = 1.280000e+14
Sent 4000000 elements to task 1 offset = 4000000
Task 1 mysum = 2.400000e+13
Sent 4000000 elements to task 2 offset = 8000000
Task 2 mysum = 4.000000e+13
Sent 4000000 elements to task 3 offset = 12000000
Task 0 mysum = 7.999998e+12
Task 3 mysum = 5.600000e+13
***Final sum = 1.280000e+14 ***