How to use shared global datasets in MPI?

Question

I am an MPI beginner. I have a large array gmat of numbers (type double, dimensions 1x14000000) which is precomputed and stored in a binary file. It will use approximately 100 MB in memory (14000000 x8 bytes /1024 /1024). I want to write a MPI code which will do some computation on this array (for example, multiply all elements of gmat by rank number of the process). This array gmat itself stays constant during run-time. The code is supposed to be something like this

#include <iostream>
#include "mpi.h"
double* gmat;
long int imax;

int main(int argc, char* argv[])
{
void performcomputation(int rank); // this function performs the computation and will be called by all processes

imax=atoi(argv[1]); // user inputs the length of gmat 

MPI::Init();
rank = MPI::COMM_WORLD.Get_rank();
size = MPI::COMM_WORLD.Get_size(); //i will use -np 16 = 4 processors x 4 cores

if rank==0 // read the gmat array using one of the processes
{
gmat = new double[imax];
// read values of gmat from a file
// next line is supposed to broadcast values of gmat to all processes which will use it
MPI::COMM_WORLD.Bcast(&gmat,imax,MPI::DOUBLE,1);     
}

MPI::COMM_WORLD.Barrier();
performcomputation(rank);  
MPI::Finalize();  

return 0;
} 

void performcomputation(int rank)
{
int i;
for (i=0;i <imax; i++)
cout << "the new value is" << gmat[i]*rank << endl; 
}    

My question is when I run this code using 16 processes (-np 16), is the gmat same for all of them ? I mean, will the code use 16 x 100 MB in memory to store gmat for each process or will it use only 100 MB since I have defined gmat to be global ? And I don't want different processes to read gmat separately from the file, since reading so many numbers takes time. What is a better way to do this ? Thanks.

Answer 1

First of all, please do not use the MPI C++ bindings. These were deprecated in MPI-2.2 and then deleted in MPI-3.0 and therefore no longer part of the specification, which means that future MPI implementations are not even required to provide C++ bindings, and that if they do, they will probably diverge in the way the interface looks like.

That said, your code contains a very common mistake:

if rank==0 // read the gmat array using one of the processes
{
   gmat = new double[imax];
   // read values of gmat from a file
   // next line is supposed to broadcast values of gmat to all processes which will use it
   MPI::COMM_WORLD.Bcast(&gmat,imax,MPI::DOUBLE,1);     
}

This won't work as there are four errors here. First, gmat is only allocated at rank 0 and not allocated in the others ranks, which is not what you want. Second, you are giving Bcast the address of the pointer gmat and not the address of the data pointed by it (i.e. you should not use the & operator). You also broadcast from rank 0 but put 1 as the broadcast root argument. But the most important error is that MPI_BCAST is a collective communication call and all ranks are required to call it with the same value of the root argument in order for it to complete successfully. The correct code (using the C bindings instead of the C++ ones) is:

gmat = new double[imax];

if (rank == 0)
{
   // read values of gmat from a file
}
MPI_Bcast(gmat, imax, MPI_DOUBLE, 0, MPI_COMM_WORLD);
//        ^^^^                   ^^^
//        no &                root == 0

Each rank has its own copy of gmat. Initially all values are different (e.g. random or all zeros, depends on the memory allocator). After the broadcast all copies will become identical to the copy of gmat at rank 0. After the call to performcomputation() each copy will be different again since each rank multiplies the elements of gmat with a different number. The answer to your question is: the code will use 100 MiB in each rank, therefore 16 x 100 MiB in total.

MPI deals with distributed memory - processes do not share variables, no matter if they are local or global ones. The only way to share data is to use MPI calls like point-to-point communication (e.g. MPI_SEND / MPI_RECV), collective calls (e.g. MPI_BCAST) or one-sided communication (e.g. MPI_PUT / MPI_GET).