Yet Another eXchange Tool  0.9.0
test_redist_collection_parallel_f.f90
1 
13 
14 !
15 ! Keywords:
16 ! Maintainer: Jörg Behrens <behrens@dkrz.de>
17 ! Moritz Hanke <hanke@dkrz.de>
18 ! Thomas Jahns <jahns@dkrz.de>
19 ! URL: https://doc.redmine.dkrz.de/yaxt/html/
20 !
21 ! Redistribution and use in source and binary forms, with or without
22 ! modification, are permitted provided that the following conditions are
23 ! met:
24 !
25 ! Redistributions of source code must retain the above copyright notice,
26 ! this list of conditions and the following disclaimer.
27 !
28 ! Redistributions in binary form must reproduce the above copyright
29 ! notice, this list of conditions and the following disclaimer in the
30 ! documentation and/or other materials provided with the distribution.
31 !
32 ! Neither the name of the DKRZ GmbH nor the names of its contributors
33 ! may be used to endorse or promote products derived from this software
34 ! without specific prior written permission.
35 !
36 ! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS
37 ! IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
38 ! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A
39 ! PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER
40 ! OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
41 ! EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
42 ! PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
43 ! PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
44 ! LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
45 ! NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
46 ! SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
47 !
48 PROGRAM test_redist_collection_parallel
49  USE mpi
50  USE ftest_common, ONLY: init_mpi, finish_mpi, test_abort, icbrt
51  USE test_idxlist_utils, ONLY: test_err_count
52  USE yaxt, ONLY: xt_initialize, xt_finalize, xt_int_kind, xi => xt_int_kind, &
53  xt_idxlist, xt_idxlist_delete, xt_stripe, xt_idxvec_new, &
58  xt_idxlist_get_indices, xt_int_mpidt, &
59  xt_request, xt_redist_a_exchange, xt_config, xt_config_delete
60  ! older PGI compilers do not handle generic interface correctly
61 #if defined __PGI && (__PGIC__ < 12 || (__PGIC__ == 12 && __PGIC_MINOR__ <= 10))
63 #endif
64  USE test_redist_common, ONLY: check_wait_request, redist_exchanger_option
65  USE iso_c_binding, ONLY: c_loc, c_ptr
66 #include "xt_slice_c_loc.inc"
67  IMPLICIT NONE
68  INTEGER :: rank, world_size, ierror
69  CHARACTER(len=*), PARAMETER :: &
70  filename = 'test_redist_collection_parallel_f.f90'
71  CHARACTER(len=*), PARAMETER :: err_msg(2) = &
72  (/ "error in xt_redist_s_exchange", "error in xt_redist_a_exchange" /)
73  TYPE(xt_config) :: config
74 
75  CALL init_mpi
76  CALL xt_initialize(mpi_comm_world)
77  config = redist_exchanger_option()
78 
79  CALL mpi_comm_rank(mpi_comm_world, rank, ierror)
80  IF (ierror /= mpi_success) &
81  CALL test_abort('mpi_comm_rank failed', filename, __line__)
82  CALL mpi_comm_size(mpi_comm_world, world_size, ierror)
83  IF (ierror /= mpi_success) &
84  CALL test_abort('mpi_comm_size failed', filename, __line__)
85 
86  IF (world_size > 1) THEN
87  CALL test_4redist(mpi_comm_world, config)
88  CALL test_rr_exchange(mpi_comm_world, config)
89  END IF
90 
91  IF (test_err_count() /= 0) &
92  CALL test_abort("non-zero error count!", filename, __line__)
93  CALL xt_finalize
94  CALL finish_mpi
95 CONTAINS
96  SUBROUTINE build_idxlists(indices_a, indices_b, indices_all)
97  TYPE(xt_idxlist), INTENT(out) :: indices_a, indices_b, indices_all
98 
99  TYPE(xt_idxlist) :: indices_a_(2)
100  INTEGER :: i
101  INTEGER(xt_int_kind), PARAMETER :: start = 0
102  INTEGER(xt_int_kind) :: global_size(2), local_start(2, 2)
103  INTEGER :: local_size(2)
104 
105  TYPE(xt_stripe) :: stripe
106 
107  global_size(1) = int(2 * world_size, xi)
108  global_size(2) = int(world_size**2, xi)
109  local_size = world_size
110  local_start = reshape((/ 0_xi, int(rank*world_size, xi), &
111  int(world_size, xi), &
112  int((world_size-(rank+1))*world_size, xi) /), (/ 2, 2 /))
113 
114  DO i = 1, 2
115  indices_a_(i) = xt_idxsection_new(start, global_size, local_size, &
116  local_start(:, i))
117  END DO
118  indices_a = xt_idxlist_collection_new(indices_a_)
119 
120  CALL xt_idxlist_delete(indices_a_(1))
121  CALL xt_idxlist_delete(indices_a_(2))
122 
123  stripe = xt_stripe(int(rank * 2 * world_size**2, xi), 1_xi, 2*world_size**2)
124  indices_b = xt_idxstripes_new(stripe)
125 
126  stripe = xt_stripe(0_xi, 1_xi, 2*world_size**3)
127  indices_all = xt_idxstripes_new(stripe)
128  END SUBROUTINE build_idxlists
129 
130  SUBROUTINE test_4redist(comm, config)
131  ! redist test with four different redists
132  INTEGER, INTENT(in) :: comm
133  TYPE(xt_config), INTENT(in) :: config
134  INTEGER, PARAMETER :: num_tx = 4
135  TYPE(xt_idxlist) :: indices_a, indices_b, indices_all
136  INTEGER(xt_int_kind), ALLOCATABLE :: index_vector_a(:), &
137  index_vector_b(:)
138  TYPE(xt_xmap) :: xmaps(num_tx)
139  TYPE(xt_redist) :: redists(num_tx), redist, redist_copy
140  INTEGER :: i, vec_size
141 
142  IF (world_size &
143  > icbrt((huge(1_xi)-mod(huge(1_xi),2_xi))/2_xi)) &
144  CALL test_abort('communicator too large for test', filename, __line__)
145 
146  vec_size = 2*world_size**2
147  ALLOCATE(index_vector_a(vec_size), index_vector_b(vec_size))
148  CALL build_idxlists(indices_a, indices_b, indices_all)
149 
150  xmaps(1) = xt_xmap_all2all_new(indices_a, indices_b, comm)
151  xmaps(2) = xt_xmap_all2all_new(indices_b, indices_a, comm)
152  xmaps(3) = xt_xmap_all2all_new(indices_a, indices_all, comm)
153  xmaps(4) = xt_xmap_all2all_new(indices_b, indices_all, comm)
154 
155  CALL xt_idxlist_get_indices(indices_a, index_vector_a)
156  CALL xt_idxlist_get_indices(indices_b, index_vector_b)
157 
158  CALL xt_idxlist_delete(indices_a)
159  CALL xt_idxlist_delete(indices_b)
160  CALL xt_idxlist_delete(indices_all)
161 
162  DO i = 1, num_tx
163  redists(i) = xt_redist_p2p_new(xmaps(i), xt_int_mpidt)
164  CALL xt_xmap_delete(xmaps(i))
165  END DO
166 
167  redist = xt_redist_collection_new(redists, num_tx, -1, comm, config)
168 
169  ! test communicator of redist
170  ! if (!test_communicator(xt_redist_get_MPI_Comm(redist), COMM))
171  ! PUT_ERR("error in xt_redist_get_MPI_Comm\n");
172 
173  CALL xt_redist_delete(redists)
174 
175  CALL exchange_4redist(redist, index_vector_a, index_vector_b)
176  redist_copy = xt_redist_copy(redist)
177  CALL xt_redist_delete(redist)
178  CALL exchange_4redist(redist_copy, index_vector_a, index_vector_b)
179 
180  ! clean up
181  CALL xt_redist_delete(redist_copy)
182  END SUBROUTINE test_4redist
183 
184  SUBROUTINE exchange_4redist(redist, index_vector_a, index_vector_b)
185  TYPE(xt_redist), INTENT(in) :: redist
186  INTEGER(xt_int_kind), INTENT(in) :: index_vector_a(2*world_size**2), &
187  index_vector_b(2*world_size**2)
188  INTEGER(xt_int_kind), TARGET, ALLOCATABLE :: buf(:)
189  INTEGER(xt_int_kind), POINTER :: results_1(:), &
190  results_2(:), results_3(:), results_4(:)
191  INTEGER :: result_sizes(4), buf_size, ofs
192  INTEGER, PARAMETER :: result_spacing(4) = (/ 2, 14, 5, 8 /)
193  INTEGER :: iexch
194 
195  result_sizes(1) = 2*world_size**2
196  result_sizes(2) = 2*world_size**2
197  result_sizes(3) = 2*world_size**3
198  result_sizes(4) = 2*world_size**3
199 
200  buf_size = sum(result_spacing) + sum(result_sizes)
201  ALLOCATE(buf(buf_size))
202  DO iexch = 1, 2
203  buf = -1_xt_int_kind
204  ofs = result_spacing(1)
205  results_1 => buf(ofs+1:ofs+result_sizes(1))
206  ofs = ofs + result_sizes(1) + result_spacing(2)
207  results_2 => buf(ofs+1:ofs+result_sizes(2))
208  ofs = ofs + result_sizes(2) + result_spacing(3)
209  results_3 => buf(ofs+1:ofs+result_sizes(3))
210  ofs = ofs + result_sizes(3) + result_spacing(4)
211  results_4 => buf(ofs+1:ofs+result_sizes(4))
212 
213  CALL do_4redist(redist, index_vector_a, index_vector_b, &
214  results_1, results_2, results_3, results_4, iexch)
215 
216  CALL check_4redist_results(results_1, results_2, results_3, results_4, &
217  index_vector_a, index_vector_b, iexch)
218  buf = -1_xt_int_kind
219  ! shift addresses around
220  IF (rank == 0) THEN
221  ofs = sum(result_spacing(1:2)) + sum(result_sizes(1:2))
222  results_3 => buf(ofs+1:ofs+result_sizes(3))
223  END IF
224 
225  CALL do_4redist(redist, index_vector_a, index_vector_b, &
226  results_1, results_2, results_3, results_4, iexch)
227 
228  CALL check_4redist_results(results_1, results_2, results_3, results_4, &
229  index_vector_a, index_vector_b, iexch)
230  ENDDO
231  ! clean up
232  DEALLOCATE(buf)
233  END SUBROUTINE exchange_4redist
234 
235  SUBROUTINE do_4redist(redist, index_vector_a, index_vector_b, &
236  results_1, results_2, results_3, results_4, iexch)
237  TYPE(xt_redist), INTENT(in) :: redist
238  INTEGER(xt_int_kind), INTENT(in), TARGET :: &
239  index_vector_a(*), index_vector_b(*)
240  INTEGER(xt_int_kind), INTENT(inout), TARGET :: &
241  results_1(*), results_2(*), results_3(*), results_4(*)
242  INTEGER, INTENT(in) :: iexch
243 
244  TYPE(c_ptr) :: results(4), input(4)
245  TYPE(xt_request) :: request
246 
247  results(1) = c_loc(results_1)
248  results(2) = c_loc(results_2)
249  results(3) = c_loc(results_3)
250  results(4) = c_loc(results_4)
251 
252  input(1) = c_loc(index_vector_a)
253  input(2) = c_loc(index_vector_b)
254  input(3) = c_loc(index_vector_a)
255  input(4) = c_loc(index_vector_b)
256  IF (iexch == 1) THEN
257  CALL xt_redist_s_exchange(redist, 4, input, results)
258  ELSE
259  CALL xt_redist_a_exchange(redist, 4, input, results, request)
260  CALL check_wait_request(request, filename, __line__)
261  ENDIF
262  END SUBROUTINE do_4redist
263 
264  SUBROUTINE check_4redist_results(results_1, results_2, results_3, results_4, &
265  index_vector_a, index_vector_b, iexch)
266  INTEGER(xt_int_kind), INTENT(in) :: index_vector_a(:), index_vector_b(:), &
267  results_1(:), results_2(:), results_3(0:), results_4(0:)
268  INTEGER, INTENT(in) :: iexch
269  INTEGER(xt_int_kind) :: i, n
270  LOGICAL :: p_3, p_4
271 
272  IF (any(results_1 /= index_vector_b)) &
273  CALL test_abort(err_msg(iexch), filename, __line__)
274 
275  IF (any(results_2 /= index_vector_a)) &
276  CALL test_abort(err_msg(iexch), filename, __line__)
277 
278  n = int(SIZE(results_3), xt_int_kind)
279  p_3 = .false.
280  p_4 = .false.
281  DO i = 0_xi, n - 1_xi
282  p_3 = p_3 .OR. results_3(i) /= i
283  p_4 = p_4 .OR. results_4(i) /= i
284  END DO
285  IF (p_3 .OR. p_4) CALL test_abort(err_msg(iexch), filename, __line__)
286  END SUBROUTINE check_4redist_results
287 
288 
289  ! redist test with two redists that do a round robin exchange in
290  ! different directions
291  SUBROUTINE test_rr_exchange(comm, config)
292  INTEGER, INTENT(in) :: comm
293  TYPE(xt_config), INTENT(in) :: config
294 
295  TYPE(xt_idxlist) :: src_indices, dst_indices(2)
296  INTEGER(xt_int_kind) :: src_indices_(5)
297  INTEGER(xt_int_kind) :: i, temp, dst_indices_(5, 2)
298  TYPE(xt_xmap) :: xmaps(2)
299  TYPE(xt_redist) :: redists(2), redist, redist_copy
300 
301  IF (world_size > (huge(1_xi)-mod(huge(1_xi),5_xi))/5_xi) &
302  CALL test_abort('communicator too large for test', filename, __line__)
303 
304  DO i = 1_xi, 5_xi
305  src_indices_(i) = int(rank, xi) * 5_xi + (i - 1_xi)
306  dst_indices_(i, 1) = mod(src_indices_(i) + 1_xi, &
307  & int(world_size, xi) * 5_xi)
308  temp = src_indices_(i) - 1_xi
309  dst_indices_(i, 2) = merge(int(world_size, xi) * 5_xi - 1_xi, &
310  & temp, temp < 0_xi)
311  END DO
312 
313  src_indices = xt_idxvec_new(src_indices_, 5)
314  dst_indices(1) = xt_idxvec_new(dst_indices_(:, 1))
315  dst_indices(2) = xt_idxvec_new(dst_indices_(:, 2))
316 
317  xmaps(1) = xt_xmap_all2all_new(src_indices, dst_indices(1), comm)
318  xmaps(2) = xt_xmap_all2all_new(src_indices, dst_indices(2), comm)
319 
320  CALL xt_idxlist_delete(src_indices)
321  CALL xt_idxlist_delete(dst_indices)
322 
323  redists(1) = xt_redist_p2p_new(xmaps(1), xt_int_mpidt)
324  redists(2) = xt_redist_p2p_new(xmaps(2), xt_int_mpidt)
325 
326  CALL xt_xmap_delete(xmaps)
327 
328  redist = xt_redist_collection_new(redists, 2, -1, comm, config)
329 
330  ! test communicator of redist
331  ! IF (!test_communicator(xt_redist_get_MPI_Comm(redist), comm))
332  ! PUT_ERR("error in xt_redist_get_MPI_Comm\n");
333 
334  CALL xt_redist_delete(redists)
335 
336  CALL rr_exchange(redist, src_indices_, dst_indices_)
337  redist_copy = xt_redist_copy(redist)
338  CALL xt_redist_delete(redist)
339  CALL rr_exchange(redist_copy, src_indices_, dst_indices_)
340 
341  ! clean up
342  CALL xt_redist_delete(redist_copy)
343  END SUBROUTINE test_rr_exchange
344 
345  SUBROUTINE rr_exchange(redist, src_indices_, ref_dst_indices_)
346 #if defined __GNUC__ && __GNUC__ >= 5 && __GNUC__ <= 8
347  ! gcc versions 5.x to 8.x have a bug that lets them evaluate the
348  ! ANY test too early if results never gets passed to some external
349  ! routine directly, 9.x is not only fixed again, but requires to
350  ! have the explicit escaping of a pointer to results via C_LOC
351  USE yaxt, ONLY: xt_slice_c_loc
352 #undef XT_SLICE_C_LOC
353 #define XT_SLICE_C_LOC(slice, cptr) CALL xt_slice_c_loc(slice, cptr)
354 #endif
355  TYPE(xt_redist), INTENT(in) :: redist
356  INTEGER, PARAMETER :: nredist = 2
357  INTEGER(xt_int_kind), TARGET, INTENT(in) :: src_indices_(5)
358  INTEGER(xt_int_kind), INTENT(in) :: ref_dst_indices_(5, nredist)
359 
360  INTEGER(xt_int_kind), TARGET :: results(5,nredist)
361  TYPE(c_ptr) :: results_p(nredist), input(nredist)
362  INTEGER :: iexch, i
363  TYPE(xt_request) :: request
364 
365  DO i = 1, nredist
366  xt_slice_c_loc(results(:,i), results_p(i))
367  input(i) = c_loc(src_indices_)
368  END DO
369 
370  DO iexch = 1, 2
371  results = -1
372 
373  IF (iexch == 1) THEN
374  CALL xt_redist_s_exchange(redist, input, results_p)
375  ELSE
376  CALL xt_redist_a_exchange(redist, input, results_p, request)
377  CALL check_wait_request(request, filename, __line__)
378  ENDIF
379 
380  ! check results
381  IF (any(results /= ref_dst_indices_)) &
382  CALL test_abort(err_msg(iexch), filename, __line__)
383  ENDDO
384  END SUBROUTINE rr_exchange
385 
386 END PROGRAM test_redist_collection_parallel
387 !
388 ! Local Variables:
389 ! f90-continuation-indent: 5
390 ! coding: utf-8
391 ! indent-tabs-mode: nil
392 ! show-trailing-whitespace: t
393 ! require-trailing-newline: t
394 ! End:
395 !
Definition: yaxt.f90:49
void xt_config_delete(Xt_config config)
Definition: xt_config.c:76
void xt_initialize(MPI_Comm default_comm)
Definition: xt_init.c:70
void xt_finalize(void)
Definition: xt_init.c:89
void xt_idxlist_get_indices(Xt_idxlist idxlist, Xt_int *indices)
Definition: xt_idxlist.c:102
void xt_idxlist_delete(Xt_idxlist idxlist)
Definition: xt_idxlist.c:74
Xt_idxlist xt_idxlist_collection_new(Xt_idxlist *idxlists, int num_idxlists)
Xt_idxlist xt_idxsection_new(Xt_int start, int num_dimensions, const Xt_int global_size[num_dimensions], const int local_size[num_dimensions], const Xt_int local_start[num_dimensions])
Xt_idxlist xt_idxstripes_new(struct Xt_stripe const *stripes, int num_stripes)
Xt_idxlist xt_idxvec_new(const Xt_int *idxlist, int num_indices)
Definition: xt_idxvec.c:163
void xt_redist_delete(Xt_redist redist)
Definition: xt_redist.c:68
Xt_redist xt_redist_copy(Xt_redist redist)
Definition: xt_redist.c:63
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
Xt_redist xt_redist_collection_new(Xt_redist *redists, int num_redists, int cache_size, MPI_Comm comm)
Xt_redist xt_redist_p2p_new(Xt_xmap xmap, MPI_Datatype datatype)
void xt_xmap_delete(Xt_xmap xmap)
Definition: xt_xmap.c:85
Xt_xmap xt_xmap_all2all_new(Xt_idxlist src_idxlist, Xt_idxlist dst_idxlist, MPI_Comm comm)