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, j, graph_layout_name);
132 for (i = 0; i < size; i++)
135 if (j != outdegree) {
136 fprintf(stderr, "outdegree does not match, expected=%d got=%d, layout='%s'\n", outdegree, j, graph_layout_name);
140 if ((indegree || outdegree) && (weighted == 0)) {
141 fprintf(stderr, "MPI_Dist_graph_neighbors_count thinks the graph is not weighted\n");
146 MPI_Dist_graph_neighbors(comm, indegree, sources, sweights, outdegree, destinations, dweights);
148 /* For each incoming and outgoing edge in the matrix, search if
149 * the query function listed it in the sources. */
150 for (i = 0; i < size; i++) {
151 if (layout[i][rank]) {
152 for (j = 0; j < indegree; j++) {
153 assert(sources[j] >= 0);
154 assert(sources[j] < size);
159 fprintf(stderr, "no edge from %d to %d specified\n", i, rank);
163 if (sweights[j] != layout[i][rank]) {
164 fprintf(stderr, "incorrect weight for edge (%d,%d): %d instead of %d\n",
165 i, rank, sweights[j], layout[i][rank]);
170 if (layout[rank][i]) {
171 for (j = 0; j < outdegree; j++) {
172 assert(destinations[j] >= 0);
173 assert(destinations[j] < size);
174 if (destinations[j] == i)
177 if (j == outdegree) {
178 fprintf(stderr, "no edge from %d to %d specified\n", rank, i);
182 if (dweights[j] != layout[rank][i]) {
183 fprintf(stderr, "incorrect weight for edge (%d,%d): %d instead of %d\n",
184 rank, i, dweights[j], layout[rank][i]);
191 /* For each incoming and outgoing edge in the sources, we should
192 * have an entry in the matrix */
193 for (i = 0; i < indegree; i++) {
194 if (layout[sources[i]][rank] != sweights[i]) {
195 fprintf(stderr, "edge (%d,%d) has a weight %d instead of %d\n", i, rank,
196 sweights[i], layout[sources[i]][rank]);
200 for (i = 0; i < outdegree; i++) {
201 if (layout[rank][destinations[i]] != dweights[i]) {
202 fprintf(stderr, "edge (%d,%d) has a weight %d instead of %d\n", rank, i,
203 dweights[i], layout[rank][destinations[i]]);
210 if (dupcomm != MPI_COMM_NULL)
211 MPI_Comm_free(&dupcomm);
216 #endif /* At least MPI 2.2 */
218 int main(int argc, char *argv[])
222 int indegree, outdegree, reorder;
223 int check_indegree, check_outdegree, check_weighted;
224 int *sources, *sweights, *destinations, *dweights, *degrees;
227 MTest_Init(&argc, &argv);
229 MPI_Comm_size(MPI_COMM_WORLD, &size);
230 MPI_Comm_rank(MPI_COMM_WORLD, &rank);
232 #if MTEST_HAVE_MIN_MPI_VERSION(2,2)
233 layout = (int **) malloc(size * sizeof(int *));
235 for (i = 0; i < size; i++) {
236 layout[i] = (int *) malloc(size * sizeof(int));
239 /* alloc size*size ints to handle the all-on-one-process case */
240 sources = (int *) malloc(size * size * sizeof(int));
241 sweights = (int *) malloc(size * size * sizeof(int));
242 destinations = (int *) malloc(size * size * sizeof(int));
243 dweights = (int *) malloc(size * size * sizeof(int));
244 degrees = (int *) malloc(size * size * sizeof(int));
246 for (i = 0; i < NUM_GRAPHS; i++) {
247 create_graph_layout(i);
249 MTestPrintfMsg( 1, "using graph layout '%s'\n", graph_layout_name );
252 /* MPI_Dist_graph_create_adjacent */
254 MTestPrintfMsg( 1, "testing MPI_Dist_graph_create_adjacent\n" );
258 for (j = 0; j < size; j++) {
259 if (layout[j][rank]) {
262 sweights[k++] = layout[j][rank];
268 for (j = 0; j < size; j++) {
269 if (layout[rank][j]) {
272 dweights[k++] = layout[rank][j];
276 for (reorder = 0; reorder <= 1; reorder++) {
277 MPI_Dist_graph_create_adjacent(MPI_COMM_WORLD, indegree, sources, sweights,
278 outdegree, destinations, dweights, MPI_INFO_NULL,
281 errs += verify_comm(comm);
282 MPI_Comm_free(&comm);
285 /* a weak check that passing MPI_UNWEIGHTED doesn't cause
286 * create_adjacent to explode */
287 MPI_Dist_graph_create_adjacent(MPI_COMM_WORLD, indegree, sources, MPI_UNWEIGHTED,
288 outdegree, destinations, MPI_UNWEIGHTED, MPI_INFO_NULL,
291 /* intentionally no verify here, weights won't match */
292 MPI_Comm_free(&comm);
295 /* MPI_Dist_graph_create() where each process specifies its
299 "testing MPI_Dist_graph_create w/ outgoing only\n" );
303 for (j = 0; j < size; j++) {
304 if (layout[rank][j]) {
306 dweights[k++] = layout[rank][j];
310 for (reorder = 0; reorder <= 1; reorder++) {
311 MPI_Dist_graph_create(MPI_COMM_WORLD, 1, sources, degrees, destinations, dweights,
312 MPI_INFO_NULL, reorder, &comm);
314 errs += verify_comm(comm);
315 MPI_Comm_free(&comm);
319 /* MPI_Dist_graph_create() where each process specifies its
323 "testing MPI_Dist_graph_create w/ incoming only\n" );
326 for (j = 0; j < size; j++) {
327 if (layout[j][rank]) {
329 sweights[k] = layout[j][rank];
331 destinations[k++] = rank;
334 for (reorder = 0; reorder <= 1; reorder++) {
335 MPI_Dist_graph_create(MPI_COMM_WORLD, k, sources, degrees, destinations, sweights,
336 MPI_INFO_NULL, reorder, &comm);
338 errs += verify_comm(comm);
339 MPI_Comm_free(&comm);
343 /* MPI_Dist_graph_create() where rank 0 specifies the entire
347 "testing MPI_Dist_graph_create w/ rank 0 specifies only\n" );
350 for (j = 0; j < size; j++) {
351 for (k = 0; k < size; k++) {
354 sweights[p] = layout[j][k];
356 destinations[p++] = k;
360 for (reorder = 0; reorder <= 1; reorder++) {
361 MPI_Dist_graph_create(MPI_COMM_WORLD, (rank == 0) ? p : 0, sources, degrees,
362 destinations, sweights, MPI_INFO_NULL, reorder, &comm);
364 errs += verify_comm(comm);
365 MPI_Comm_free(&comm);
368 /* MPI_Dist_graph_create() where rank 0 specifies the entire
369 * graph and all other ranks pass NULL. Can catch implementation
370 * problems when MPI_UNWEIGHTED==NULL. */
373 "testing MPI_Dist_graph_create w/ rank 0 specifies only -- NULLs\n");
376 for (j = 0; j < size; j++) {
377 for (k = 0; k < size; k++) {
380 sweights[p] = layout[j][k];
382 destinations[p++] = k;
386 for (reorder = 0; reorder <= 1; reorder++) {
388 MPI_Dist_graph_create(MPI_COMM_WORLD, p, sources, degrees,
389 destinations, sweights, MPI_INFO_NULL, reorder, &comm);
392 MPI_Dist_graph_create(MPI_COMM_WORLD, 0, NULL, NULL,
393 NULL, NULL, MPI_INFO_NULL, reorder, &comm);
396 errs += verify_comm(comm);
397 MPI_Comm_free(&comm);
402 /* now tests that don't depend on the layout[][] array */
404 /* The MPI-2.2 standard recommends implementations set
405 * MPI_UNWEIGHTED==NULL, but this leads to an ambiguity. The draft
406 * MPI-3.0 standard specifically recommends _not_ setting it equal
408 if (MPI_UNWEIGHTED == NULL) {
409 fprintf(stderr, "MPI_UNWEIGHTED should not be NULL\n");
413 /* MPI_Dist_graph_create() with no graph */
415 MTestPrintfMsg( 1, "testing MPI_Dist_graph_create w/ no graph\n" );
417 for (reorder = 0; reorder <= 1; reorder++) {
418 MPI_Dist_graph_create(MPI_COMM_WORLD, 0, sources, degrees,
419 destinations, sweights, MPI_INFO_NULL, reorder, &comm);
420 MPI_Dist_graph_neighbors_count(comm, &check_indegree, &check_outdegree, &check_weighted);
421 if (!check_weighted) {
422 fprintf(stderr, "expected weighted == TRUE for the \"no graph\" case\n");
425 MPI_Comm_free(&comm);
428 /* MPI_Dist_graph_create() with no graph -- passing MPI_WEIGHTS_EMPTY
430 /* NOTE that MPI_WEIGHTS_EMPTY was added in MPI-3 and does not
431 appear before then. This part of the test thus requires a check
432 on the MPI major version */
435 MTestPrintfMsg( 1, "testing MPI_Dist_graph_create w/ no graph\n" );
437 for (reorder = 0; reorder <= 1; reorder++) {
438 MPI_Dist_graph_create(MPI_COMM_WORLD, 0, sources, degrees,
439 destinations, MPI_WEIGHTS_EMPTY, MPI_INFO_NULL, reorder, &comm);
440 MPI_Dist_graph_neighbors_count(comm, &check_indegree, &check_outdegree, &check_weighted);
441 if (!check_weighted) {
442 fprintf(stderr, "expected weighted == TRUE for the \"no graph -- MPI_WEIGHTS_EMPTY\" case\n");
445 MPI_Comm_free(&comm);
449 /* MPI_Dist_graph_create() with no graph -- passing NULLs instead */
452 "testing MPI_Dist_graph_create w/ no graph -- NULLs\n" );
454 for (reorder = 0; reorder <= 1; reorder++) {
455 MPI_Dist_graph_create(MPI_COMM_WORLD, 0, NULL, NULL,
456 NULL, NULL, MPI_INFO_NULL, reorder, &comm);
457 MPI_Dist_graph_neighbors_count(comm, &check_indegree, &check_outdegree, &check_weighted);
458 /* ambiguous if they are equal, only check when they are distinct values. */
459 if (MPI_UNWEIGHTED != NULL) {
460 if (!check_weighted) {
461 fprintf(stderr, "expected weighted == TRUE for the \"no graph -- NULLs\" case\n");
465 MPI_Comm_free(&comm);
468 /* MPI_Dist_graph_create() with no graph -- passing NULLs+MPI_UNWEIGHTED instead */
471 "testing MPI_Dist_graph_create w/ no graph -- NULLs+MPI_UNWEIGHTED\n" );
473 for (reorder = 0; reorder <= 1; reorder++) {
474 MPI_Dist_graph_create(MPI_COMM_WORLD, 0, NULL, NULL,
475 NULL, MPI_UNWEIGHTED, MPI_INFO_NULL, reorder, &comm);
476 MPI_Dist_graph_neighbors_count(comm, &check_indegree, &check_outdegree, &check_weighted);
477 /* ambiguous if they are equal, only check when they are distinct values. */
478 if (MPI_UNWEIGHTED != NULL) {
479 if (check_weighted) {
480 fprintf(stderr, "expected weighted == FALSE for the \"no graph -- NULLs+MPI_UNWEIGHTED\" case\n");
484 MPI_Comm_free(&comm);
487 /* MPI_Dist_graph_create_adjacent() with no graph */
490 "testing MPI_Dist_graph_create_adjacent w/ no graph\n" );
492 for (reorder = 0; reorder <= 1; reorder++) {
493 MPI_Dist_graph_create_adjacent(MPI_COMM_WORLD, 0, sources, sweights,
494 0, destinations, dweights, MPI_INFO_NULL, reorder, &comm);
495 MPI_Dist_graph_neighbors_count(comm, &check_indegree, &check_outdegree, &check_weighted);
496 if (!check_weighted) {
497 fprintf(stderr, "expected weighted == TRUE for the \"no graph\" case\n");
500 MPI_Comm_free(&comm);
503 /* MPI_Dist_graph_create_adjacent() with no graph -- passing MPI_WEIGHTS_EMPTY instead */
504 /* NOTE that MPI_WEIGHTS_EMPTY was added in MPI-3 and does not
505 appear before then. This part of the test thus requires a check
506 on the MPI major version */
510 "testing MPI_Dist_graph_create_adjacent w/ no graph -- MPI_WEIGHTS_EMPTY\n" );
512 for (reorder = 0; reorder <= 1; reorder++) {
513 MPI_Dist_graph_create_adjacent(MPI_COMM_WORLD, 0, sources, MPI_WEIGHTS_EMPTY,
514 0, destinations, MPI_WEIGHTS_EMPTY, MPI_INFO_NULL, reorder, &comm);
515 MPI_Dist_graph_neighbors_count(comm, &check_indegree, &check_outdegree, &check_weighted);
516 if (!check_weighted) {
517 fprintf(stderr, "expected weighted == TRUE for the \"no graph -- MPI_WEIGHTS_EMPTY\" case\n");
520 MPI_Comm_free(&comm);
524 /* MPI_Dist_graph_create_adjacent() with no graph -- passing NULLs instead */
527 "testing MPI_Dist_graph_create_adjacent w/ no graph -- NULLs\n" );
529 for (reorder = 0; reorder <= 1; reorder++) {
530 MPI_Dist_graph_create_adjacent(MPI_COMM_WORLD, 0, NULL, NULL,
531 0, NULL, NULL, MPI_INFO_NULL, reorder, &comm);
532 MPI_Dist_graph_neighbors_count(comm, &check_indegree, &check_outdegree, &check_weighted);
533 /* ambiguous if they are equal, only check when they are distinct values. */
534 if (MPI_UNWEIGHTED != NULL) {
535 if (!check_weighted) {
536 fprintf(stderr, "expected weighted == TRUE for the \"no graph -- NULLs\" case\n");
540 MPI_Comm_free(&comm);
543 /* MPI_Dist_graph_create_adjacent() with no graph -- passing NULLs+MPI_UNWEIGHTED instead */
546 "testing MPI_Dist_graph_create_adjacent w/ no graph -- NULLs+MPI_UNWEIGHTED\n");
548 for (reorder = 0; reorder <= 1; reorder++) {
549 MPI_Dist_graph_create_adjacent(MPI_COMM_WORLD, 0, NULL, MPI_UNWEIGHTED,
550 0, NULL, MPI_UNWEIGHTED, MPI_INFO_NULL, reorder, &comm);
551 MPI_Dist_graph_neighbors_count(comm, &check_indegree, &check_outdegree, &check_weighted);
552 /* ambiguous if they are equal, only check when they are distinct values. */
553 if (MPI_UNWEIGHTED != NULL) {
554 if (check_weighted) {
555 fprintf(stderr, "expected weighted == FALSE for the \"no graph -- NULLs+MPI_UNWEIGHTED\" case\n");
559 MPI_Comm_free(&comm);
563 for (i = 0; i < size; i++)
568 MTest_Finalize(errs);