1 /* -*- Mode: C; c-basic-offset:4 ; indent-tabs-mode:nil ; -*- */
3 * (C) 2010 by Argonne National Laboratory.
4 * See COPYRIGHT in top-level directory.
7 * Test of reduce scatter with large data (needed in MPICH to trigger the
10 * Each processor contributes its rank + the index to the reduction,
11 * then receives the ith sum
13 * Can be called with any number of processors.
21 /* Limit the number of error reports */
24 int main( int argc, char **argv )
27 int *sendbuf, *recvbuf, *recvcounts;
28 int size, rank, i, j, idx, mycount, sumval;
32 MTest_Init( &argc, &argv );
33 comm = MPI_COMM_WORLD;
35 MPI_Comm_size( comm, &size );
36 MPI_Comm_rank( comm, &rank );
37 recvcounts = (int *)malloc( size * sizeof(int) );
39 fprintf( stderr, "Could not allocate %d ints for recvcounts\n",
41 MPI_Abort( MPI_COMM_WORLD, 1 );
44 mycount = (1024 * 1024) / size;
45 for (i=0; i<size; i++)
46 recvcounts[i] = mycount;
47 sendbuf = (int *) malloc( mycount * size * sizeof(int) );
49 fprintf( stderr, "Could not allocate %d ints for sendbuf\n",
51 MPI_Abort( MPI_COMM_WORLD, 1 );
55 for (i=0; i<size; i++) {
56 for (j=0; j<mycount; j++) {
57 sendbuf[idx++] = rank + i;
60 recvbuf = (int *)malloc( mycount * sizeof(int) );
62 fprintf( stderr, "Could not allocate %d ints for recvbuf\n",
64 MPI_Abort( MPI_COMM_WORLD, 1 );
67 for (i=0; i<mycount; i++) {
71 MPI_Reduce_scatter( sendbuf, recvbuf, recvcounts, MPI_INT, MPI_SUM, comm );
73 sumval = size * rank + ((size - 1) * size)/2;
74 /* recvbuf should be size * (rank + i) */
75 for (i=0; i<mycount; i++) {
76 if (recvbuf[i] != sumval) {
78 if (err < MAX_ERRORS) {
79 fprintf( stdout, "Did not get expected value for reduce scatter\n" );
80 fprintf( stdout, "[%d] Got recvbuf[%d] = %d expected %d\n",
81 rank, i, recvbuf[i], sumval );
86 #if MTEST_HAVE_MIN_MPI_VERSION(2,2)
87 MPI_Reduce_scatter( MPI_IN_PLACE, sendbuf, recvcounts, MPI_INT, MPI_SUM,
90 sumval = size * rank + ((size - 1) * size)/2;
91 /* recv'ed values for my process should be size * (rank + i) */
92 for (i=0; i<mycount; i++) {
93 if (sendbuf[i] != sumval) {
95 if (err < MAX_ERRORS) {
96 fprintf( stdout, "Did not get expected value for reduce scatter (in place)\n" );
97 fprintf( stdout, "[%d] Got buf[%d] = %d expected %d\n",
98 rank, i, sendbuf[i], sumval );
103 MPI_Comm_set_errhandler(MPI_COMM_WORLD, MPI_ERRORS_RETURN);
104 if (MPI_SUCCESS == MPI_Reduce_scatter(sendbuf, sendbuf, recvcounts, MPI_INT, MPI_SUM, comm))
112 MTest_Finalize( err );