Yet Another eXchange Tool  0.9.0
xt_redist.c
Go to the documentation of this file.
1 
12 /*
13  * Keywords:
14  * Maintainer: Jörg Behrens <behrens@dkrz.de>
15  * Moritz Hanke <hanke@dkrz.de>
16  * Thomas Jahns <jahns@dkrz.de>
17  * URL: https://doc.redmine.dkrz.de/yaxt/html/
18  *
19  * Redistribution and use in source and binary forms, with or without
20  * modification, are permitted provided that the following conditions are
21  * met:
22  *
23  * Redistributions of source code must retain the above copyright notice,
24  * this list of conditions and the following disclaimer.
25  *
26  * Redistributions in binary form must reproduce the above copyright
27  * notice, this list of conditions and the following disclaimer in the
28  * documentation and/or other materials provided with the distribution.
29  *
30  * Neither the name of the DKRZ GmbH nor the names of its contributors
31  * may be used to endorse or promote products derived from this software
32  * without specific prior written permission.
33  *
34  * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS
35  * IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
36  * TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A
37  * PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER
38  * OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
39  * EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
40  * PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
41  * PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
42  * LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
43  * NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
44  * SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
45  */
46 #ifdef HAVE_CONFIG_H
47 #include <config.h>
48 #endif
49 
50 #include <limits.h>
51 #include <stdbool.h>
52 #include <stdlib.h>
53 
54 #include "core/core.h"
55 #include "xt/xt_core.h"
56 #include "xt/xt_redist.h"
57 #include "xt/xt_mpi.h"
58 #include "xt/xt_request.h"
59 #include "xt/xt_sort.h"
60 #include "core/ppm_xfuncs.h"
61 #include "xt_redist_internal.h"
62 
64 
65  return redist->vtable->copy(redist);
66 }
67 
69 
70  redist->vtable->delete(redist);
71 }
72 
73 void xt_redist_s_exchange(Xt_redist redist, int num_arrays,
74  const void **src_data, void **dst_data) {
75 
76  redist->vtable->s_exchange(redist, num_arrays, src_data, dst_data);
77 }
78 
79 void xt_redist_a_exchange(Xt_redist redist, int num_arrays,
80  const void **src_data, void **dst_data,
81  Xt_request *request) {
82 
83  redist->vtable->a_exchange(redist, num_arrays, src_data, dst_data, request);
84 }
85 
86 void xt_redist_s_exchange1(Xt_redist redist, const void *src_data, void *dst_data) {
87 
88  redist->vtable->s_exchange1(redist, src_data, dst_data);
89 }
90 
91 void xt_redist_a_exchange1(Xt_redist redist, const void *src_data,
92  void *dst_data, Xt_request *request) {
93 
94  redist->vtable->a_exchange1(redist, src_data, dst_data, request);
95 }
96 
98 
99  return redist->vtable->get_num_msg(redist, SEND);
100 }
101 
103 
104  return redist->vtable->get_num_msg(redist, RECV);
105 }
106 
107 MPI_Datatype xt_redist_get_send_MPI_Datatype(Xt_redist redist, int rank) {
108 
109  return redist->vtable->get_msg_MPI_Datatype(redist, rank, SEND);
110 }
111 
112 MPI_Datatype xt_redist_get_recv_MPI_Datatype(Xt_redist redist, int rank) {
113 
114  return redist->vtable->get_msg_MPI_Datatype(redist, rank, RECV);
115 }
116 
117 MPI_Datatype xt_redist_get_MPI_Datatype(Xt_redist redist, int rank,
118  enum xt_msg_direction direction)
119 {
120  return redist->vtable->get_msg_MPI_Datatype(redist, rank, direction);
121 }
122 
124 
125  return redist->vtable->get_MPI_Comm(redist);
126 }
127 
129  int *restrict *ranks)
130 {
131  return redist->vtable->get_msg_ranks(redist, direction, ranks);
132 }
133 
134 
135 void
137  const struct Xt_redist_msg *restrict src,
138  size_t src_stride,
139  struct Xt_redist_msg *restrict dst,
140  size_t dst_stride,
141  MPI_Comm comm) {
142 
143 
144  const unsigned char *restrict src_store = (const unsigned char *)src;
145  unsigned char *restrict dst_store = (unsigned char *)dst;
146  for (size_t i = 0; i < n; ++i) {
147  const struct Xt_redist_msg *restrict src_msg
148  = (const struct Xt_redist_msg *)(const void *)(src_store + i * src_stride);
149  struct Xt_redist_msg *dst_msg
150  = (struct Xt_redist_msg *)(void *)(dst_store + i * dst_stride);
151  dst_msg->rank = src_msg->rank;
152  xt_mpi_call(MPI_Type_dup(src_msg->datatype, &(dst_msg->datatype)), comm);
153  }
154 }
155 
156 void xt_redist_msgs_strided_destruct(size_t n, struct Xt_redist_msg *msgs,
157  MPI_Comm comm, size_t ofs_stride) {
158 
159  unsigned char *restrict msgs_store = (unsigned char *)msgs;
160  for (size_t i = 0; i < n; ++i) {
161  MPI_Datatype *dt
162  = &(((struct Xt_redist_msg *)(void *)(msgs_store + i * ofs_stride))->datatype);
163  if (*dt != MPI_DATATYPE_NULL)
164  xt_mpi_call(MPI_Type_free(dt), comm);
165  }
166 }
167 
168 void
169 xt_redist_check_comms(Xt_redist *redists, int num_redists, MPI_Comm comm) {
170  int result;
171 
172  for (int i = 0; i < num_redists; ++i) {
173 
174  if (redists[i] == NULL)
175  Xt_abort(comm, "ERROR: invalid redist; cannot build "
176  "redistribution collection\n", __FILE__, __LINE__);
177 
178  xt_mpi_call(MPI_Comm_compare(xt_redist_get_MPI_Comm(redists[i]),
179  comm, &result), comm);
180 
181  if ((result != MPI_IDENT) && (result != MPI_CONGRUENT))
182  Xt_abort(comm, "ERROR: MPI communicators do not match; cannot build "
183  "redistribution collection\n", __FILE__, __LINE__);
184  }
185 }
186 
187 static size_t
188 xt_ranks_uniq_count(size_t num_rank_sets,
189  const size_t *restrict num_ranks,
190  const int *const ranks[num_rank_sets])
191 {
192  size_t rank_pos[num_rank_sets];
193  for (size_t j = 0; j < num_rank_sets; ++j)
194  rank_pos[j] = 0;
195  bool ranks_left;
196  size_t num_messages = 0;
197  do {
198  int min_rank = INT_MAX;
199  /* find minimal rank in list, guaranteed to be smaller than comm_size */
200  for (size_t j = 0; j < num_rank_sets; ++j)
201  if (rank_pos[j] < num_ranks[j] && ranks[j][rank_pos[j]] < min_rank)
202  min_rank = ranks[j][rank_pos[j]];
203  ranks_left = false;
204  /* increment list index for all redists matching minimal rank and
205  * see if any ranks are left */
206  for (size_t j = 0; j < num_rank_sets; ++j) {
207  rank_pos[j]
208  += (rank_pos[j] < num_ranks[j] && ranks[j][rank_pos[j]] == min_rank);
209  ranks_left |= (rank_pos[j] < num_ranks[j]);
210  }
211  ++num_messages;
212  } while (ranks_left);
213  return num_messages;
214 }
215 
216 unsigned
217 xt_redist_agg_msg_count(size_t num_redists, enum xt_msg_direction direction,
218  const Xt_redist redists[num_redists],
219  size_t num_ranks[num_redists],
220  int *restrict ranks[num_redists])
221 {
222  bool ranks_left = false;
223  /* get lists of ranks to send/receive message to/from */
224  for (size_t j = 0; j < num_redists; ++j) {
225  size_t nranks = num_ranks[j]
226  = (size_t)xt_redist_get_msg_ranks(redists[j], direction, ranks + j);
227  ranks_left |= (nranks > 0);
228  /* sort list */
229  xt_sort_int(ranks[j], nranks);
230  }
231  /* count number of different ranks to send/receive message to/from */
232  size_t num_messages = ranks_left
233  ? xt_ranks_uniq_count(num_redists, num_ranks, (const int *const *)ranks)
234  : 0;
235  return (unsigned)num_messages;
236 }
237 
238 MPI_Datatype
240  const MPI_Aint displacements[count],
241  const MPI_Datatype datatypes[count],
242  const int block_lengths[count],
243  MPI_Comm comm)
244 {
245  size_t num_datatypes = 0;
246  /* allocate more than max_auto_dt datatype items from heap */
247  enum { max_auto_dt = 8 };
248  for (size_t i = 0; i < count; ++i)
249  num_datatypes += (datatypes[i] != MPI_DATATYPE_NULL);
250  MPI_Datatype *datatypes_, dt_auto[max_auto_dt];
251  MPI_Aint *displacements_, disp_auto[max_auto_dt];
252  int *block_lengths_, bl_auto[max_auto_dt];
253 
254  if (num_datatypes != count) {
255  if (num_datatypes > max_auto_dt) {
256  datatypes_ = xmalloc(num_datatypes * sizeof(*datatypes_));
257  displacements_ = xmalloc(num_datatypes * sizeof(*displacements_));
258  block_lengths_ = xmalloc(num_datatypes * sizeof(*block_lengths_));
259  } else {
260  datatypes_ = dt_auto;
261  displacements_ = disp_auto;
262  block_lengths_ = bl_auto;
263  }
264  num_datatypes = 0;
265 
266  for (size_t i = 0; i < count; ++i) {
267  if (datatypes[i] != MPI_DATATYPE_NULL) {
268 
269  datatypes_[num_datatypes] = datatypes[i];
270  displacements_[num_datatypes] = displacements[i];
271  block_lengths_[num_datatypes] = block_lengths[i];
272  ++num_datatypes;
273  }
274  }
275  } else {
276  datatypes_ = (MPI_Datatype *)datatypes;
277  displacements_ = (MPI_Aint *)displacements;
278  block_lengths_ = (int *)block_lengths;
279  }
280  MPI_Datatype datatype;
281  if (num_datatypes > 1)
282  xt_mpi_call(MPI_Type_create_struct((int)num_datatypes, block_lengths_,
283  displacements_, datatypes_, &datatype),
284  comm);
285  else if (displacements_[0] == 0)
286  xt_mpi_call(MPI_Type_dup(datatypes_[0], &datatype), comm);
287  else
288  xt_mpi_call(MPI_Type_create_hindexed(1, (int [1]){1}, displacements_,
289  datatypes_[0], &datatype), comm);
290 
291  xt_mpi_call(MPI_Type_commit(&datatype), comm);
292 
293  if (num_datatypes != count && num_datatypes > max_auto_dt) {
294  free(datatypes_);
295  free(displacements_);
296  }
297 
298  return datatype;
299 }
300 
301 /*
302  * Local Variables:
303  * c-basic-offset: 2
304  * coding: utf-8
305  * indent-tabs-mode: nil
306  * show-trailing-whitespace: t
307  * require-trailing-newline: t
308  * End:
309  */
int MPI_Comm
Definition: core.h:64
add versions of standard API functions not returning on error
#define xmalloc(size)
Definition: ppm_xfuncs.h:70
const struct xt_redist_vtable * vtable
MPI_Datatype datatype
Definition: xt_redist.h:69
Xt_redist(* copy)(Xt_redist)
void(* a_exchange)(Xt_redist, int, const void **, void **, Xt_request *)
void(* s_exchange1)(Xt_redist, const void *, void *)
MPI_Datatype(* get_msg_MPI_Datatype)(Xt_redist, int, enum xt_msg_direction)
void(* a_exchange1)(Xt_redist, const void *, void *, Xt_request *)
int(* get_num_msg)(Xt_redist, enum xt_msg_direction)
int(* get_msg_ranks)(Xt_redist, enum xt_msg_direction, int *restrict *)
void(* s_exchange)(Xt_redist, int, const void **, void **)
void(* delete)(Xt_redist)
MPI_Comm(* get_MPI_Comm)(Xt_redist)
base definitions header file
utility routines for MPI
#define xt_mpi_call(call, comm)
Definition: xt_mpi.h:68
MPI_Datatype xt_redist_get_MPI_Datatype(Xt_redist redist, int rank, enum xt_msg_direction direction)
Definition: xt_redist.c:117
int xt_redist_get_msg_ranks(Xt_redist redist, enum xt_msg_direction direction, int *restrict *ranks)
Definition: xt_redist.c:128
MPI_Datatype xt_redist_get_recv_MPI_Datatype(Xt_redist redist, int rank)
Definition: xt_redist.c:112
void xt_redist_delete(Xt_redist redist)
Definition: xt_redist.c:68
void xt_redist_msgs_strided_copy(size_t n, const struct Xt_redist_msg *restrict src, size_t src_stride, struct Xt_redist_msg *restrict dst, size_t dst_stride, MPI_Comm comm)
Definition: xt_redist.c:136
int xt_redist_get_num_recv_msg(Xt_redist redist)
Definition: xt_redist.c:102
static size_t xt_ranks_uniq_count(size_t num_rank_sets, const size_t *restrict num_ranks, const int *const ranks[num_rank_sets])
Definition: xt_redist.c:188
unsigned xt_redist_agg_msg_count(size_t num_redists, enum xt_msg_direction direction, const Xt_redist redists[num_redists], size_t num_ranks[num_redists], int *restrict ranks[num_redists])
Definition: xt_redist.c:217
void xt_redist_check_comms(Xt_redist *redists, int num_redists, MPI_Comm comm)
Definition: xt_redist.c:169
int xt_redist_get_num_send_msg(Xt_redist redist)
Definition: xt_redist.c:97
MPI_Datatype xt_create_compound_datatype(size_t count, const MPI_Aint displacements[count], const MPI_Datatype datatypes[count], const int block_lengths[count], MPI_Comm comm)
Definition: xt_redist.c:239
void xt_redist_a_exchange1(Xt_redist redist, const void *src_data, void *dst_data, Xt_request *request)
Definition: xt_redist.c:91
Xt_redist xt_redist_copy(Xt_redist redist)
Definition: xt_redist.c:63
void xt_redist_msgs_strided_destruct(size_t n, struct Xt_redist_msg *msgs, MPI_Comm comm, size_t ofs_stride)
Definition: xt_redist.c:156
MPI_Comm xt_redist_get_MPI_Comm(Xt_redist redist)
Definition: xt_redist.c:123
MPI_Datatype xt_redist_get_send_MPI_Datatype(Xt_redist redist, int rank)
Definition: xt_redist.c:107
void xt_redist_a_exchange(Xt_redist redist, int num_arrays, const void **src_data, void **dst_data, Xt_request *request)
Definition: xt_redist.c:79
void xt_redist_s_exchange(Xt_redist redist, int num_arrays, const void **src_data, void **dst_data)
Definition: xt_redist.c:73
void xt_redist_s_exchange1(Xt_redist redist, const void *src_data, void *dst_data)
Definition: xt_redist.c:86
redistribution of data
redistribution of data, non-public declarations
xt_msg_direction
@ RECV
@ SEND
void(* xt_sort_int)(int *a, size_t n)
Definition: xt_sort.c:53