1 /*
2  * Copyright (C) by Argonne National Laboratory
3  *     See COPYRIGHT in top-level directory
4  */
5 
6 /* This test tries to overlap multiple Comm_idups with other communicator
7   creation functions, either intracomm or intercomm.
8  */
9 
10 #include <stdio.h>
11 #include <stdlib.h>
12 #include <mpi.h>
13 #include "mpitest.h"
14 #include "mpithreadtest.h"
15 
16 #define NUM_THREADS 4
17 #define NUM_IDUPS   5
18 
19 MPI_Comm comms[NUM_THREADS];
20 int errs[NUM_THREADS] = { 0 };
21 
22 int verbose = 0;
23 
test_idup(void * arg)24 MTEST_THREAD_RETURN_TYPE test_idup(void *arg)
25 {
26     int i;
27     int size, rank;
28     int ranges[1][3];
29     int rleader, isLeft;
30     int *excl = NULL;
31     int tid = *(int *) arg;
32 
33     MPI_Group ingroup, high_group, even_group;
34     MPI_Comm local_comm, inter_comm;
35     MPI_Comm idupcomms[NUM_IDUPS];
36     MPI_Request reqs[NUM_IDUPS];
37 
38     MPI_Comm outcomm;
39     MPI_Comm incomm = comms[tid];
40 
41     MPI_Comm_size(incomm, &size);
42     MPI_Comm_rank(incomm, &rank);
43     MPI_Comm_group(incomm, &ingroup);
44 
45     /* Idup incomm multiple times */
46     for (i = 0; i < NUM_IDUPS; i++) {
47         MPI_Comm_idup(incomm, &idupcomms[i], &reqs[i]);
48     }
49 
50     /* Overlap pending idups with various comm generation functions */
51     /* Comm_dup */
52     MPI_Comm_dup(incomm, &outcomm);
53     errs[tid] += MTestTestComm(outcomm);
54     MTestFreeComm(&outcomm);
55 
56     /* Comm_split */
57     MPI_Comm_split(incomm, rank % 2, size - rank, &outcomm);
58     errs[tid] += MTestTestComm(outcomm);
59     MTestFreeComm(&outcomm);
60 
61     /* Comm_create, high half of incomm */
62     ranges[0][0] = size / 2;
63     ranges[0][1] = size - 1;
64     ranges[0][2] = 1;
65     MPI_Group_range_incl(ingroup, 1, ranges, &high_group);
66     MPI_Comm_create(incomm, high_group, &outcomm);
67     MPI_Group_free(&high_group);
68     errs[tid] += MTestTestComm(outcomm);
69     MTestFreeComm(&outcomm);
70 
71     /* Comm_create_group, even ranks of incomm */
72     /* exclude the odd ranks */
73     excl = malloc((size / 2) * sizeof(int));
74     for (i = 0; i < size / 2; i++)
75         excl[i] = (2 * i) + 1;
76 
77     MPI_Group_excl(ingroup, size / 2, excl, &even_group);
78     free(excl);
79 
80     if (rank % 2 == 0) {
81         MPI_Comm_create_group(incomm, even_group, 0, &outcomm);
82     } else {
83         outcomm = MPI_COMM_NULL;
84     }
85     MPI_Group_free(&even_group);
86     errs[tid] += MTestTestComm(outcomm);
87     MTestFreeComm(&outcomm);
88 
89     /* Intercomm_create & Intercomm_merge */
90     MPI_Comm_split(incomm, (rank < size / 2), rank, &local_comm);
91     if (rank == 0) {
92         rleader = size / 2;
93     } else if (rank == size / 2) {
94         rleader = 0;
95     } else {
96         rleader = -1;
97     }
98     isLeft = rank < size / 2;
99 
100     MPI_Intercomm_create(local_comm, 0, incomm, rleader, 99, &inter_comm);
101     MPI_Intercomm_merge(inter_comm, isLeft, &outcomm);
102     MPI_Comm_free(&local_comm);
103 
104     errs[tid] += MTestTestComm(inter_comm);
105     MTestFreeComm(&inter_comm);
106     errs[tid] += MTestTestComm(outcomm);
107     MTestFreeComm(&outcomm);
108 
109     MPI_Waitall(NUM_IDUPS, reqs, MPI_STATUSES_IGNORE);
110     for (i = 0; i < NUM_IDUPS; i++) {
111         errs[tid] += MTestTestComm(idupcomms[i]);
112         MPI_Comm_free(&idupcomms[i]);
113     }
114     MPI_Group_free(&ingroup);
115     return NULL;
116 }
117 
main(int argc,char ** argv)118 int main(int argc, char **argv)
119 {
120     int thread_args[NUM_THREADS];
121     int i, provided;
122     int toterrs = 0;
123     int size;
124 
125     MTest_Init_thread(&argc, &argv, MPI_THREAD_MULTIPLE, &provided);
126     MPI_Comm_size(MPI_COMM_WORLD, &size);
127 
128     if (provided < MPI_THREAD_MULTIPLE) {
129         printf("MPI_THREAD_MULTIPLE for the test\n");
130         MPI_Abort(MPI_COMM_WORLD, 1);
131     }
132     if (size < 2) {
133         printf("This test requires at least two processes\n");
134         MPI_Abort(MPI_COMM_WORLD, 1);
135     }
136 
137     for (i = 0; i < NUM_THREADS; i++) {
138         MPI_Comm_dup(MPI_COMM_WORLD, &comms[i]);
139     }
140     for (i = 0; i < NUM_THREADS; i++) {
141         thread_args[i] = i;
142         MTest_Start_thread(test_idup, (void *) &thread_args[i]);
143     }
144     MTest_Join_threads();
145 
146     for (i = 0; i < NUM_THREADS; i++) {
147         MPI_Comm_free(&comms[i]);
148         toterrs += errs[i];
149     }
150     MTest_Finalize(toterrs);
151     return MTestReturnValue(toterrs);
152 }
153