Add build scripts for Jenkins test jobs
[mpich.git] / test / mpi / topo / distgraph1.c
1 /* -*- Mode: C; c-basic-offset:4 ; indent-tabs-mode:nil ; -*- */
2 /*
3  *
4  *  (C) 2003 by Argonne National Laboratory.
5  *      See COPYRIGHT in top-level directory.
6  */
7
8 #include "mpi.h"
9 #include <stdio.h>
10 #include <stdlib.h>
11 #include <string.h>
12 #include <assert.h>
13 #include "mpitest.h"
14
15 #define NUM_GRAPHS 10
16 #define MAX_WEIGHT 100
17
18 /* convenience globals */
19 int size, rank;
20
21 /* We need MPI 2.2 to be able to compile the following routines. */
22 #if MTEST_HAVE_MIN_MPI_VERSION(2,2)
23
24 /* Maybe use a bit vector instead? */
25 int **layout;
26
27 #define MAX_LAYOUT_NAME_LEN 256
28 char graph_layout_name[MAX_LAYOUT_NAME_LEN] = { '\0' };
29
30 static void create_graph_layout(int graph_num)
31 {
32     int i, j;
33
34     if (rank == 0) {
35         switch (graph_num) {
36         case 0:
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);
41             break;
42         case 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);
47             break;
48         case 2:
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);
54                     else
55                         layout[i][j] = 0;
56                 }
57             }
58             break;
59         case 3:
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++)
63                     layout[i][j] = 0;
64             break;
65         default:
66             strncpy(graph_layout_name, "a random incomplete graph", MAX_LAYOUT_NAME_LEN);
67             srand(graph_num);
68
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)
75                         layout[i][j] = 0;
76                     else
77                         layout[i][j] = rand() % MAX_WEIGHT;
78                 }
79             }
80             break;
81         }
82     }
83
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);
89     }
90 }
91
92 static int verify_comm(MPI_Comm comm)
93 {
94     int local_errs = 0;
95     int i, j;
96     int indegree, outdegree, weighted;
97     int *sources, *sweights, *destinations, *dweights;
98     int use_dup;
99     int topo_type = MPI_UNDEFINED;
100     MPI_Comm dupcomm = MPI_COMM_NULL;
101
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));
106
107     for (use_dup = 0; use_dup <= 1; ++use_dup) {
108         if (!use_dup) {
109             MPI_Dist_graph_neighbors_count(comm, &indegree, &outdegree, &weighted);
110         }
111         else {
112             MPI_Comm_dup(comm, &dupcomm);
113             comm = dupcomm;     /* caller retains original comm value */
114         }
115
116         MPI_Topo_test(comm, &topo_type);
117         if (topo_type != MPI_DIST_GRAPH) {
118             fprintf(stderr, "topo_type != MPI_DIST_GRAPH\n");
119             ++local_errs;
120         }
121
122         j = 0;
123         for (i = 0; i < size; i++)
124             if (layout[i][rank])
125                 j++;
126         if (j != indegree) {
127             fprintf(stderr, "indegree does not match, expected=%d got=%d, layout='%s'\n", indegree,
128                     j, graph_layout_name);
129             ++local_errs;
130         }
131
132         j = 0;
133         for (i = 0; i < size; i++)
134             if (layout[rank][i])
135                 j++;
136         if (j != outdegree) {
137             fprintf(stderr, "outdegree does not match, expected=%d got=%d, layout='%s'\n",
138                     outdegree, j, graph_layout_name);
139             ++local_errs;
140         }
141
142         if ((indegree || outdegree) && (weighted == 0)) {
143             fprintf(stderr, "MPI_Dist_graph_neighbors_count thinks the graph is not weighted\n");
144             ++local_errs;
145         }
146
147
148         MPI_Dist_graph_neighbors(comm, indegree, sources, sweights, outdegree, destinations,
149                                  dweights);
150
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);
158                     if (sources[j] == i)
159                         break;
160                 }
161                 if (j == indegree) {
162                     fprintf(stderr, "no edge from %d to %d specified\n", i, rank);
163                     ++local_errs;
164                 }
165                 else {
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]);
169                         ++local_errs;
170                     }
171                 }
172             }
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)
178                         break;
179                 }
180                 if (j == outdegree) {
181                     fprintf(stderr, "no edge from %d to %d specified\n", rank, i);
182                     ++local_errs;
183                 }
184                 else {
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]);
188                         ++local_errs;
189                     }
190                 }
191             }
192         }
193
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]);
200                 ++local_errs;
201             }
202         }
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]]);
207                 ++local_errs;
208             }
209         }
210
211     }
212
213     if (dupcomm != MPI_COMM_NULL)
214         MPI_Comm_free(&dupcomm);
215
216     return local_errs;
217 }
218
219 #endif /* At least MPI 2.2 */
220
221 int main(int argc, char *argv[])
222 {
223     int errs = 0;
224     int i, j, k, p;
225     int indegree, outdegree, reorder;
226     int check_indegree, check_outdegree, check_weighted;
227     int *sources, *sweights, *destinations, *dweights, *degrees;
228     MPI_Comm comm;
229
230     MTest_Init(&argc, &argv);
231
232     MPI_Comm_size(MPI_COMM_WORLD, &size);
233     MPI_Comm_rank(MPI_COMM_WORLD, &rank);
234
235 #if MTEST_HAVE_MIN_MPI_VERSION(2,2)
236     layout = (int **) malloc(size * sizeof(int *));
237     assert(layout);
238     for (i = 0; i < size; i++) {
239         layout[i] = (int *) malloc(size * sizeof(int));
240         assert(layout[i]);
241     }
242     /* alloc size*size ints to handle the all-on-one-process case */
243     sources = (int *) malloc(size * size * sizeof(int));
244     sweights = (int *) malloc(size * size * sizeof(int));
245     destinations = (int *) malloc(size * size * sizeof(int));
246     dweights = (int *) malloc(size * size * sizeof(int));
247     degrees = (int *) malloc(size * size * sizeof(int));
248
249     for (i = 0; i < NUM_GRAPHS; i++) {
250         create_graph_layout(i);
251         if (rank == 0) {
252             MTestPrintfMsg(1, "using graph layout '%s'\n", graph_layout_name);
253         }
254
255         /* MPI_Dist_graph_create_adjacent */
256         if (rank == 0) {
257             MTestPrintfMsg(1, "testing MPI_Dist_graph_create_adjacent\n");
258         }
259         indegree = 0;
260         k = 0;
261         for (j = 0; j < size; j++) {
262             if (layout[j][rank]) {
263                 indegree++;
264                 sources[k] = j;
265                 sweights[k++] = layout[j][rank];
266             }
267         }
268
269         outdegree = 0;
270         k = 0;
271         for (j = 0; j < size; j++) {
272             if (layout[rank][j]) {
273                 outdegree++;
274                 destinations[k] = j;
275                 dweights[k++] = layout[rank][j];
276             }
277         }
278
279         for (reorder = 0; reorder <= 1; reorder++) {
280             MPI_Dist_graph_create_adjacent(MPI_COMM_WORLD, indegree, sources, sweights,
281                                            outdegree, destinations, dweights, MPI_INFO_NULL,
282                                            reorder, &comm);
283             MPI_Barrier(comm);
284             errs += verify_comm(comm);
285             MPI_Comm_free(&comm);
286         }
287
288         /* a weak check that passing MPI_UNWEIGHTED doesn't cause
289          * create_adjacent to explode */
290         MPI_Dist_graph_create_adjacent(MPI_COMM_WORLD, indegree, sources, MPI_UNWEIGHTED,
291                                        outdegree, destinations, MPI_UNWEIGHTED, MPI_INFO_NULL,
292                                        reorder, &comm);
293         MPI_Barrier(comm);
294         /* intentionally no verify here, weights won't match */
295         MPI_Comm_free(&comm);
296
297
298         /* MPI_Dist_graph_create() where each process specifies its
299          * outgoing edges */
300         if (rank == 0) {
301             MTestPrintfMsg(1, "testing MPI_Dist_graph_create w/ outgoing only\n");
302         }
303         sources[0] = rank;
304         k = 0;
305         for (j = 0; j < size; j++) {
306             if (layout[rank][j]) {
307                 destinations[k] = j;
308                 dweights[k++] = layout[rank][j];
309             }
310         }
311         degrees[0] = k;
312         for (reorder = 0; reorder <= 1; reorder++) {
313             MPI_Dist_graph_create(MPI_COMM_WORLD, 1, sources, degrees, destinations, dweights,
314                                   MPI_INFO_NULL, reorder, &comm);
315             MPI_Barrier(comm);
316             errs += verify_comm(comm);
317             MPI_Comm_free(&comm);
318         }
319
320
321         /* MPI_Dist_graph_create() where each process specifies its
322          * incoming edges */
323         if (rank == 0) {
324             MTestPrintfMsg(1, "testing MPI_Dist_graph_create w/ incoming only\n");
325         }
326         k = 0;
327         for (j = 0; j < size; j++) {
328             if (layout[j][rank]) {
329                 sources[k] = j;
330                 sweights[k] = layout[j][rank];
331                 degrees[k] = 1;
332                 destinations[k++] = rank;
333             }
334         }
335         for (reorder = 0; reorder <= 1; reorder++) {
336             MPI_Dist_graph_create(MPI_COMM_WORLD, k, sources, degrees, destinations, sweights,
337                                   MPI_INFO_NULL, reorder, &comm);
338             MPI_Barrier(comm);
339             errs += verify_comm(comm);
340             MPI_Comm_free(&comm);
341         }
342
343
344         /* MPI_Dist_graph_create() where rank 0 specifies the entire
345          * graph */
346         if (rank == 0) {
347             MTestPrintfMsg(1, "testing MPI_Dist_graph_create w/ rank 0 specifies only\n");
348         }
349         p = 0;
350         for (j = 0; j < size; j++) {
351             for (k = 0; k < size; k++) {
352                 if (layout[j][k]) {
353                     sources[p] = j;
354                     sweights[p] = layout[j][k];
355                     degrees[p] = 1;
356                     destinations[p++] = k;
357                 }
358             }
359         }
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);
363             MPI_Barrier(comm);
364             errs += verify_comm(comm);
365             MPI_Comm_free(&comm);
366         }
367
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. */
371         if (rank == 0) {
372             MTestPrintfMsg(1, "testing MPI_Dist_graph_create w/ rank 0 specifies only -- NULLs\n");
373         }
374         p = 0;
375         for (j = 0; j < size; j++) {
376             for (k = 0; k < size; k++) {
377                 if (layout[j][k]) {
378                     sources[p] = j;
379                     sweights[p] = layout[j][k];
380                     degrees[p] = 1;
381                     destinations[p++] = k;
382                 }
383             }
384         }
385         for (reorder = 0; reorder <= 1; reorder++) {
386             if (rank == 0) {
387                 MPI_Dist_graph_create(MPI_COMM_WORLD, p, sources, degrees,
388                                       destinations, sweights, MPI_INFO_NULL, reorder, &comm);
389             }
390             else {
391                 MPI_Dist_graph_create(MPI_COMM_WORLD, 0, NULL, NULL,
392                                       NULL, NULL, MPI_INFO_NULL, reorder, &comm);
393             }
394             MPI_Barrier(comm);
395             errs += verify_comm(comm);
396             MPI_Comm_free(&comm);
397         }
398
399     }
400
401     /* now tests that don't depend on the layout[][] array */
402
403     /* The MPI-2.2 standard recommends implementations set
404      * MPI_UNWEIGHTED==NULL, but this leads to an ambiguity.  The draft
405      * MPI-3.0 standard specifically recommends _not_ setting it equal
406      * to NULL. */
407     if (MPI_UNWEIGHTED == NULL) {
408         fprintf(stderr, "MPI_UNWEIGHTED should not be NULL\n");
409         ++errs;
410     }
411
412     /* MPI_Dist_graph_create() with no graph */
413     if (rank == 0) {
414         MTestPrintfMsg(1, "testing MPI_Dist_graph_create w/ no graph\n");
415     }
416     for (reorder = 0; reorder <= 1; reorder++) {
417         MPI_Dist_graph_create(MPI_COMM_WORLD, 0, sources, degrees,
418                               destinations, sweights, MPI_INFO_NULL, reorder, &comm);
419         MPI_Dist_graph_neighbors_count(comm, &check_indegree, &check_outdegree, &check_weighted);
420         if (!check_weighted) {
421             fprintf(stderr, "expected weighted == TRUE for the \"no graph\" case\n");
422             ++errs;
423         }
424         MPI_Comm_free(&comm);
425     }
426
427     /* MPI_Dist_graph_create() with no graph -- passing MPI_WEIGHTS_EMPTY
428      * instead */
429     /* NOTE that MPI_WEIGHTS_EMPTY was added in MPI-3 and does not
430      * appear before then.  This part of the test thus requires a check
431      * on the MPI major version */
432 #if MPI_VERSION >= 3
433     if (rank == 0) {
434         MTestPrintfMsg(1, "testing MPI_Dist_graph_create w/ no graph\n");
435     }
436     for (reorder = 0; reorder <= 1; reorder++) {
437         MPI_Dist_graph_create(MPI_COMM_WORLD, 0, sources, degrees,
438                               destinations, MPI_WEIGHTS_EMPTY, MPI_INFO_NULL, reorder, &comm);
439         MPI_Dist_graph_neighbors_count(comm, &check_indegree, &check_outdegree, &check_weighted);
440         if (!check_weighted) {
441             fprintf(stderr,
442                     "expected weighted == TRUE for the \"no graph -- MPI_WEIGHTS_EMPTY\" case\n");
443             ++errs;
444         }
445         MPI_Comm_free(&comm);
446     }
447 #endif
448
449     /* MPI_Dist_graph_create() with no graph -- passing NULLs instead */
450     if (rank == 0) {
451         MTestPrintfMsg(1, "testing MPI_Dist_graph_create w/ no graph -- NULLs\n");
452     }
453     for (reorder = 0; reorder <= 1; reorder++) {
454         MPI_Dist_graph_create(MPI_COMM_WORLD, 0, NULL, NULL,
455                               NULL, NULL, MPI_INFO_NULL, reorder, &comm);
456         MPI_Dist_graph_neighbors_count(comm, &check_indegree, &check_outdegree, &check_weighted);
457         /* ambiguous if they are equal, only check when they are distinct values. */
458         if (MPI_UNWEIGHTED != NULL) {
459             if (!check_weighted) {
460                 fprintf(stderr, "expected weighted == TRUE for the \"no graph -- NULLs\" case\n");
461                 ++errs;
462             }
463         }
464         MPI_Comm_free(&comm);
465     }
466
467     /* MPI_Dist_graph_create() with no graph -- passing NULLs+MPI_UNWEIGHTED instead */
468     if (rank == 0) {
469         MTestPrintfMsg(1, "testing MPI_Dist_graph_create w/ no graph -- NULLs+MPI_UNWEIGHTED\n");
470     }
471     for (reorder = 0; reorder <= 1; reorder++) {
472         MPI_Dist_graph_create(MPI_COMM_WORLD, 0, NULL, NULL,
473                               NULL, MPI_UNWEIGHTED, MPI_INFO_NULL, reorder, &comm);
474         MPI_Dist_graph_neighbors_count(comm, &check_indegree, &check_outdegree, &check_weighted);
475         /* ambiguous if they are equal, only check when they are distinct values. */
476         if (MPI_UNWEIGHTED != NULL) {
477             if (check_weighted) {
478                 fprintf(stderr,
479                         "expected weighted == FALSE for the \"no graph -- NULLs+MPI_UNWEIGHTED\" case\n");
480                 ++errs;
481             }
482         }
483         MPI_Comm_free(&comm);
484     }
485
486     /* MPI_Dist_graph_create_adjacent() with no graph */
487     if (rank == 0) {
488         MTestPrintfMsg(1, "testing MPI_Dist_graph_create_adjacent w/ no graph\n");
489     }
490     for (reorder = 0; reorder <= 1; reorder++) {
491         MPI_Dist_graph_create_adjacent(MPI_COMM_WORLD, 0, sources, sweights,
492                                        0, destinations, dweights, MPI_INFO_NULL, reorder, &comm);
493         MPI_Dist_graph_neighbors_count(comm, &check_indegree, &check_outdegree, &check_weighted);
494         if (!check_weighted) {
495             fprintf(stderr, "expected weighted == TRUE for the \"no graph\" case\n");
496             ++errs;
497         }
498         MPI_Comm_free(&comm);
499     }
500
501     /* MPI_Dist_graph_create_adjacent() with no graph -- passing MPI_WEIGHTS_EMPTY instead */
502     /* NOTE that MPI_WEIGHTS_EMPTY was added in MPI-3 and does not
503      * appear before then.  This part of the test thus requires a check
504      * on the MPI major version */
505 #if MPI_VERSION >= 3
506     if (rank == 0) {
507         MTestPrintfMsg(1,
508                        "testing MPI_Dist_graph_create_adjacent w/ no graph -- MPI_WEIGHTS_EMPTY\n");
509     }
510     for (reorder = 0; reorder <= 1; reorder++) {
511         MPI_Dist_graph_create_adjacent(MPI_COMM_WORLD, 0, sources, MPI_WEIGHTS_EMPTY,
512                                        0, destinations, MPI_WEIGHTS_EMPTY, MPI_INFO_NULL, reorder,
513                                        &comm);
514         MPI_Dist_graph_neighbors_count(comm, &check_indegree, &check_outdegree, &check_weighted);
515         if (!check_weighted) {
516             fprintf(stderr,
517                     "expected weighted == TRUE for the \"no graph -- MPI_WEIGHTS_EMPTY\" case\n");
518             ++errs;
519         }
520         MPI_Comm_free(&comm);
521     }
522 #endif
523
524     /* MPI_Dist_graph_create_adjacent() with no graph -- passing NULLs instead */
525     if (rank == 0) {
526         MTestPrintfMsg(1, "testing MPI_Dist_graph_create_adjacent w/ no graph -- NULLs\n");
527     }
528     for (reorder = 0; reorder <= 1; reorder++) {
529         MPI_Dist_graph_create_adjacent(MPI_COMM_WORLD, 0, NULL, NULL,
530                                        0, NULL, NULL, MPI_INFO_NULL, reorder, &comm);
531         MPI_Dist_graph_neighbors_count(comm, &check_indegree, &check_outdegree, &check_weighted);
532         /* ambiguous if they are equal, only check when they are distinct values. */
533         if (MPI_UNWEIGHTED != NULL) {
534             if (!check_weighted) {
535                 fprintf(stderr, "expected weighted == TRUE for the \"no graph -- NULLs\" case\n");
536                 ++errs;
537             }
538         }
539         MPI_Comm_free(&comm);
540     }
541
542     /* MPI_Dist_graph_create_adjacent() with no graph -- passing NULLs+MPI_UNWEIGHTED instead */
543     if (rank == 0) {
544         MTestPrintfMsg(1,
545                        "testing MPI_Dist_graph_create_adjacent w/ no graph -- NULLs+MPI_UNWEIGHTED\n");
546     }
547     for (reorder = 0; reorder <= 1; reorder++) {
548         MPI_Dist_graph_create_adjacent(MPI_COMM_WORLD, 0, NULL, MPI_UNWEIGHTED,
549                                        0, NULL, MPI_UNWEIGHTED, MPI_INFO_NULL, reorder, &comm);
550         MPI_Dist_graph_neighbors_count(comm, &check_indegree, &check_outdegree, &check_weighted);
551         /* ambiguous if they are equal, only check when they are distinct values. */
552         if (MPI_UNWEIGHTED != NULL) {
553             if (check_weighted) {
554                 fprintf(stderr,
555                         "expected weighted == FALSE for the \"no graph -- NULLs+MPI_UNWEIGHTED\" case\n");
556                 ++errs;
557             }
558         }
559         MPI_Comm_free(&comm);
560     }
561
562
563     for (i = 0; i < size; i++)
564         free(layout[i]);
565     free(layout);
566 #endif
567
568     MTest_Finalize(errs);
569     MPI_Finalize();
570
571     return 0;
572 }