2 This source code file is part of thread_mpi.
3 Written by Sander Pronk, Erik Lindahl, and possibly others.
5 Copyright (c) 2009, Sander Pronk, Erik Lindahl.
8 Redistribution and use in source and binary forms, with or without
9 modification, are permitted provided that the following conditions are met:
10 1) Redistributions of source code must retain the above copyright
11 notice, this list of conditions and the following disclaimer.
12 2) Redistributions in binary form must reproduce the above copyright
13 notice, this list of conditions and the following disclaimer in the
14 documentation and/or other materials provided with the distribution.
15 3) Neither the name of the copyright holders nor the
16 names of its contributors may be used to endorse or promote products
17 derived from this software without specific prior written permission.
19 THIS SOFTWARE IS PROVIDED BY US ''AS IS'' AND ANY
20 EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
21 WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
22 DISCLAIMED. IN NO EVENT SHALL WE BE LIABLE FOR ANY
23 DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
24 (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
25 LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
26 ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
27 (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
28 SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
30 If you want to redistribute modifications, please consider that
31 scientific software is very special. Version control is crucial -
32 bugs must be traceable. We will be happy to consider code for
33 inclusion in the official distribution, but derived work should not
34 be called official thread_mpi. Details are found in the README & COPYING
38 #ifdef HAVE_TMPI_CONFIG_H
39 #include "tmpi_config.h"
58 #include "collective.h"
62 int tMPI_Scan(void* sendbuf, void* recvbuf, int count,
63 tMPI_Datatype datatype, tMPI_Op op, tMPI_Comm comm)
65 struct tmpi_thread *cur=tMPI_Get_current();
66 int myrank=tMPI_Comm_seek_rank(comm, cur);
67 int N=tMPI_Comm_N(comm);
68 int prev=myrank - 1; /* my previous neighbor */
69 int next=myrank + 1; /* my next neighbor */
72 tMPI_Profile_count_start(cur);
75 tMPI_Trace_print("tMPI_Scan(%p, %p, %d, %p, %p, %p)",
76 sendbuf, recvbuf, count, datatype, op, comm);
82 return tMPI_Error(comm, TMPI_ERR_BUF);
84 if (sendbuf==TMPI_IN_PLACE)
89 /* we set our send and recv buffers */
90 tMPI_Atomic_ptr_set(&(comm->reduce_sendbuf[myrank]),sendbuf);
91 tMPI_Atomic_ptr_set(&(comm->reduce_recvbuf[myrank]),recvbuf);
93 /* now wait for the previous rank to finish */
99 #if defined(TMPI_PROFILE) && defined(TMPI_CYCLE_COUNT)
100 tMPI_Profile_wait_start(cur);
102 /* wait for the previous neighbor's data to be ready */
103 tMPI_Event_wait( &(comm->csync[myrank].events[prev]) );
104 tMPI_Event_process( &(comm->csync[myrank].events[prev]), 1);
105 #if defined(TMPI_PROFILE) && defined(TMPI_CYCLE_COUNT)
106 tMPI_Profile_wait_stop(cur, TMPIWAIT_Reduce);
109 printf("%d: scanning with %d \n", myrank, prev, iteration);
112 /* now do the reduction */
115 a = (void*)tMPI_Atomic_ptr_get(&(comm->reduce_recvbuf[prev]));
119 a = (void*)tMPI_Atomic_ptr_get(&(comm->reduce_sendbuf[prev]));
123 if ((ret=tMPI_Reduce_run_op(recvbuf, a, b, datatype,
124 count, op, comm)) != TMPI_SUCCESS)
129 /* signal to my previous neighbor that I'm done with the data */
130 tMPI_Event_signal( &(comm->csync[prev].events[prev]) );
134 if (sendbuf != recvbuf)
136 /* copy the data if this is rank 0, and not MPI_IN_PLACE */
137 memcpy(recvbuf, sendbuf, count*datatype->size);
143 /* signal to my next neighbor that I have the data */
144 tMPI_Event_signal( &(comm->csync[next].events[myrank]) );
145 /* and wait for my next neighbor to finish */
146 tMPI_Event_wait( &(comm->csync[myrank].events[myrank]) );
147 tMPI_Event_process( &(comm->csync[myrank].events[myrank]), 1);
151 #if defined(TMPI_PROFILE) && defined(TMPI_CYCLE_COUNT)
152 tMPI_Profile_wait_start(cur);
154 /*tMPI_Barrier_wait( &(comm->barrier));*/
155 #if defined(TMPI_PROFILE)
156 /*tMPI_Profile_wait_stop(cur, TMPIWAIT_Reduce);*/
157 tMPI_Profile_count_stop(cur, TMPIFN_Scan);