Alibemz
Alibemz

Reputation: 115

mpi_allreduce sum over a derived datatype vector

I'm trying to reduce (sum) a derived data type which is created by MPI_type_vector. When I'm running the code it crashes and complains about the reduction MPI_SUM is not defined for non-instrinsic datatype. I wrote a peice of simple code to show my problem. The code tries to reduce the diagonal elements of a 3*3 matrix:

#include "mpi.h"
#include <stdio.h>

int main(int argc, char *argv[]) {
int rank, size, i, j;
double a[3][3] ;

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

MPI_Datatype diag3;
MPI_Type_vector(3,1,4,MPI_DOUBLE,&diag3);
MPI_Type_commit(&diag3);

if(rank==0)
    for(i=0; i < 3 ; i++)
       for(j=0; j < 3 ; j++)
           a[i][j]=1;

if(rank==1)
    for(i=0; i < 3 ; i++)
       for(j=0; j < 3 ; j++)
           a[i][j]=-1;

MPI_Allreduce( MPI_IN_PLACE, &a[0][0], 1, diag3, MPI_SUM, MPI_COMM_WORLD );


for(i=0; i < 3 ; i++)
       for(j=0; j < 3 ; j++)
           printf("rank=%d\ta[%d][%d]=%f\n",rank,i,j,a[i][j]);

MPI_Finalize();

}

The error after running is something like this:

*** An error occurred in MPI_Allreduce: the reduction operation MPI_SUM is not defined for non-intrinsic datatypes
*** reported by process [140130307538945,1]
*** on communicator MPI_COMM_WORLD
*** MPI_ERR_OP: invalid reduce operation

I thought the Reduce and MPI_SUM could be performed on a derived datatype as the MPI documentations say. So, what's the problem in the code?

Upvotes: 4

Views: 1819

Answers (2)

Jonathan Dursi
Jonathan Dursi

Reputation: 50927

Ed Smith is right, you'll need to define your own operation; but it needs to be a little more complicated for non-contiguous types than the version he listed. Below we have an add_double_vector function which will decode any double_vector type and operate on it; it's relatively straightforwardly extended to len > 1.

#include "mpi.h"
#include <stdio.h>

void add_double_vector(void *in, void *inout, int *len, MPI_Datatype *dtype)
{
    double *invec = in;
    double *inoutvec = inout;
    int nints, naddresses, ntypes;
    int combiner;

    if (*len != 1) {
        fprintf(stderr,"my_add: len>1 not implemented.\n");
        return;
    } 

    MPI_Type_get_envelope(*dtype, &nints, &naddresses, &ntypes, &combiner); 
    if (combiner != MPI_COMBINER_VECTOR) {
        fprintf(stderr,"my_add: do not understand composite datatype.\n");
        return;
    } 

    int vecargs [nints];
    MPI_Aint vecaddrs[naddresses];
    MPI_Datatype vectypes[ntypes];

    MPI_Type_get_contents(*dtype, nints, naddresses, ntypes, 
            vecargs, vecaddrs, vectypes);

    if (vectypes[0] != MPI_DOUBLE) {
        fprintf(stderr,"my_add: not a vector of DOUBLEs.\n");
    }

    int count    = vecargs[0];
    int blocklen = vecargs[1];
    int stride   = vecargs[2];

    for ( int i=0; i<count; i++ ) {
        for ( int j=0; j<blocklen; j++) {
            inoutvec[i*stride+j] += invec[i*stride+j]; 
        } 
    }
}

int main(int argc, char *argv[]) {
    int rank, size, i, j;
    const int n=3;
    double a[n][n] ;

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

    MPI_Datatype diag3;
    MPI_Type_vector(n,1,n+1,MPI_DOUBLE,&diag3);
    MPI_Type_commit(&diag3);

    if(rank==0)
        for(i=0; i < n ; i++)
            for(j=0; j < n ; j++)
                a[i][j]=1;

    if(rank==1)
        for(i=0; i < n ; i++)
            for(j=0; j < n ; j++)
                a[i][j]=-1;

    MPI_Op vector_add;
    MPI_Op_create( add_double_vector, 1, &vector_add );
    MPI_Allreduce( MPI_IN_PLACE, &a[0][0], 1, diag3, vector_add, MPI_COMM_WORLD );
    MPI_Op_free( &vector_add );

    for(i=0; i < n ; i++)
        for(j=0; j < n ; j++)
            printf("rank=%d\ta[%d][%d]=%f\n",rank,i,j,a[i][j]);

    MPI_Finalize();

}

Compiling and running gives the correct answer:

$ mpicc -o foo foo.c -std=c99
$ mpirun -np 2 ./foo 
rank=1  a[0][0]=0.000000
rank=1  a[0][1]=-1.000000
rank=1  a[0][2]=-1.000000
rank=1  a[1][0]=-1.000000
rank=1  a[1][1]=0.000000
rank=1  a[1][2]=-1.000000
rank=1  a[2][0]=-1.000000
rank=1  a[2][1]=-1.000000
rank=1  a[2][2]=0.000000
rank=0  a[0][0]=0.000000
rank=0  a[0][1]=1.000000
rank=0  a[0][2]=1.000000
rank=0  a[1][0]=1.000000
rank=0  a[1][1]=0.000000
rank=0  a[1][2]=1.000000
rank=0  a[2][0]=1.000000
rank=0  a[2][1]=1.000000
rank=0  a[2][2]=0.000000

Upvotes: 3

Ed Smith
Ed Smith

Reputation: 13206

I think the error is because there is no defined way to add the vectors you have created. If you define you're own sum operation:

#include "mpi.h"
#include <stdio.h>

void mySum ( int *, int *, int *, MPI_Datatype * );

void mySum(int *invec, int *inoutvec, int *len, MPI_Datatype *dtype)
{
    int i;
    for ( i=0; i<*len; i++ ) 
        inoutvec[i] += invec[i];
}
int main(int argc, char *argv[]) {
int rank, size, i, j;
double a[3][3] ;

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

MPI_Datatype diag3;
MPI_Type_vector(3,1,4,MPI_DOUBLE,&diag3);
MPI_Type_commit(&diag3);

MPI_Op diagSum;
MPI_Op_create( (MPI_User_function *)mySum, 1, &diagSum );

if(rank==0)
    for(i=0; i < 3 ; i++)
       for(j=0; j < 3 ; j++)
           a[i][j]=i+j;

if(rank==1)
    for(i=0; i < 3 ; i++)
       for(j=0; j < 3 ; j++)
           a[i][j]=-1;


MPI_Allreduce( MPI_IN_PLACE, &a[0][0], 1, diag3, diagSum, MPI_COMM_WORLD );

for(i=0; i < 3 ; i++)
       for(j=0; j < 3 ; j++)
           printf("rank=%d\ta[%d][%d]=%f\n",rank,i,j,a[i][j]);

MPI_Op_free( &diagSum );
MPI_Finalize();

}

Upvotes: 0

Related Questions