1 /* Copyright (c) 2013-2017. The SimGrid Team.
2 * All rights reserved. */
4 /* This program is free software; you can redistribute it and/or modify it
5 * under the terms of the license (GNU LGPL) which comes with this package. */
7 /*****************************************************************************
9 * Function: alltoall_bruck
14 send_buff: send input buffer
15 send_count: number of elements to send
16 send_type: data type of elements being sent
17 recv_buff: receive output buffer
18 recv_count: number of elements to received
19 recv_type: data type of elements being received
22 * Descrp: Function realizes the alltoall operation using the bruck algorithm.
24 * Auther: MPICH / modified by Ahmad Faraj
26 ****************************************************************************/
28 #include "../colls_private.hpp"
35 Coll_alltoall_bruck::alltoall(void *send_buff, int send_count,
36 MPI_Datatype send_type, void *recv_buff,
37 int recv_count, MPI_Datatype recv_type,
42 MPI_Datatype new_type;
44 int *blocks_length, *disps;
45 int i, src, dst, rank, num_procs, count, block, position;
46 int pack_size, tag = COLL_TAG_ALLTOALL, pof2 = 1;
50 char *send_ptr = (char *) send_buff;
51 char *recv_ptr = (char *) recv_buff;
53 num_procs = comm->size();
56 extent = recv_type->get_extent();
58 tmp_buff = (char *) smpi_get_tmp_sendbuffer(num_procs * recv_count * extent);
59 disps = (int *) xbt_malloc(sizeof(int) * num_procs);
60 blocks_length = (int *) xbt_malloc(sizeof(int) * num_procs);
62 Request::sendrecv(send_ptr + rank * send_count * extent,
63 (num_procs - rank) * send_count, send_type, rank, tag,
64 recv_ptr, (num_procs - rank) * recv_count, recv_type, rank,
67 Request::sendrecv(send_ptr, rank * send_count, send_type, rank, tag,
68 recv_ptr + (num_procs - rank) * recv_count * extent,
69 rank * recv_count, recv_type, rank, tag, comm, &status);
73 MPI_Pack_size(send_count * num_procs, send_type, comm, &pack_size);
75 while (pof2 < num_procs) {
76 dst = (rank + pof2) % num_procs;
77 src = (rank - pof2 + num_procs) % num_procs;
81 for (block = 1; block < num_procs; block++)
83 blocks_length[count] = send_count;
84 disps[count] = block * send_count;
88 MPI_Type_indexed(count, blocks_length, disps, recv_type, &new_type);
92 MPI_Pack(recv_buff, 1, new_type, tmp_buff, pack_size, &position, comm);
94 Request::sendrecv(tmp_buff, position, MPI_PACKED, dst, tag, recv_buff, 1,
95 new_type, src, tag, comm, &status);
96 Datatype::unref(new_type);
104 Request::sendrecv(recv_ptr + (rank + 1) * recv_count * extent,
105 (num_procs - rank - 1) * recv_count, send_type,
106 rank, tag, tmp_buff, (num_procs - rank - 1) * recv_count,
107 recv_type, rank, tag, comm, &status);
109 Request::sendrecv(recv_ptr, (rank + 1) * recv_count, send_type, rank, tag,
110 tmp_buff + (num_procs - rank - 1) * recv_count * extent,
111 (rank + 1) * recv_count, recv_type, rank, tag, comm, &status);
114 for (i = 0; i < num_procs; i++)
115 Request::sendrecv(tmp_buff + i * recv_count * extent, recv_count, send_type,
117 recv_ptr + (num_procs - i - 1) * recv_count * extent,
118 recv_count, recv_type, rank, tag, comm, &status);
120 smpi_free_tmp_buffer(tmp_buff);