test/dtpools: Fix autoconf to allow separated autotools
[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         } else {
111             MPI_Comm_dup(comm, &dupcomm);
112             comm = dupcomm;     /* caller retains original comm value */
113         }
114
115         MPI_Topo_test(comm, &topo_type);
116         if (topo_type != MPI_DIST_GRAPH) {
117             fprintf(stderr, "topo_type != MPI_DIST_GRAPH\n");
118             ++local_errs;
119         }
120
121         j = 0;
122         for (i = 0; i < size; i++)
123             if (layout[i][rank])
124                 j++;
125         if (j != indegree) {
126             fprintf(stderr, "indegree does not match, expected=%d got=%d, layout='%s'\n", indegree,
127                     j, graph_layout_name);
128             ++local_errs;
129         }
130
131         j = 0;
132         for (i = 0; i < size; i++)
133             if (layout[rank][i])
134                 j++;
135         if (j != outdegree) {
136             fprintf(stderr, "outdegree does not match, expected=%d got=%d, layout='%s'\n",
137                     outdegree, j, graph_layout_name);
138             ++local_errs;
139         }
140
141         if ((indegree || outdegree) && (weighted == 0)) {
142             fprintf(stderr, "MPI_Dist_graph_neighbors_count thinks the graph is not weighted\n");
143             ++local_errs;
144         }
145
146
147         MPI_Dist_graph_neighbors(comm, indegree, sources, sweights, outdegree, destinations,
148                                  dweights);
149
150         /* For each incoming and outgoing edge in the matrix, search if
151          * the query function listed it in the sources. */
152         for (i = 0; i < size; i++) {
153             if (layout[i][rank]) {
154                 for (j = 0; j < indegree; j++) {
155                     assert(sources[j] >= 0);
156                     assert(sources[j] < size);
157                     if (sources[j] == i)
158                         break;
159                 }
160                 if (j == indegree) {
161                     fprintf(stderr, "no edge from %d to %d specified\n", i, rank);
162                     ++local_errs;
163                 } else {
164                     if (sweights[j] != layout[i][rank]) {
165                         fprintf(stderr, "incorrect weight for edge (%d,%d): %d instead of %d\n",
166                                 i, rank, sweights[j], layout[i][rank]);
167                         ++local_errs;
168                     }
169                 }
170             }
171             if (layout[rank][i]) {
172                 for (j = 0; j < outdegree; j++) {
173                     assert(destinations[j] >= 0);
174                     assert(destinations[j] < size);
175                     if (destinations[j] == i)
176                         break;
177                 }
178                 if (j == outdegree) {
179                     fprintf(stderr, "no edge from %d to %d specified\n", rank, i);
180                     ++local_errs;
181                 } else {
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]);
185                         ++local_errs;
186                     }
187                 }
188             }
189         }
190
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]);
197                 ++local_errs;
198             }
199         }
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]]);
204                 ++local_errs;
205             }
206         }
207     }
208
209     if (dupcomm != MPI_COMM_NULL)
210         MPI_Comm_free(&dupcomm);
211
212     free(sources);
213     free(sweights);
214     free(destinations);
215     free(dweights);
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             } else {
390                 MPI_Dist_graph_create(MPI_COMM_WORLD, 0, NULL, NULL,
391                                       NULL, NULL, MPI_INFO_NULL, reorder, &comm);
392             }
393             MPI_Barrier(comm);
394             errs += verify_comm(comm);
395             MPI_Comm_free(&comm);
396         }
397     }
398
399     /* now tests that don't depend on the layout[][] array */
400
401     /* The MPI-2.2 standard recommends implementations set
402      * MPI_UNWEIGHTED==NULL, but this leads to an ambiguity.  The draft
403      * MPI-3.0 standard specifically recommends _not_ setting it equal
404      * to NULL. */
405     if (MPI_UNWEIGHTED == NULL) {
406         fprintf(stderr, "MPI_UNWEIGHTED should not be NULL\n");
407         ++errs;
408     }
409
410     /* MPI_Dist_graph_create() with no graph */
411     if (rank == 0) {
412         MTestPrintfMsg(1, "testing MPI_Dist_graph_create w/ no graph\n");
413     }
414     for (reorder = 0; reorder <= 1; reorder++) {
415         MPI_Dist_graph_create(MPI_COMM_WORLD, 0, sources, degrees,
416                               destinations, sweights, MPI_INFO_NULL, reorder, &comm);
417         MPI_Dist_graph_neighbors_count(comm, &check_indegree, &check_outdegree, &check_weighted);
418         if (!check_weighted) {
419             fprintf(stderr, "expected weighted == TRUE for the \"no graph\" case\n");
420             ++errs;
421         }
422         MPI_Comm_free(&comm);
423     }
424
425     /* MPI_Dist_graph_create() with no graph -- passing MPI_WEIGHTS_EMPTY
426      * instead */
427     /* NOTE that MPI_WEIGHTS_EMPTY was added in MPI-3 and does not
428      * appear before then.  This part of the test thus requires a check
429      * on the MPI major version */
430 #if MPI_VERSION >= 3
431     if (rank == 0) {
432         MTestPrintfMsg(1, "testing MPI_Dist_graph_create w/ no graph\n");
433     }
434     for (reorder = 0; reorder <= 1; reorder++) {
435         MPI_Dist_graph_create(MPI_COMM_WORLD, 0, sources, degrees,
436                               destinations, MPI_WEIGHTS_EMPTY, MPI_INFO_NULL, reorder, &comm);
437         MPI_Dist_graph_neighbors_count(comm, &check_indegree, &check_outdegree, &check_weighted);
438         if (!check_weighted) {
439             fprintf(stderr,
440                     "expected weighted == TRUE for the \"no graph -- MPI_WEIGHTS_EMPTY\" case\n");
441             ++errs;
442         }
443         MPI_Comm_free(&comm);
444     }
445 #endif
446
447     /* MPI_Dist_graph_create() with no graph -- passing NULLs instead */
448     if (rank == 0) {
449         MTestPrintfMsg(1, "testing MPI_Dist_graph_create w/ no graph -- NULLs\n");
450     }
451     for (reorder = 0; reorder <= 1; reorder++) {
452         MPI_Dist_graph_create(MPI_COMM_WORLD, 0, NULL, NULL,
453                               NULL, NULL, MPI_INFO_NULL, reorder, &comm);
454         MPI_Dist_graph_neighbors_count(comm, &check_indegree, &check_outdegree, &check_weighted);
455         /* ambiguous if they are equal, only check when they are distinct values. */
456         if (MPI_UNWEIGHTED != NULL) {
457             if (!check_weighted) {
458                 fprintf(stderr, "expected weighted == TRUE for the \"no graph -- NULLs\" case\n");
459                 ++errs;
460             }
461         }
462         MPI_Comm_free(&comm);
463     }
464
465     /* MPI_Dist_graph_create() with no graph -- passing NULLs+MPI_UNWEIGHTED instead */
466     if (rank == 0) {
467         MTestPrintfMsg(1, "testing MPI_Dist_graph_create w/ no graph -- NULLs+MPI_UNWEIGHTED\n");
468     }
469     for (reorder = 0; reorder <= 1; reorder++) {
470         MPI_Dist_graph_create(MPI_COMM_WORLD, 0, NULL, NULL,
471                               NULL, MPI_UNWEIGHTED, MPI_INFO_NULL, reorder, &comm);
472         MPI_Dist_graph_neighbors_count(comm, &check_indegree, &check_outdegree, &check_weighted);
473         /* ambiguous if they are equal, only check when they are distinct values. */
474         if (MPI_UNWEIGHTED != NULL) {
475             if (check_weighted) {
476                 fprintf(stderr,
477                         "expected weighted == FALSE for the \"no graph -- NULLs+MPI_UNWEIGHTED\" case\n");
478                 ++errs;
479             }
480         }
481         MPI_Comm_free(&comm);
482     }
483
484     /* MPI_Dist_graph_create_adjacent() with no graph */
485     if (rank == 0) {
486         MTestPrintfMsg(1, "testing MPI_Dist_graph_create_adjacent w/ no graph\n");
487     }
488     for (reorder = 0; reorder <= 1; reorder++) {
489         MPI_Dist_graph_create_adjacent(MPI_COMM_WORLD, 0, sources, sweights,
490                                        0, destinations, dweights, MPI_INFO_NULL, reorder, &comm);
491         MPI_Dist_graph_neighbors_count(comm, &check_indegree, &check_outdegree, &check_weighted);
492         if (!check_weighted) {
493             fprintf(stderr, "expected weighted == TRUE for the \"no graph\" case\n");
494             ++errs;
495         }
496         MPI_Comm_free(&comm);
497     }
498
499     /* MPI_Dist_graph_create_adjacent() with no graph -- passing MPI_WEIGHTS_EMPTY instead */
500     /* NOTE that MPI_WEIGHTS_EMPTY was added in MPI-3 and does not
501      * appear before then.  This part of the test thus requires a check
502      * on the MPI major version */
503 #if MPI_VERSION >= 3
504     if (rank == 0) {
505         MTestPrintfMsg(1,
506                        "testing MPI_Dist_graph_create_adjacent w/ no graph -- MPI_WEIGHTS_EMPTY\n");
507     }
508     for (reorder = 0; reorder <= 1; reorder++) {
509         MPI_Dist_graph_create_adjacent(MPI_COMM_WORLD, 0, sources, MPI_WEIGHTS_EMPTY,
510                                        0, destinations, MPI_WEIGHTS_EMPTY, MPI_INFO_NULL, reorder,
511                                        &comm);
512         MPI_Dist_graph_neighbors_count(comm, &check_indegree, &check_outdegree, &check_weighted);
513         if (!check_weighted) {
514             fprintf(stderr,
515                     "expected weighted == TRUE for the \"no graph -- MPI_WEIGHTS_EMPTY\" case\n");
516             ++errs;
517         }
518         MPI_Comm_free(&comm);
519     }
520 #endif
521
522     /* MPI_Dist_graph_create_adjacent() with no graph -- passing NULLs instead */
523     if (rank == 0) {
524         MTestPrintfMsg(1, "testing MPI_Dist_graph_create_adjacent w/ no graph -- NULLs\n");
525     }
526     for (reorder = 0; reorder <= 1; reorder++) {
527         MPI_Dist_graph_create_adjacent(MPI_COMM_WORLD, 0, NULL, NULL,
528                                        0, NULL, NULL, MPI_INFO_NULL, reorder, &comm);
529         MPI_Dist_graph_neighbors_count(comm, &check_indegree, &check_outdegree, &check_weighted);
530         /* ambiguous if they are equal, only check when they are distinct values. */
531         if (MPI_UNWEIGHTED != NULL) {
532             if (!check_weighted) {
533                 fprintf(stderr, "expected weighted == TRUE for the \"no graph -- NULLs\" case\n");
534                 ++errs;
535             }
536         }
537         MPI_Comm_free(&comm);
538     }
539
540     /* MPI_Dist_graph_create_adjacent() with no graph -- passing NULLs+MPI_UNWEIGHTED instead */
541     if (rank == 0) {
542         MTestPrintfMsg(1,
543                        "testing MPI_Dist_graph_create_adjacent w/ no graph -- NULLs+MPI_UNWEIGHTED\n");
544     }
545     for (reorder = 0; reorder <= 1; reorder++) {
546         MPI_Dist_graph_create_adjacent(MPI_COMM_WORLD, 0, NULL, MPI_UNWEIGHTED,
547                                        0, NULL, MPI_UNWEIGHTED, MPI_INFO_NULL, reorder, &comm);
548         MPI_Dist_graph_neighbors_count(comm, &check_indegree, &check_outdegree, &check_weighted);
549         /* ambiguous if they are equal, only check when they are distinct values. */
550         if (MPI_UNWEIGHTED != NULL) {
551             if (check_weighted) {
552                 fprintf(stderr,
553                         "expected weighted == FALSE for the \"no graph -- NULLs+MPI_UNWEIGHTED\" case\n");
554                 ++errs;
555             }
556         }
557         MPI_Comm_free(&comm);
558     }
559
560
561     for (i = 0; i < size; i++)
562         free(layout[i]);
563     free(layout);
564     free(sources);
565     free(sweights);
566     free(destinations);
567     free(dweights);
568     free(degrees);
569 #endif
570
571     MTest_Finalize(errs);
572
573     return MTestReturnValue(errs);
574 }