3 XBT_LOG_NEW_DEFAULT_SUBCATEGORY(smpi_colls, smpi,
4 "Logging specific to SMPI collectives");
6 /*****************************************************************************
8 * Function: alltoall_2dmesh_shoot
13 send_buff: send input buffer
14 send_count: number of elements to send
15 send_type: data type of elements being sent
16 recv_buff: receive output buffer
17 recv_count: number of elements to received
18 recv_type: data type of elements being received
21 * Descrp: Function realizes the alltoall operation using the 2dmesh
22 algorithm. It actually performs allgather operation in x dimension
23 then in the y dimension. Each node then extracts the needed data.
24 The communication in each dimension follows "simple."
28 ****************************************************************************/
29 int alltoall_check_is_2dmesh(int num, int * i, int * j)
56 smpi_coll_tuned_alltoall_2dmesh(void * send_buff, int send_count, MPI_Datatype send_type,
57 void * recv_buff, int recv_count, MPI_Datatype recv_type,
60 MPI_Status * statuses, s;
61 MPI_Request * reqs, * req_ptr;;
64 char * tmp_buff1, * tmp_buff2;
65 int i, j, src, dst, rank, num_procs, count, num_reqs;
66 int rows, cols, my_row, my_col, X, Y, send_offset, recv_offset;
67 int two_dsize, my_row_base, my_col_base, src_row_base, block_size;
68 int tag = 1, failure = 0, success = 1;
70 MPI_Comm_rank(comm, &rank);
71 MPI_Comm_size(comm, &num_procs);
72 MPI_Type_extent(send_type, &extent);
74 if (!alltoall_check_is_2dmesh(num_procs, &X, &Y))
79 my_row_base = (rank / Y) * Y;
80 my_col_base = rank % Y;
82 block_size = extent * send_count;
84 tmp_buff1 =(char *) malloc(block_size * num_procs * Y);
87 XBT_DEBUG("alltoall-2dmesh_shoot.c:88: cannot allocate memory");
92 tmp_buff2 =(char *) malloc(block_size * Y);
95 XBT_WARN("alltoall-2dmesh_shoot.c:88: cannot allocate memory");
103 if (Y > X) num_reqs = Y;
105 statuses = (MPI_Status *) malloc(num_reqs * sizeof(MPI_Status));
106 reqs = (MPI_Request *) malloc(num_reqs * sizeof(MPI_Request));
109 XBT_WARN("alltoall-2dmesh_shoot.c:88: cannot allocate memory");
116 send_offset = recv_offset = (rank % Y) * block_size * num_procs;
118 count = send_count * num_procs;
120 for (i = 0; i < Y; i++)
122 src = i + my_row_base;
126 recv_offset = (src % Y) * block_size * num_procs;
127 MPI_Irecv(tmp_buff1 + recv_offset, count, recv_type, src, tag, comm,
131 for (i = 0; i < Y; i++)
133 dst = i + my_row_base;
136 MPI_Send(send_buff, count, send_type, dst, tag, comm);
139 MPI_Waitall(Y - 1, reqs, statuses);
142 for (i = 0; i < Y; i++)
144 send_offset = (rank * block_size) + (i * block_size * num_procs);
145 recv_offset = (my_row_base * block_size) + (i * block_size);
147 if (i + my_row_base == rank)
148 MPI_Sendrecv (send_buff + recv_offset, send_count, send_type,
149 rank, tag, recv_buff + recv_offset, recv_count,
150 recv_type, rank, tag, comm, &s);
153 MPI_Sendrecv (tmp_buff1 + send_offset, send_count, send_type,
155 recv_buff + recv_offset, recv_count, recv_type,
156 rank, tag, comm, &s);
160 for (i = 0; i < X; i++)
162 src = (i * Y + my_col_base);
165 src_row_base = (src / Y) * Y;
167 MPI_Irecv(recv_buff + src_row_base * block_size, recv_count * Y,
168 recv_type, src, tag, comm, req_ptr++);
171 for (i = 0; i < X; i++)
173 dst = (i * Y + my_col_base);
178 for (j = 0; j < Y; j++)
180 send_offset = (dst + j * num_procs) * block_size;
182 if (j + my_row_base == rank)
183 MPI_Sendrecv (send_buff + dst * block_size, send_count, send_type,
185 tmp_buff2 + recv_offset, recv_count, recv_type,
186 rank, tag, comm, &s);
188 MPI_Sendrecv (tmp_buff1 + send_offset, send_count, send_type,
190 tmp_buff2 + recv_offset, recv_count, recv_type,
191 rank, tag, comm, &s);
193 recv_offset += block_size;
196 MPI_Send(tmp_buff2, send_count * Y, send_type, dst, tag, comm);
198 MPI_Waitall(X - 1, reqs, statuses);