1 /* -*- Mode: C; c-basic-offset:4 ; indent-tabs-mode:nil ; -*- */
4 * (C) 2003 by Argonne National Laboratory.
5 * See COPYRIGHT in top-level directory.
16 #define MAX_WEIGHT 100
18 /* convenience globals */
21 /* We need MPI 2.2 to be able to compile the following routines. */
22 #if MTEST_HAVE_MIN_MPI_VERSION(2,2)
24 /* Maybe use a bit vector instead? */
27 #define MAX_LAYOUT_NAME_LEN 256
28 char graph_layout_name[MAX_LAYOUT_NAME_LEN] = { '\0' };
30 static void create_graph_layout(int graph_num)
37 strncpy(graph_layout_name, "deterministic complete graph", MAX_LAYOUT_NAME_LEN);
38 for (i = 0; i < size; i++)
39 for (j = 0; j < size; j++)
40 layout[i][j] = (i + 2) * (j + 1);
43 strncpy(graph_layout_name, "every other edge deleted", MAX_LAYOUT_NAME_LEN);
44 for (i = 0; i < size; i++)
45 for (j = 0; j < size; j++)
46 layout[i][j] = (j % 2 ? (i + 2) * (j + 1) : 0);
49 strncpy(graph_layout_name, "only self-edges", MAX_LAYOUT_NAME_LEN);
50 for (i = 0; i < size; i++) {
51 for (j = 0; j < size; j++) {
52 if (i == rank && j == rank)
53 layout[i][j] = 10 * (i + 1);
60 strncpy(graph_layout_name, "no edges", MAX_LAYOUT_NAME_LEN);
61 for (i = 0; i < size; i++)
62 for (j = 0; j < size; j++)
66 strncpy(graph_layout_name, "a random incomplete graph", MAX_LAYOUT_NAME_LEN);
69 /* Create a connectivity graph; layout[i,j]==w represents an outward
70 * connectivity from i to j with weight w, w==0 is no edge. */
71 for (i = 0; i < size; i++) {
72 for (j = 0; j < size; j++) {
73 /* disable about a third of the edges */
74 if (((rand() * 1.0) / RAND_MAX) < 0.33)
77 layout[i][j] = rand() % MAX_WEIGHT;
84 /* because of the randomization we must determine the graph on rank 0 and
85 * send the layout to all other processes */
86 MPI_Bcast(graph_layout_name, MAX_LAYOUT_NAME_LEN, MPI_CHAR, 0, MPI_COMM_WORLD);
87 for (i = 0; i < size; ++i) {
88 MPI_Bcast(layout[i], size, MPI_INT, 0, MPI_COMM_WORLD);
92 static int verify_comm(MPI_Comm comm)
96 int indegree, outdegree, weighted;
97 int *sources, *sweights, *destinations, *dweights;
99 int topo_type = MPI_UNDEFINED;
100 MPI_Comm dupcomm = MPI_COMM_NULL;
102 sources = (int *) malloc(size * sizeof(int));
103 sweights = (int *) malloc(size * sizeof(int));
104 destinations = (int *) malloc(size * sizeof(int));
105 dweights = (int *) malloc(size * sizeof(int));
107 for (use_dup = 0; use_dup <= 1; ++use_dup) {
109 MPI_Dist_graph_neighbors_count(comm, &indegree, &outdegree, &weighted);
112 MPI_Comm_dup(comm, &dupcomm);
113 comm = dupcomm; /* caller retains original comm value */
116 MPI_Topo_test(comm, &topo_type);
117 if (topo_type != MPI_DIST_GRAPH) {
118 fprintf(stderr, "topo_type != MPI_DIST_GRAPH\n");
123 for (i = 0; i < size; i++)
127 fprintf(stderr, "indegree does not match, expected=%d got=%d, layout='%s'\n", indegree,
128 j, graph_layout_name);
133 for (i = 0; i < size; i++)
136 if (j != outdegree) {
137 fprintf(stderr, "outdegree does not match, expected=%d got=%d, layout='%s'\n",
138 outdegree, j, graph_layout_name);
142 if ((indegree || outdegree) && (weighted == 0)) {
143 fprintf(stderr, "MPI_Dist_graph_neighbors_count thinks the graph is not weighted\n");
148 MPI_Dist_graph_neighbors(comm, indegree, sources, sweights, outdegree, destinations,
151 /* For each incoming and outgoing edge in the matrix, search if
152 * the query function listed it in the sources. */
153 for (i = 0; i < size; i++) {
154 if (layout[i][rank]) {
155 for (j = 0; j < indegree; j++) {
156 assert(sources[j] >= 0);
157 assert(sources[j] < size);
162 fprintf(stderr, "no edge from %d to %d specified\n", i, rank);
166 if (sweights[j] != layout[i][rank]) {
167 fprintf(stderr, "incorrect weight for edge (%d,%d): %d instead of %d\n",
168 i, rank, sweights[j], layout[i][rank]);
173 if (layout[rank][i]) {
174 for (j = 0; j < outdegree; j++) {
175 assert(destinations[j] >= 0);
176 assert(destinations[j] < size);
177 if (destinations[j] == i)
180 if (j == outdegree) {
181 fprintf(stderr, "no edge from %d to %d specified\n", rank, i);
185 if (dweights[j] != layout[rank][i]) {
186 fprintf(stderr, "incorrect weight for edge (%d,%d): %d instead of %d\n",
187 rank, i, dweights[j], layout[rank][i]);
194 /* For each incoming and outgoing edge in the sources, we should
195 * have an entry in the matrix */
196 for (i = 0; i < indegree; i++) {
197 if (layout[sources[i]][rank] != sweights[i]) {
198 fprintf(stderr, "edge (%d,%d) has a weight %d instead of %d\n", i, rank,
199 sweights[i], layout[sources[i]][rank]);
203 for (i = 0; i < outdegree; i++) {
204 if (layout[rank][destinations[i]] != dweights[i]) {
205 fprintf(stderr, "edge (%d,%d) has a weight %d instead of %d\n", rank, i,
206 dweights[i], layout[rank][destinations[i]]);
213 if (dupcomm != MPI_COMM_NULL)
214 MPI_Comm_free(&dupcomm);
223 #endif /* At least MPI 2.2 */
225 int main(int argc, char *argv[])
229 int indegree, outdegree, reorder;
230 int check_indegree, check_outdegree, check_weighted;
231 int *sources, *sweights, *destinations, *dweights, *degrees;
234 MTest_Init(&argc, &argv);
236 MPI_Comm_size(MPI_COMM_WORLD, &size);
237 MPI_Comm_rank(MPI_COMM_WORLD, &rank);
239 #if MTEST_HAVE_MIN_MPI_VERSION(2,2)
240 layout = (int **) malloc(size * sizeof(int *));
242 for (i = 0; i < size; i++) {
243 layout[i] = (int *) malloc(size * sizeof(int));
246 /* alloc size*size ints to handle the all-on-one-process case */
247 sources = (int *) malloc(size * size * sizeof(int));
248 sweights = (int *) malloc(size * size * sizeof(int));
249 destinations = (int *) malloc(size * size * sizeof(int));
250 dweights = (int *) malloc(size * size * sizeof(int));
251 degrees = (int *) malloc(size * size * sizeof(int));
253 for (i = 0; i < NUM_GRAPHS; i++) {
254 create_graph_layout(i);
256 MTestPrintfMsg(1, "using graph layout '%s'\n", graph_layout_name);
259 /* MPI_Dist_graph_create_adjacent */
261 MTestPrintfMsg(1, "testing MPI_Dist_graph_create_adjacent\n");
265 for (j = 0; j < size; j++) {
266 if (layout[j][rank]) {
269 sweights[k++] = layout[j][rank];
275 for (j = 0; j < size; j++) {
276 if (layout[rank][j]) {
279 dweights[k++] = layout[rank][j];
283 for (reorder = 0; reorder <= 1; reorder++) {
284 MPI_Dist_graph_create_adjacent(MPI_COMM_WORLD, indegree, sources, sweights,
285 outdegree, destinations, dweights, MPI_INFO_NULL,
288 errs += verify_comm(comm);
289 MPI_Comm_free(&comm);
292 /* a weak check that passing MPI_UNWEIGHTED doesn't cause
293 * create_adjacent to explode */
294 MPI_Dist_graph_create_adjacent(MPI_COMM_WORLD, indegree, sources, MPI_UNWEIGHTED,
295 outdegree, destinations, MPI_UNWEIGHTED, MPI_INFO_NULL,
298 /* intentionally no verify here, weights won't match */
299 MPI_Comm_free(&comm);
302 /* MPI_Dist_graph_create() where each process specifies its
305 MTestPrintfMsg(1, "testing MPI_Dist_graph_create w/ outgoing only\n");
309 for (j = 0; j < size; j++) {
310 if (layout[rank][j]) {
312 dweights[k++] = layout[rank][j];
316 for (reorder = 0; reorder <= 1; reorder++) {
317 MPI_Dist_graph_create(MPI_COMM_WORLD, 1, sources, degrees, destinations, dweights,
318 MPI_INFO_NULL, reorder, &comm);
320 errs += verify_comm(comm);
321 MPI_Comm_free(&comm);
325 /* MPI_Dist_graph_create() where each process specifies its
328 MTestPrintfMsg(1, "testing MPI_Dist_graph_create w/ incoming only\n");
331 for (j = 0; j < size; j++) {
332 if (layout[j][rank]) {
334 sweights[k] = layout[j][rank];
336 destinations[k++] = rank;
339 for (reorder = 0; reorder <= 1; reorder++) {
340 MPI_Dist_graph_create(MPI_COMM_WORLD, k, sources, degrees, destinations, sweights,
341 MPI_INFO_NULL, reorder, &comm);
343 errs += verify_comm(comm);
344 MPI_Comm_free(&comm);
348 /* MPI_Dist_graph_create() where rank 0 specifies the entire
351 MTestPrintfMsg(1, "testing MPI_Dist_graph_create w/ rank 0 specifies only\n");
354 for (j = 0; j < size; j++) {
355 for (k = 0; k < size; k++) {
358 sweights[p] = layout[j][k];
360 destinations[p++] = k;
364 for (reorder = 0; reorder <= 1; reorder++) {
365 MPI_Dist_graph_create(MPI_COMM_WORLD, (rank == 0) ? p : 0, sources, degrees,
366 destinations, sweights, MPI_INFO_NULL, reorder, &comm);
368 errs += verify_comm(comm);
369 MPI_Comm_free(&comm);
372 /* MPI_Dist_graph_create() where rank 0 specifies the entire
373 * graph and all other ranks pass NULL. Can catch implementation
374 * problems when MPI_UNWEIGHTED==NULL. */
376 MTestPrintfMsg(1, "testing MPI_Dist_graph_create w/ rank 0 specifies only -- NULLs\n");
379 for (j = 0; j < size; j++) {
380 for (k = 0; k < size; k++) {
383 sweights[p] = layout[j][k];
385 destinations[p++] = k;
389 for (reorder = 0; reorder <= 1; reorder++) {
391 MPI_Dist_graph_create(MPI_COMM_WORLD, p, sources, degrees,
392 destinations, sweights, MPI_INFO_NULL, reorder, &comm);
395 MPI_Dist_graph_create(MPI_COMM_WORLD, 0, NULL, NULL,
396 NULL, NULL, MPI_INFO_NULL, reorder, &comm);
399 errs += verify_comm(comm);
400 MPI_Comm_free(&comm);
405 /* now tests that don't depend on the layout[][] array */
407 /* The MPI-2.2 standard recommends implementations set
408 * MPI_UNWEIGHTED==NULL, but this leads to an ambiguity. The draft
409 * MPI-3.0 standard specifically recommends _not_ setting it equal
411 if (MPI_UNWEIGHTED == NULL) {
412 fprintf(stderr, "MPI_UNWEIGHTED should not be NULL\n");
416 /* MPI_Dist_graph_create() with no graph */
418 MTestPrintfMsg(1, "testing MPI_Dist_graph_create w/ no graph\n");
420 for (reorder = 0; reorder <= 1; reorder++) {
421 MPI_Dist_graph_create(MPI_COMM_WORLD, 0, sources, degrees,
422 destinations, sweights, MPI_INFO_NULL, reorder, &comm);
423 MPI_Dist_graph_neighbors_count(comm, &check_indegree, &check_outdegree, &check_weighted);
424 if (!check_weighted) {
425 fprintf(stderr, "expected weighted == TRUE for the \"no graph\" case\n");
428 MPI_Comm_free(&comm);
431 /* MPI_Dist_graph_create() with no graph -- passing MPI_WEIGHTS_EMPTY
433 /* NOTE that MPI_WEIGHTS_EMPTY was added in MPI-3 and does not
434 * appear before then. This part of the test thus requires a check
435 * on the MPI major version */
438 MTestPrintfMsg(1, "testing MPI_Dist_graph_create w/ no graph\n");
440 for (reorder = 0; reorder <= 1; reorder++) {
441 MPI_Dist_graph_create(MPI_COMM_WORLD, 0, sources, degrees,
442 destinations, MPI_WEIGHTS_EMPTY, MPI_INFO_NULL, reorder, &comm);
443 MPI_Dist_graph_neighbors_count(comm, &check_indegree, &check_outdegree, &check_weighted);
444 if (!check_weighted) {
446 "expected weighted == TRUE for the \"no graph -- MPI_WEIGHTS_EMPTY\" case\n");
449 MPI_Comm_free(&comm);
453 /* MPI_Dist_graph_create() with no graph -- passing NULLs instead */
455 MTestPrintfMsg(1, "testing MPI_Dist_graph_create w/ no graph -- NULLs\n");
457 for (reorder = 0; reorder <= 1; reorder++) {
458 MPI_Dist_graph_create(MPI_COMM_WORLD, 0, NULL, NULL,
459 NULL, NULL, MPI_INFO_NULL, reorder, &comm);
460 MPI_Dist_graph_neighbors_count(comm, &check_indegree, &check_outdegree, &check_weighted);
461 /* ambiguous if they are equal, only check when they are distinct values. */
462 if (MPI_UNWEIGHTED != NULL) {
463 if (!check_weighted) {
464 fprintf(stderr, "expected weighted == TRUE for the \"no graph -- NULLs\" case\n");
468 MPI_Comm_free(&comm);
471 /* MPI_Dist_graph_create() with no graph -- passing NULLs+MPI_UNWEIGHTED instead */
473 MTestPrintfMsg(1, "testing MPI_Dist_graph_create w/ no graph -- NULLs+MPI_UNWEIGHTED\n");
475 for (reorder = 0; reorder <= 1; reorder++) {
476 MPI_Dist_graph_create(MPI_COMM_WORLD, 0, NULL, NULL,
477 NULL, MPI_UNWEIGHTED, MPI_INFO_NULL, reorder, &comm);
478 MPI_Dist_graph_neighbors_count(comm, &check_indegree, &check_outdegree, &check_weighted);
479 /* ambiguous if they are equal, only check when they are distinct values. */
480 if (MPI_UNWEIGHTED != NULL) {
481 if (check_weighted) {
483 "expected weighted == FALSE for the \"no graph -- NULLs+MPI_UNWEIGHTED\" case\n");
487 MPI_Comm_free(&comm);
490 /* MPI_Dist_graph_create_adjacent() with no graph */
492 MTestPrintfMsg(1, "testing MPI_Dist_graph_create_adjacent w/ no graph\n");
494 for (reorder = 0; reorder <= 1; reorder++) {
495 MPI_Dist_graph_create_adjacent(MPI_COMM_WORLD, 0, sources, sweights,
496 0, destinations, dweights, MPI_INFO_NULL, reorder, &comm);
497 MPI_Dist_graph_neighbors_count(comm, &check_indegree, &check_outdegree, &check_weighted);
498 if (!check_weighted) {
499 fprintf(stderr, "expected weighted == TRUE for the \"no graph\" case\n");
502 MPI_Comm_free(&comm);
505 /* MPI_Dist_graph_create_adjacent() with no graph -- passing MPI_WEIGHTS_EMPTY instead */
506 /* NOTE that MPI_WEIGHTS_EMPTY was added in MPI-3 and does not
507 * appear before then. This part of the test thus requires a check
508 * on the MPI major version */
512 "testing MPI_Dist_graph_create_adjacent w/ no graph -- MPI_WEIGHTS_EMPTY\n");
514 for (reorder = 0; reorder <= 1; reorder++) {
515 MPI_Dist_graph_create_adjacent(MPI_COMM_WORLD, 0, sources, MPI_WEIGHTS_EMPTY,
516 0, destinations, MPI_WEIGHTS_EMPTY, MPI_INFO_NULL, reorder,
518 MPI_Dist_graph_neighbors_count(comm, &check_indegree, &check_outdegree, &check_weighted);
519 if (!check_weighted) {
521 "expected weighted == TRUE for the \"no graph -- MPI_WEIGHTS_EMPTY\" case\n");
524 MPI_Comm_free(&comm);
528 /* MPI_Dist_graph_create_adjacent() with no graph -- passing NULLs instead */
530 MTestPrintfMsg(1, "testing MPI_Dist_graph_create_adjacent w/ no graph -- NULLs\n");
532 for (reorder = 0; reorder <= 1; reorder++) {
533 MPI_Dist_graph_create_adjacent(MPI_COMM_WORLD, 0, NULL, NULL,
534 0, NULL, NULL, MPI_INFO_NULL, reorder, &comm);
535 MPI_Dist_graph_neighbors_count(comm, &check_indegree, &check_outdegree, &check_weighted);
536 /* ambiguous if they are equal, only check when they are distinct values. */
537 if (MPI_UNWEIGHTED != NULL) {
538 if (!check_weighted) {
539 fprintf(stderr, "expected weighted == TRUE for the \"no graph -- NULLs\" case\n");
543 MPI_Comm_free(&comm);
546 /* MPI_Dist_graph_create_adjacent() with no graph -- passing NULLs+MPI_UNWEIGHTED instead */
549 "testing MPI_Dist_graph_create_adjacent w/ no graph -- NULLs+MPI_UNWEIGHTED\n");
551 for (reorder = 0; reorder <= 1; reorder++) {
552 MPI_Dist_graph_create_adjacent(MPI_COMM_WORLD, 0, NULL, MPI_UNWEIGHTED,
553 0, NULL, MPI_UNWEIGHTED, MPI_INFO_NULL, reorder, &comm);
554 MPI_Dist_graph_neighbors_count(comm, &check_indegree, &check_outdegree, &check_weighted);
555 /* ambiguous if they are equal, only check when they are distinct values. */
556 if (MPI_UNWEIGHTED != NULL) {
557 if (check_weighted) {
559 "expected weighted == FALSE for the \"no graph -- NULLs+MPI_UNWEIGHTED\" case\n");
563 MPI_Comm_free(&comm);
567 for (i = 0; i < size; i++)
577 MTest_Finalize(errs);