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++) {
53 for (j = 0; j < 3; j++) {
54 for (i = 0; i < 3; i++) {
56 for (k = 0; k < 3; k++) {
57 /* col[i] += cin(i,k) * cout(k,j) */
60 tempcol[i] += cin[offset1] * cout[offset2];
63 for (i = 0; i < 3; i++) {
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[])
83 for (i = 0; i < 9; i++) {
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;
108 case 0: /* Identity */
126 /* Compare a matrix with the known result */
127 static int checkResult(int nmat, int mat[], const char *msg)
129 int n, k, errs = 0, wrank;
130 static int solution[9] = { 0, 1, 0,
135 MPI_Comm_rank(MPI_COMM_WORLD, &wrank);
137 for (n = 0; n < nmat; n++) {
138 for (k = 0; k < 9; k++) {
139 if (mat[k] != solution[k]) {
142 printf("Errors for communicators %s\n", MTestGetIntracommName());
147 printf("[%d]matrix #%d(%s): Expected mat[%d,%d] = %d, got %d\n",
148 wrank, n, msg, k / 3, k % 3, solution[k], mat[k]);
153 /* Advance to the next matrix */
159 int main(int argc, char *argv[])
163 int minsize = 2, count;
167 MPI_Datatype mattype;
170 MTest_Init(&argc, &argv);
172 MPI_Op_create(matmult, 0, &op);
174 /* A single rotation matrix (3x3, stored as 9 consequetive elements) */
175 MPI_Type_contiguous(9, MPI_INT, &mattype);
176 MPI_Type_commit(&mattype);
178 /* Sanity check: test that our routines work properly */
181 buf = (int *) malloc(4 * 9 * sizeof(int));
182 initMat(0, 4, 0, &buf[0]);
183 initMat(1, 4, 0, &buf[9]);
184 initMat(2, 4, 0, &buf[18]);
185 initMat(3, 4, 0, &buf[27]);
186 matmult(&buf[0], &buf[9], &one, &mattype);
187 matmult(&buf[9], &buf[18], &one, &mattype);
188 matmult(&buf[18], &buf[27], &one, &mattype);
189 checkResult(1, &buf[27], "Sanity Check");
193 while (MTestGetIntracommGeneral(&comm, minsize, 1)) {
194 if (comm == MPI_COMM_NULL)
197 MPI_Comm_size(comm, &size);
198 MPI_Comm_rank(comm, &rank);
200 for (count = 1; count < size; count++) {
202 /* Allocate the matrices */
203 buf = (int *) malloc(count * 9 * sizeof(int));
205 MPI_Abort(MPI_COMM_WORLD, 1);
208 bufout = (int *) malloc(count * 9 * sizeof(int));
210 MPI_Abort(MPI_COMM_WORLD, 1);
213 for (i = 0; i < count; i++) {
214 initMat(rank, size, i, &buf[i * 9]);
217 MPI_Allreduce(buf, bufout, count, mattype, op, comm);
218 errs += checkResult(count, bufout, "");
220 /* Try the same test, but using MPI_IN_PLACE */
221 for (i = 0; i < count; i++) {
222 initMat(rank, size, i, &bufout[i * 9]);
224 MPI_Allreduce(MPI_IN_PLACE, bufout, count, mattype, op, comm);
225 errs += checkResult(count, bufout, "IN_PLACE");
230 MTestFreeComm(&comm);
234 MPI_Type_free(&mattype);
236 MTest_Finalize(errs);