1 /* -*- Mode: C; c-basic-offset:4 ; indent-tabs-mode:nil ; -*- */
3 * (C) 2004 by Argonne National Laboratory.
4 * See COPYRIGHT in top-level directory.
13 static char MTEST_Descrip[] = "Test MPI_Allreduce with non-commutative user-defined operations using matrix rotations";
16 /* This example is similar to allred3.c, but uses only 3x3 matrics with
17 integer-valued entries. This is an associative but not commutative
19 The number of matrices is the count argument. The matrix is stored
23 Three different matrices are used:
31 I^k A I^(p-2-k-j) B I^j
39 for all values of k, p, and j.
42 void matmult( void *cinPtr, void *coutPtr, int *count, MPI_Datatype *dtype );
44 void matmult( void *cinPtr, void *coutPtr, int *count, MPI_Datatype *dtype )
46 const int *cin = (const int *)cinPtr;
47 int *cout = (int *)coutPtr;
52 for (nmat = 0; nmat < *count; nmat++) {
57 /* col[i] += cin(i,k) * cout(k,j) */
60 tempcol[i] += cin[offset1] * cout[offset2];
65 cout[offset1] = tempcol[i];
68 /* Advance to the next matrix */
74 /* Initialize the integer matrix as one of the
75 above matrix entries, as a function of count.
76 We guarantee that both the A and B matrices are included.
78 static void initMat( int rank, int size, int nmat, int mat[] )
87 /* Decide which matrix to create (I, A, or B) */
89 /* rank 0 is A, 1 is B */
94 /* Most ranks are identity matrices */
96 /* Make sure exactly one rank gets the A matrix
97 and one the B matrix */
99 tmpB = (3 * size) / 4;
101 if (rank == tmpA) kind = 1;
102 if (rank == tmpB) kind = 2;
106 case 0: /* Identity */
124 /* Compare a matrix with the known result */
125 static int checkResult( int nmat, int mat[], const char *msg )
127 int n, k, errs = 0, wrank;
128 static int solution[9] = { 0, 1, 0,
132 MPI_Comm_rank( MPI_COMM_WORLD, &wrank );
134 for (n=0; n<nmat; n++) {
135 for (k=0; k<9; k++) {
136 if (mat[k] != solution[k]) {
139 printf( "Errors for communicators %s\n",
140 MTestGetIntracommName() ); fflush(stdout);
144 printf( "[%d]matrix #%d(%s): Expected mat[%d,%d] = %d, got %d\n",
145 wrank, n, msg, k / 3, k % 3, solution[k], mat[k] );
150 /* Advance to the next matrix */
156 int main( int argc, char *argv[] )
160 int minsize = 2, count;
164 MPI_Datatype mattype;
167 MTest_Init( &argc, &argv );
169 MPI_Op_create( matmult, 0, &op );
171 /* A single rotation matrix (3x3, stored as 9 consequetive elements) */
172 MPI_Type_contiguous( 9, MPI_INT, &mattype );
173 MPI_Type_commit( &mattype );
175 /* Sanity check: test that our routines work properly */
177 buf = (int *)malloc( 4*9 * sizeof(int) );
178 initMat( 0, 4, 0, &buf[0] );
179 initMat( 1, 4, 0, &buf[9] );
180 initMat( 2, 4, 0, &buf[18] );
181 initMat( 3, 4, 0, &buf[27] );
182 matmult( &buf[0], &buf[9], &one, &mattype );
183 matmult( &buf[9], &buf[18], &one, &mattype );
184 matmult( &buf[18], &buf[27], &one, &mattype );
185 checkResult( 1, &buf[27], "Sanity Check" );
189 while (MTestGetIntracommGeneral( &comm, minsize, 1 )) {
190 if (comm == MPI_COMM_NULL) continue;
192 MPI_Comm_size( comm, &size );
193 MPI_Comm_rank( comm, &rank );
195 for (count = 1; count < size; count ++ ) {
197 /* Allocate the matrices */
198 buf = (int *)malloc( count * 9 * sizeof(int) );
200 MPI_Abort( MPI_COMM_WORLD, 1 );
204 bufout = (int *)malloc( count * 9 * sizeof(int) );
206 MPI_Abort( MPI_COMM_WORLD, 1 );
210 for (i=0; i < count; i++) {
211 initMat( rank, size, i, &buf[i*9] );
214 MPI_Allreduce( buf, bufout, count, mattype, op, comm );
215 errs += checkResult( count, bufout, "" );
217 /* Try the same test, but using MPI_IN_PLACE */
218 for (i=0; i < count; i++) {
219 initMat( rank, size, i, &bufout[i*9] );
221 MPI_Allreduce( MPI_IN_PLACE, bufout, count, mattype, op, comm );
222 errs += checkResult( count, bufout, "IN_PLACE" );
227 MTestFreeComm( &comm );
231 MPI_Type_free( &mattype );
233 MTest_Finalize( errs );