matrix multiplication using Mpi_Scatter and Mpi_Gather

Question

I newbie to mpi programming. I was trying to write matrix multiplication. Went through the post MPI Matrix Multiplication with scatter gather about matrix multiplication using scatter and gather routine. I tried modifying the code available on above post as below...

#define N 4
#include <stdio.h>
#include <math.h>
#include <sys/time.h>
#include <stdlib.h>
#include <stddef.h>
#include "mpi.h"


void print_results(char *prompt, int a[N][N]);

int main(int argc, char *argv[])
{
    int i, j, k, rank, size, tag = 99, blksz, sum = 0;
    int a[N][N]={{1,2,3,4},{5,6,7,8},{9,1,2,3},{4,5,6,7,}};
    int b[N][N]={{1,2,3,4},{5,6,7,8},{9,1,2,3},{4,5,6,7,}};
    int c[N][N];
    int aa[N],cc[N];

    MPI_Init(&argc, &argv);
    MPI_Comm_size(MPI_COMM_WORLD, &size);
    MPI_Comm_rank(MPI_COMM_WORLD, &rank);

    //scatter rows of first matrix to different processes     
    MPI_Scatter(a, N*N/size, MPI_INT, aa, N*N/size, MPI_INT,0,MPI_COMM_WORLD);

    //broadcast second matrix to all processes
    MPI_Bcast(b, N*N, MPI_INT, 0, MPI_COMM_WORLD);

    MPI_Barrier(MPI_COMM_WORLD);

          //perform vector multiplication by all processes
          for (i = 0; i < N; i++)
            {
                    for (j = 0; j < N; j++)
                    {
                            sum = sum + aa[j] * b[i][j];                
                    }
                    cc[i] = sum;
                    sum = 0;
            }

    MPI_Gather(cc, N*N/size, MPI_INT, c, N*N/size, MPI_INT, 0, MPI_COMM_WORLD);

    MPI_Barrier(MPI_COMM_WORLD);        
    MPI_Finalize();
    print_results("C = ", c);
}

void print_results(char *prompt, int a[N][N])
{
    int i, j;

    printf ("\n\n%s\n", prompt);
    for (i = 0; i < N; i++) {
            for (j = 0; j < N; j++) {
                    printf(" %d", a[i][j]);
            }
            printf ("\n");
    }
    printf ("\n\n");
}

I ran above program as

$mpirun -np 4 ./a.out

For above program I am getting following incorrect output..

C = 
 0 0 -562242168 32766
 1 0 4197933 0
 -562242176 32766 0 0
 4197856 0 4196672 0

C = 
 0 0 -1064802792 32765
 1 0 4197933 0
 -1064802800 32765 0 0
 4197856 0 4196672 0

C = 
 30 70 29 60
 70 174 89 148
 29 89 95 74
 60 148 74 126

C = 
 0 0 -1845552920 32765
 1 0 4197933 0
 -1845552928 32765 0 0
 4197856 0 4196672 0

I have following queries 1. Why result matrix C is getting printed by all processes. It is supposed to be printed by only main process. 2. Why incorrect result is being printed?

Corrections and help in this regard will be appreciated.

Solution

The result matrix c is getting printed by all processes because every process executes the function void print_results(char *prompt, int a[N][N]). Since you are gathering at the process having rank 0, add a statement if (rank == 0) before calling the print_results(...) function. Further, the result is incorrect because of a wrong loop logic in :

                for (j = 0; j < N; j++)
                {
                        sum = sum + aa[j] * b[i][j];                
                }

This should be :

                for (j = 0; j < N; j++)
                {
                        sum = sum + aa[j] * b[j][i];                
                }

Also there is no need to broadcast b as all processes already already have a copy of it and you can avoid MPI_Barrier(). The complete program then becomes :

#define N 4
#include <stdio.h>
#include <math.h>
#include <sys/time.h>
#include <stdlib.h>
#include <stddef.h>
#include "mpi.h"


void print_results(char *prompt, int a[N][N]);

int main(int argc, char *argv[])
{
    int i, j, k, rank, size, tag = 99, blksz, sum = 0;
    int a[N][N]={{1,2,3,4},{5,6,7,8},{9,1,2,3},{4,5,6,7,}};
    int b[N][N]={{1,2,3,4},{5,6,7,8},{9,1,2,3},{4,5,6,7,}};
    int c[N][N];
    int aa[N],cc[N];

    MPI_Init(&argc, &argv);
    MPI_Comm_size(MPI_COMM_WORLD, &size);
    MPI_Comm_rank(MPI_COMM_WORLD, &rank);

    //scatter rows of first matrix to different processes     
    MPI_Scatter(a, N*N/size, MPI_INT, aa, N*N/size, MPI_INT,0,MPI_COMM_WORLD);

    //broadcast second matrix to all processes
    MPI_Bcast(b, N*N, MPI_INT, 0, MPI_COMM_WORLD);

    MPI_Barrier(MPI_COMM_WORLD);

          //perform vector multiplication by all processes
          for (i = 0; i < N; i++)
            {
                    for (j = 0; j < N; j++)
                    {
                            sum = sum + aa[j] * b[j][i];  //MISTAKE_WAS_HERE               
                    }
                    cc[i] = sum;
                    sum = 0;
            }

    MPI_Gather(cc, N*N/size, MPI_INT, c, N*N/size, MPI_INT, 0, MPI_COMM_WORLD);

    MPI_Barrier(MPI_COMM_WORLD);        
    MPI_Finalize();
    if (rank == 0)                         //I_ADDED_THIS
        print_results("C = ", c);
}

void print_results(char *prompt, int a[N][N])
{
    int i, j;

    printf ("\n\n%s\n", prompt);
    for (i = 0; i < N; i++) {
            for (j = 0; j < N; j++) {
                    printf(" %d", a[i][j]);
            }
            printf ("\n");
    }
    printf ("\n\n");
}

Then c =

C = 
 54 37 47 57

 130 93 119 145

 44 41 56 71

 111 79 101 123 

Comments

[Gilles Gouaillardet]

If you choose not to broadcast b, then you have the option not to scatter a (and declare int *aa and have it point to the right row of a instead). Also, c only has to be declared/allocated on rank 0.