72 size_t *num_src_intersections,
74 size_t *num_dst_intersections,
90 int comm_size, rank, is_inter;
92 xt_mpi_call(MPI_Comm_test_inter(comm, &is_inter), comm);
93 int (*get_comm_size)(
MPI_Comm,
int *)
94 = is_inter ? MPI_Comm_remote_size : MPI_Comm_size;
104 size_t size_sum = src_pack_size + dst_pack_size;
106 if (size_sum >= INT_MAX || size_sum < src_pack_size
107 || size_sum < dst_pack_size)
108 die(
"local src+dst index lists are too large");
110 int send_buffer_size = (int)size_sum;
113 int *restrict pack_sizes
114 =
xmalloc((
size_t)comm_size *
sizeof(*pack_sizes) * 2);
116 xt_mpi_call(MPI_Allgather(&send_buffer_size, 1, MPI_INT,
117 pack_sizes, 1, MPI_INT, comm), comm);
119 int *restrict displ = pack_sizes + comm_size;
120 unsigned recv_buffer_size = 0, size_overflow = 0;
121 for (
size_t i = 0; i < (size_t)comm_size; ++i) {
122 displ[i] = (int)recv_buffer_size;
123 recv_buffer_size += (unsigned)pack_sizes[i];
124 size_overflow |= recv_buffer_size & (1U << (
sizeof(int) * CHAR_BIT - 1));
127 die(
"accumulated buffer sizes too big,"
128 " use distributed directory (xt_xmap_dist_dir_new)!");
129 void *recv_buffer =
xmalloc((
size_t)recv_buffer_size + size_sum),
130 *send_buffer = (
unsigned char *)recv_buffer + (
size_t)recv_buffer_size;
140 xt_mpi_call(MPI_Allgatherv(send_buffer, send_buffer_size, MPI_PACKED,
141 recv_buffer, pack_sizes, displ, MPI_PACKED,
144 size_t dst_isect_count = 0, src_isect_count = 0;
145 int large_list_seen = 0;
147 for (
int i = 0; i < comm_size; ++i) {
152 if (is_inter || i !=
rank) {
154 pack_sizes[i], &position, comm);
156 pack_sizes[i], &position, comm);
158 src = src_idxlist_local;
159 dst = dst_idxlist_local;
161 large_list_seen = large_list_seen
168 dsti[dst_isect_count].
list = intersect;
169 dsti[dst_isect_count].
rank = i;
179 srci[src_isect_count].
list = intersect;
180 srci[src_isect_count].
rank = i;
185 if (is_inter || i !=
rank) {
191 *stripify = large_list_seen;
197 *num_src_intersections = src_isect_count;
198 srci =
xrealloc(srci, src_isect_count *
sizeof (**src_intersections));
199 *src_intersections = srci;
201 *num_dst_intersections = dst_isect_count;
202 dsti =
xrealloc(dsti, dst_isect_count *
sizeof (**dst_intersections));
203 *dst_intersections = dsti;
217 INSTR_DEF(t_xt_xmap_all2all_new,
"xt_xmap_all2all_new")
226 struct Xt_com_list * src_intersections = NULL, * dst_intersections = NULL;
227 size_t num_src_intersections, num_dst_intersections;
232 &dst_intersections, &num_dst_intersections,
233 &stripify, src_idxlist, dst_idxlist, newcomm,
236 Xt_xmap (*xmap_new)(
int num_src_intersections,
238 int num_dst_intersections,
244 Xt_xmap xmap = xmap_new((
int)num_src_intersections, src_intersections,
245 (
int)num_dst_intersections, dst_intersections,
246 src_idxlist, dst_idxlist, newcomm);
249 for (
size_t i = 0; i < num_src_intersections; ++i)
250 if (src_intersections[i].
list != NULL)
252 for (
size_t i = 0; i < num_dst_intersections; ++i)
253 if (dst_intersections[i].
list != NULL)
255 free(src_intersections);
256 free(dst_intersections);
static void exchange_idxlists(struct Xt_com_list **src_intersections, size_t *num_src_intersections, struct Xt_com_list **dst_intersections, size_t *num_dst_intersections, int *stripify, Xt_idxlist src_idxlist_local, Xt_idxlist dst_idxlist_local, MPI_Comm comm, Xt_config config)
Xt_xmap xt_xmap_intersection_ext_new(int num_src_intersections, const struct Xt_com_list src_com[num_src_intersections], int num_dst_intersections, const struct Xt_com_list dst_com[num_dst_intersections], Xt_idxlist src_idxlist, Xt_idxlist dst_idxlist, MPI_Comm comm)
Xt_xmap xt_xmap_intersection_new(int num_src_intersections, const struct Xt_com_list src_com[num_src_intersections], int num_dst_intersections, const struct Xt_com_list dst_com[num_dst_intersections], Xt_idxlist src_idxlist, Xt_idxlist dst_idxlist, MPI_Comm comm)