[svn-r8400] intercomm Ibarrier modernization and bugfix
[mpich.git] / src / mpi / coll / ibarrier.c
1 /* -*- Mode: C; c-basic-offset:4 ; indent-tabs-mode:nil ; -*- */
2 /*
3  *  (C) 2010 by Argonne National Laboratory.
4  *      See COPYRIGHT in top-level directory.
5  */
6
7 #include "mpiimpl.h"
8
9 /* -- Begin Profiling Symbol Block for routine MPIX_Ibarrier */
10 #if defined(HAVE_PRAGMA_WEAK)
11 #pragma weak MPIX_Ibarrier = PMPIX_Ibarrier
12 #elif defined(HAVE_PRAGMA_HP_SEC_DEF)
13 #pragma _HP_SECONDARY_DEF PMPIX_Ibarrier  MPIX_Ibarrier
14 #elif defined(HAVE_PRAGMA_CRI_DUP)
15 #pragma _CRI duplicate MPIX_Ibarrier as PMPIX_Ibarrier
16 #endif
17 /* -- End Profiling Symbol Block */
18
19 /* Define MPICH_MPI_FROM_PMPI if weak symbols are not supported to build
20    the MPI routines */
21 #ifndef MPICH_MPI_FROM_PMPI
22 #undef MPIX_Ibarrier
23 #define MPIX_Ibarrier PMPIX_Ibarrier
24
25 /* any non-MPI functions go here, especially non-static ones */
26
27 /* This is the default implementation of the barrier operation.  The
28    algorithm is:
29
30    Algorithm: MPI_Ibarrier
31
32    We use the dissemination algorithm described in:
33    Debra Hensgen, Raphael Finkel, and Udi Manber, "Two Algorithms for
34    Barrier Synchronization," International Journal of Parallel
35    Programming, 17(1):1-17, 1988.
36
37    It uses ceiling(lgp) steps. In step k, 0 <= k <= (ceiling(lgp)-1),
38    process i sends to process (i + 2^k) % p and receives from process
39    (i - 2^k + p) % p.
40
41    Possible improvements:
42
43    End Algorithm: MPI_Ibarrier
44
45    This is an intracommunicator barrier only!
46 */
47 /* Provides a generic "flat" barrier that doesn't know anything about hierarchy.
48  * It will choose between several different algorithms based on the given
49  * parameters. */
50 #undef FUNCNAME
51 #define FUNCNAME MPIR_Ibarrier_intra
52 #undef FCNAME
53 #define FCNAME MPIU_QUOTE(FUNCNAME)
54 int MPIR_Ibarrier_intra(MPID_Comm *comm_ptr, MPID_Sched_t s)
55 {
56     int mpi_errno = MPI_SUCCESS;
57     int size, rank, src, dst, mask;
58
59     MPIU_Assert(comm_ptr->comm_kind == MPID_INTRACOMM);
60
61     size = comm_ptr->local_size;
62     rank = comm_ptr->rank;
63
64     /* Trivial barriers return immediately */
65     if (size == 1) goto fn_exit;
66
67     mask = 0x1;
68     while (mask < size) {
69         dst = (rank + mask) % size;
70         src = (rank - mask + size) % size;
71
72         mpi_errno = MPID_Sched_send(NULL, 0, MPI_BYTE, dst, comm_ptr, s);
73         if (mpi_errno) MPIU_ERR_POP(mpi_errno);
74
75         mpi_errno = MPID_Sched_recv(NULL, 0, MPI_BYTE, src, comm_ptr, s);
76         if (mpi_errno) MPIU_ERR_POP(mpi_errno);
77
78         mpi_errno = MPID_Sched_barrier(s);
79         if (mpi_errno) MPIU_ERR_POP(mpi_errno);
80
81         mask <<= 1;
82     }
83
84 fn_exit:
85     return mpi_errno;
86 fn_fail:
87     goto fn_exit;
88 }
89
90 /* Provides a generic "flat" barrier that doesn't know anything about hierarchy.
91  * It will choose between several different algorithms based on the given
92  * parameters. */
93 #undef FUNCNAME
94 #define FUNCNAME MPIR_Ibarrier_inter
95 #undef FCNAME
96 #define FCNAME MPIU_QUOTE(FUNCNAME)
97 int MPIR_Ibarrier_inter(MPID_Comm *comm_ptr, MPID_Sched_t s)
98 {
99     int mpi_errno = MPI_SUCCESS;
100     int size, rank, root;
101     MPIR_SCHED_CHKPMEM_DECL(1);
102     char *buf = NULL;
103
104     MPIU_Assert(comm_ptr->comm_kind == MPID_INTERCOMM);
105
106     size = comm_ptr->local_size;
107     rank = comm_ptr->rank;
108
109     /* Get the local intracommunicator */
110     if (!comm_ptr->local_comm) {
111         mpi_errno = MPIR_Setup_intercomm_localcomm(comm_ptr);
112         if (mpi_errno) MPIU_ERR_POP(mpi_errno);
113     }
114
115     /* do a barrier on the local intracommunicator */
116     MPIU_Assert(comm_ptr->local_comm->coll_fns && comm_ptr->local_comm->coll_fns->Ibarrier);
117     mpi_errno = comm_ptr->local_comm->coll_fns->Ibarrier(comm_ptr->local_comm, s);
118     if (mpi_errno) MPIU_ERR_POP(mpi_errno);
119
120     MPID_SCHED_BARRIER(s);
121
122     /* rank 0 on each group does an intercommunicator broadcast to the
123        remote group to indicate that all processes in the local group
124        have reached the barrier. We do a 1-byte bcast because a 0-byte
125        bcast will just return without doing anything. */
126
127     MPIR_SCHED_CHKPMEM_MALLOC(buf, char *, 1, mpi_errno, "bcast buf");
128     buf[0] = 'D'; /* avoid valgrind warnings */
129
130     /* first broadcast from left to right group, then from right to
131        left group */
132     MPIU_Assert(comm_ptr->coll_fns && comm_ptr->coll_fns->Ibcast);
133     if (comm_ptr->is_low_group) {
134         root = (rank == 0) ? MPI_ROOT : MPI_PROC_NULL;
135         mpi_errno = comm_ptr->coll_fns->Ibcast(buf, 1, MPI_BYTE, root, comm_ptr, s);
136         if (mpi_errno) MPIU_ERR_POP(mpi_errno);
137
138         MPID_SCHED_BARRIER(s);
139
140         /* receive bcast from right */
141         root = 0;
142         mpi_errno = comm_ptr->coll_fns->Ibcast(buf, 1, MPI_BYTE, root, comm_ptr, s);
143         if (mpi_errno) MPIU_ERR_POP(mpi_errno);
144     }
145     else {
146         /* receive bcast from left */
147         root = 0;
148         mpi_errno = comm_ptr->coll_fns->Ibcast(buf, 1, MPI_BYTE, root, comm_ptr, s);
149         if (mpi_errno) MPIU_ERR_POP(mpi_errno);
150
151         MPID_SCHED_BARRIER(s);
152
153         /* bcast to left */
154         root = (rank == 0) ? MPI_ROOT : MPI_PROC_NULL;
155         mpi_errno = comm_ptr->coll_fns->Ibcast(buf, 1, MPI_BYTE, root, comm_ptr, s);
156         if (mpi_errno) MPIU_ERR_POP(mpi_errno);
157     }
158
159     MPIR_SCHED_CHKPMEM_COMMIT(s);
160 fn_exit:
161     return mpi_errno;
162 fn_fail:
163     MPIR_SCHED_CHKPMEM_REAP(s);
164     goto fn_exit;
165 }
166
167 #undef FUNCNAME
168 #define FUNCNAME MPIR_Ibarrier_impl
169 #undef FCNAME
170 #define FCNAME MPIU_QUOTE(FUNCNAME)
171 int MPIR_Ibarrier_impl(MPID_Comm *comm_ptr, MPI_Request *request)
172 {
173     int mpi_errno = MPI_SUCCESS;
174     int tag = -1;
175     MPID_Request *reqp = NULL;
176     MPID_Sched_t s = MPID_SCHED_NULL;
177
178     *request = MPI_REQUEST_NULL;
179
180     mpi_errno = MPID_Sched_next_tag(comm_ptr, &tag);
181     if (mpi_errno) MPIU_ERR_POP(mpi_errno);
182     mpi_errno = MPID_Sched_create(&s);
183     if (mpi_errno) MPIU_ERR_POP(mpi_errno);
184
185     MPIU_Assert(comm_ptr->coll_fns != NULL);
186     MPIU_Assert(comm_ptr->coll_fns->Ibarrier != NULL);
187     mpi_errno = comm_ptr->coll_fns->Ibarrier(comm_ptr, s);
188     if (mpi_errno) MPIU_ERR_POP(mpi_errno);
189
190     mpi_errno = MPID_Sched_start(&s, comm_ptr, tag, &reqp);
191     if (reqp)
192         *request = reqp->handle;
193     if (mpi_errno) MPIU_ERR_POP(mpi_errno);
194
195 fn_exit:
196     return mpi_errno;
197 fn_fail:
198     goto fn_exit;
199 }
200
201 #endif /* MPICH_MPI_FROM_PMPI */
202
203 #undef FUNCNAME
204 #define FUNCNAME MPIX_Ibarrier
205 #undef FCNAME
206 #define FCNAME MPIU_QUOTE(FUNCNAME)
207 /*@
208 MPIX_Ibarrier - XXX description here
209
210 Input Parameters:
211 . comm - communicator (handle)
212
213 Output Parameters:
214 . request - communication request (handle)
215
216 .N ThreadSafe
217
218 .N Fortran
219
220 .N Errors
221 @*/
222 int MPIX_Ibarrier(MPI_Comm comm, MPI_Request *request)
223 {
224     int mpi_errno = MPI_SUCCESS;
225     MPID_Comm *comm_ptr = NULL;
226     MPID_MPI_STATE_DECL(MPID_STATE_MPIX_IBARRIER);
227
228     MPIU_THREAD_CS_ENTER(ALLFUNC,);
229     MPID_MPI_FUNC_ENTER(MPID_STATE_MPIX_IBARRIER);
230
231     /* Validate parameters, especially handles needing to be converted */
232 #   ifdef HAVE_ERROR_CHECKING
233     {
234         MPID_BEGIN_ERROR_CHECKS
235         {
236             MPIR_ERRTEST_COMM(comm, mpi_errno);
237
238             /* TODO more checks may be appropriate */
239             if (mpi_errno != MPI_SUCCESS) goto fn_fail;
240         }
241         MPID_END_ERROR_CHECKS
242     }
243 #   endif /* HAVE_ERROR_CHECKING */
244
245     /* Convert MPI object handles to object pointers */
246     MPID_Comm_get_ptr(comm, comm_ptr);
247
248     /* Validate parameters and objects (post conversion) */
249 #   ifdef HAVE_ERROR_CHECKING
250     {
251         MPID_BEGIN_ERROR_CHECKS
252         {
253             MPID_Comm_valid_ptr(comm_ptr, mpi_errno);
254             MPIR_ERRTEST_ARGNULL(request,"request", mpi_errno);
255             /* TODO more checks may be appropriate (counts, in_place, buffer aliasing, etc) */
256             if (mpi_errno != MPI_SUCCESS) goto fn_fail;
257         }
258         MPID_END_ERROR_CHECKS
259     }
260 #   endif /* HAVE_ERROR_CHECKING */
261
262     /* ... body of routine ...  */
263
264     mpi_errno = MPIR_Ibarrier_impl(comm_ptr, request);
265     if (mpi_errno) MPIU_ERR_POP(mpi_errno);
266
267     /* ... end of body of routine ... */
268
269 fn_exit:
270     MPID_MPI_FUNC_EXIT(MPID_STATE_MPIX_IBARRIER);
271     MPIU_THREAD_CS_EXIT(ALLFUNC,);
272     return mpi_errno;
273
274 fn_fail:
275     /* --BEGIN ERROR HANDLING-- */
276 #   ifdef HAVE_ERROR_CHECKING
277     {
278         mpi_errno = MPIR_Err_create_code(
279             mpi_errno, MPIR_ERR_RECOVERABLE, FCNAME, __LINE__, MPI_ERR_OTHER,
280             "**mpix_ibarrier", "**mpix_ibarrier %C %p", comm, request);
281     }
282 #   endif
283     mpi_errno = MPIR_Err_return_comm(comm_ptr, FCNAME, mpi_errno);
284     goto fn_exit;
285     /* --END ERROR HANDLING-- */
286     goto fn_exit;
287 }