1 /* -*- Mode: C; c-basic-offset:4 ; indent-tabs-mode:nil ; -*- */
2 /*
3 * (C) 2012 by Argonne National Laboratory.
4 * See COPYRIGHT in top-level directory.
5 */
6 #include <stdio.h>
7 #include <stdlib.h>
8 #include <string.h>
9
10 #include <mpi.h>
11 #include "mpitest.h"
12
13 #if !defined(USE_STRICT_MPI) && defined(MPICH2)
14 #define TEST_NEIGHB_COLL 1
15 #endif
16
17 /* assert-like macro that bumps the err count and emits a message */
18 #define check(x_) \
19 do { \
20 if (!(x_)) { \
21 ++errs; \
22 if (errs < 10) { \
23 fprintf(stderr, "check failed: (%s), line %d\n", #x_, __LINE__); \
24 } \
25 } \
26 } while (0)
27
main(int argc,char * argv[])28 int main(int argc, char *argv[])
29 {
30 int errs = 0;
31 int wrank, wsize;
32 int periods[1] = { 0 };
33 MPI_Comm cart, dgraph, graph;
34
35 MPI_Init(&argc, &argv);
36 MPI_Comm_rank(MPI_COMM_WORLD, &wrank);
37 MPI_Comm_size(MPI_COMM_WORLD, &wsize);
38
39 #if defined(TEST_NEIGHB_COLL)
40 /* a basic test for the 10 (5 patterns x {blocking,nonblocking}) MPI-3
41 * neighborhood collective routines */
42
43 /* (wrap)--> 0 <--> 1 <--> ... <--> p-1 <--(wrap) */
44 MPI_Cart_create(MPI_COMM_WORLD, 1, &wsize, periods, /*reorder=*/0, &cart);
45
46 /* allgather */
47 {
48 int sendbuf[1] = { wrank };
49 int recvbuf[2] = { 0xdeadbeef, 0xdeadbeef };
50
51 /* should see one send to each neighbor (rank-1 and rank+1) and one receive
52 * each from same */
53 MPIX_Neighbor_allgather(sendbuf, 1, MPI_INT, recvbuf, 1, MPI_INT, cart);
54
55 if (wrank == 0)
56 check(recvbuf[0] == 0xdeadbeef);
57 else
58 check(recvbuf[0] == wrank - 1);
59
60 if (wrank == wsize - 1)
61 check(recvbuf[1] == 0xdeadbeef);
62 else
63 check(recvbuf[1] == wrank + 1);
64 }
65
66 /* allgatherv */
67 {
68 int sendbuf[1] = { wrank };
69 int recvbuf[2] = { 0xdeadbeef, 0xdeadbeef };
70 int recvcounts[2] = { 1, 1 };
71 int displs[2] = { 1, 0};
72
73 /* should see one send to each neighbor (rank-1 and rank+1) and one receive
74 * each from same, but put them in opposite slots in the buffer */
75 MPIX_Neighbor_allgatherv(sendbuf, 1, MPI_INT, recvbuf, recvcounts, displs, MPI_INT, cart);
76
77 if (wrank == 0)
78 check(recvbuf[1] == 0xdeadbeef);
79 else
80 check(recvbuf[1] == wrank - 1);
81
82 if (wrank == wsize - 1)
83 check(recvbuf[0] == 0xdeadbeef);
84 else
85 check(recvbuf[0] == wrank + 1);
86 }
87
88 /* alltoall */
89 {
90 int sendbuf[2] = { -(wrank+1), wrank+1 };
91 int recvbuf[2] = { 0xdeadbeef, 0xdeadbeef };
92
93 /* should see one send to each neighbor (rank-1 and rank+1) and one
94 * receive each from same */
95 MPIX_Neighbor_alltoall(sendbuf, 1, MPI_INT, recvbuf, 1, MPI_INT, cart);
96
97 if (wrank == 0)
98 check(recvbuf[0] == 0xdeadbeef);
99 else
100 check(recvbuf[0] == wrank);
101
102 if (wrank == wsize - 1)
103 check(recvbuf[1] == 0xdeadbeef);
104 else
105 check(recvbuf[1] == -(wrank + 2));
106 }
107
108 /* alltoallv */
109 {
110 int sendbuf[2] = { -(wrank+1), wrank+1 };
111 int recvbuf[2] = { 0xdeadbeef, 0xdeadbeef };
112 int sendcounts[2] = { 1, 1 };
113 int recvcounts[2] = { 1, 1 };
114 int sdispls[2] = { 0, 1 };
115 int rdispls[2] = { 1, 0 };
116
117 /* should see one send to each neighbor (rank-1 and rank+1) and one receive
118 * each from same, but put them in opposite slots in the buffer */
119 MPIX_Neighbor_alltoallv(sendbuf, sendcounts, sdispls, MPI_INT,
120 recvbuf, recvcounts, rdispls, MPI_INT,
121 cart);
122
123 if (wrank == 0)
124 check(recvbuf[1] == 0xdeadbeef);
125 else
126 check(recvbuf[1] == wrank);
127
128 if (wrank == wsize - 1)
129 check(recvbuf[0] == 0xdeadbeef);
130 else
131 check(recvbuf[0] == -(wrank + 2));
132 }
133
134 /* alltoallw */
135 {
136 int sendbuf[2] = { -(wrank+1), wrank+1 };
137 int recvbuf[2] = { 0xdeadbeef, 0xdeadbeef };
138 int sendcounts[2] = { 1, 1 };
139 int recvcounts[2] = { 1, 1 };
140 MPI_Aint sdispls[2] = { 0, sizeof(int) };
141 MPI_Aint rdispls[2] = { sizeof(int), 0 };
142 MPI_Datatype sendtypes[2] = { MPI_INT, MPI_INT };
143 MPI_Datatype recvtypes[2] = { MPI_INT, MPI_INT };
144
145 /* should see one send to each neighbor (rank-1 and rank+1) and one receive
146 * each from same, but put them in opposite slots in the buffer */
147 MPIX_Neighbor_alltoallw(sendbuf, sendcounts, sdispls, sendtypes,
148 recvbuf, recvcounts, rdispls, recvtypes,
149 cart);
150
151 if (wrank == 0)
152 check(recvbuf[1] == 0xdeadbeef);
153 else
154 check(recvbuf[1] == wrank);
155
156 if (wrank == wsize - 1)
157 check(recvbuf[0] == 0xdeadbeef);
158 else
159 check(recvbuf[0] == -(wrank + 2));
160 }
161
162
163 MPI_Comm_free(&cart);
164 #endif /* defined(TEST_NEIGHB_COLL) */
165
166 MPI_Reduce((wrank == 0 ? MPI_IN_PLACE : &errs), &errs, 1, MPI_INT, MPI_SUM, 0, MPI_COMM_WORLD);
167 if (wrank == 0) {
168 if (errs) {
169 printf("found %d errors\n", errs);
170 }
171 else {
172 printf(" No errors\n");
173 }
174 }
175 MPI_Finalize();
176
177 return 0;
178 }
179
180