Yet Another eXchange Tool  0.9.0
test_redist_p2p_parallel_f.f90
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 PROGRAM test_redist_p2p_parallel
47  USE ftest_common, ONLY: init_mpi, finish_mpi, test_abort
48  USE mpi
49  USE yaxt, ONLY: xt_initialize, xt_finalize, &
50  xt_int_kind, xi => xt_int_kind, &
51  xt_idxlist, xt_idxvec_new, xt_idxlist_delete, &
53  xt_redist, xt_redist_p2p_new, xt_redist_get_mpi_comm, &
56  xt_config, xt_config_delete
57  USE test_idxlist_utils, ONLY: test_err_count
58  USE test_redist_common, ONLY: communicators_are_congruent, &
59  check_redist, redist_exchanger_option
60  IMPLICIT NONE
61 
62  CHARACTER(len=*), PARAMETER :: filename = 'test_redist_p2p_parallel_f.f90'
63  TYPE(xt_config) :: config
64  INTEGER :: comm_rank, comm_size, ierror
65 
66  CALL init_mpi
67  CALL xt_initialize(mpi_comm_world)
68  config = redist_exchanger_option()
69 
70  CALL mpi_comm_rank(mpi_comm_world, comm_rank, ierror)
71  IF (ierror /= mpi_success) &
72  CALL test_abort("MPI error!", filename, __line__)
73 
74  CALL mpi_comm_size(mpi_comm_world, comm_size, ierror)
75  IF (ierror /= mpi_success) &
76  CALL test_abort("MPI error!", filename, __line__)
77 
78  CALL simple_test
79  CALL nonuniform_test
80  CALL block_redist_test
81 
82  IF (test_err_count() /= 0) &
83  CALL test_abort("non-zero error count!", filename, __line__)
84  CALL xt_config_delete(config)
85  CALL xt_finalize
86  CALL finish_mpi
87 
88 CONTAINS
89  SUBROUTINE simple_test
90  INTEGER, PARAMETER :: data_size = 10
91  INTEGER, PARAMETER :: src_num_indices = data_size, &
92  dst_num_indices = data_size
93  INTEGER(xt_int_kind) :: src_index_list(data_size), &
94  dst_index_list(data_size)
95  TYPE(xt_idxlist) :: src_idxlist, dst_idxlist
96  TYPE(xt_xmap) :: xmap
97  TYPE(xt_redist) :: redist
98  DOUBLE PRECISION :: src_data(data_size), dst_data(data_size)
99  INTEGER :: i
100 
101  ! source index list
102  DO i = 1, src_num_indices
103  src_index_list(i) = int(comm_rank * data_size + (i - 1), xi)
104  END DO
105 
106  src_idxlist = xt_idxvec_new(src_index_list)
107  ! destination index list
108  DO i = 1, dst_num_indices
109  dst_index_list(i) &
110  = int(mod(comm_rank * data_size + i + 1, comm_size * data_size), xi)
111  END DO
112  dst_idxlist = xt_idxvec_new(dst_index_list)
113  ! xmap
114  xmap = xt_xmap_all2all_new(src_idxlist, dst_idxlist, mpi_comm_world)
115  ! redist_p2p
116  redist = xt_redist_p2p_new(xmap, mpi_double_precision, config)
117 
118  ! test communicator of redist
119  IF (.NOT. communicators_are_congruent(xt_redist_get_mpi_comm(redist), &
120  mpi_comm_world)) &
121  CALL test_abort("error in xt_redist_get_mpi_comm", filename, __line__)
122 
123  DO i = 1, src_num_indices
124  src_data(i) = dble(comm_rank * data_size + i - 1)
125  END DO
126 
127  CALL check_redist(redist, src_data, dst_data, dst_index_list)
128 
129  ! clean up
130  CALL xt_redist_delete(redist)
131  CALL xt_xmap_delete(xmap)
132  CALL xt_idxlist_delete(src_idxlist)
133  CALL xt_idxlist_delete(dst_idxlist)
134  END SUBROUTINE simple_test
135 
136  ! test nonuniform numbers of send and receive partners
137  SUBROUTINE nonuniform_test
138  ! source index list
139  INTEGER(xt_int_kind), ALLOCATABLE :: src_index_list(:), dst_index_list(:)
140  DOUBLE PRECISION, ALLOCATABLE :: src_data(:), dst_data(:), ref_dst_data(:)
141  TYPE(xt_idxlist) :: src_idxlist, dst_idxlist
142  TYPE(xt_xmap) :: xmap
143  TYPE(Xt_redist) :: redist
144  INTEGER :: i, src_num_indices, dst_num_indices
145 
146  ALLOCATE(src_index_list(comm_size), dst_index_list(comm_size), &
147  src_data(comm_size), dst_data(comm_size), ref_dst_data(comm_size))
148  src_num_indices = merge(comm_size, 0, comm_rank == 0)
149  DO i = 1, src_num_indices
150  src_index_list(i) = int(i - 1, xi)
151  END DO
152 
153  src_idxlist = xt_idxvec_new(src_index_list, src_num_indices)
154 
155  ! destination index list
156  dst_num_indices = comm_size
157  DO i = 1, dst_num_indices
158  dst_index_list(i) = int(i - 1, xi)
159  END DO
160 
161  dst_idxlist = xt_idxvec_new(dst_index_list, dst_num_indices)
162 
163  ! xmap
164  xmap = xt_xmap_all2all_new(src_idxlist, dst_idxlist, mpi_comm_world)
165 
166  ! redist_p2p
167  redist = xt_redist_p2p_new(xmap, mpi_double_precision, config)
168 
169  ! test communicator of redist
170  IF (.NOT. communicators_are_congruent(xt_redist_get_mpi_comm(redist), &
171  mpi_comm_world)) &
172  CALL test_abort("error in xt_redist_get_mpi_comm", filename, __line__)
173 
174  ! test exchange
175  IF (comm_rank == 0) THEN
176  DO i = 1, comm_size
177  src_data(i) = dble(i - 1)
178  END DO
179  ELSE
180  src_data(:) = -2.0d0
181  END IF
182 
183  DO i = 1, comm_size
184  ref_dst_data(i) = dble(i-1)
185  END DO
186  CALL check_redist(redist, src_data, dst_data, ref_dst_data)
187 
188  ! clean up
189  CALL xt_redist_delete(redist)
190  CALL xt_xmap_delete(xmap)
191  CALL xt_idxlist_delete(src_idxlist)
192  CALL xt_idxlist_delete(dst_idxlist)
193  END SUBROUTINE nonuniform_test
194 
195  ! test redist with blocks
196  SUBROUTINE block_redist_test
197  ! gvol_size: volume of deep ocean
198  INTEGER :: ngdom, gvol_size, i, nwin, ig0, ig, j, p, qa, qb
199  ! gdepth: ocean depth of an one dim. ocean
200  INTEGER, ALLOCATABLE :: gdoma(:), gdomb(:), gsurfdata(:), &
201  gdepth(:), ig2col_off(:), b_surfdata_ref(:), gvoldata(:), &
202  src_block_offsets(:), src_block_sizes(:), dst_block_offsets(:), &
203  dst_block_sizes(:), b_voldata_ref(:)
204  INTEGER, ALLOCATABLE :: a_surfdata(:), b_surfdata(:), &
205  a_voldata(:), b_voldata(:)
206  INTEGER(xi), ALLOCATABLE :: iveca(:), ivecb(:)
207  INTEGER :: ia, ib, blk_ofs_accum, gdepth_i
208  TYPE(Xt_idxlist) :: idxlist_a, idxlist_b
209  TYPE(xt_xmap) :: xmap
210  TYPE(Xt_redist) :: redist, block_redist, block_redist2
211 
212  IF (2 * comm_size > huge(1_xt_int_kind)) &
213  CALL test_abort('too large number of tasks', filename, __line__)
214  ! the global index domain (1dim problem):
215  ngdom = 2 * comm_size
216  ! start state (index distribution) of global domain
217  ALLOCATE(gdoma(ngdom), gdomb(ngdom))
218  ! end state ""
219  ALLOCATE(gsurfdata(ngdom), gdepth(ngdom))
220  ALLOCATE(ig2col_off(ngdom)) ! offset of surface DATA within vol
221  gvol_size = 0
222  DO i = 1, ngdom
223  gdoma(i) = i - 1
224  gdomb(i) = ngdom - i
225  gsurfdata(i) = 99 + i
226  gdepth(i) = i
227  ig2col_off(i) = gvol_size
228  gvol_size = gvol_size + gdepth(i)
229  END DO
230 
231  nwin = ngdom / comm_size ! my local window size of the global surface domain
232  ! start of my window within global index domain (== global offset)
233  ig0 = comm_rank * nwin
234  IF (nwin * comm_size /= ngdom) &
235  CALL test_abort("internal error", filename, __line__)
236 
237  ! local index
238  ALLOCATE(iveca(nwin), ivecb(nwin))
239  DO i = 1, nwin
240  ig = ig0 + i
241  iveca(i) = int(gdoma(ig), xi)
242  ivecb(i) = int(gdomb(ig), xi)
243  END DO
244 
245  idxlist_a = xt_idxvec_new(iveca, nwin)
246  idxlist_b = xt_idxvec_new(ivecb, nwin)
247 
248  xmap = xt_xmap_all2all_new(idxlist_a, idxlist_b, mpi_comm_world)
249 
250  ! simple redist
251  redist = xt_redist_p2p_new(xmap, mpi_integer, config)
252 
253  ! test communicator of redist
254  IF (.NOT. communicators_are_congruent(xt_redist_get_mpi_comm(redist), &
255  mpi_comm_world)) &
256  CALL test_abort("error in xt_redist_get_mpi_comm", filename, __line__)
257 
258  ALLOCATE(a_surfdata(nwin), b_surfdata(nwin), b_surfdata_ref(nwin))
259  DO i = 1, nwin
260  a_surfdata(i) = gsurfdata(iveca(i) + 1)
261  b_surfdata(i) = -1
262  b_surfdata_ref(i) = gsurfdata(ivecb(i) + 1)
263  END DO
264 
265  CALL check_redist(redist, a_surfdata, b_surfdata, b_surfdata_ref)
266  CALL xt_redist_delete(redist)
267 
268  ! generate global volume data
269  ALLOCATE(gvoldata(gvol_size))
270  DO i = 1, ngdom
271  DO j = 1, gdepth(i)
272  p = ig2col_off(i) + j
273  gvoldata(p) = (i - 1) * 100 + j - 1
274  END DO
275  END DO
276 
277  ! generate blocks
278  ALLOCATE(src_block_offsets(nwin), src_block_sizes(nwin), &
279  dst_block_offsets(nwin), dst_block_sizes(nwin))
280  ! we only need local size but simply oversize here
281  ALLOCATE(a_voldata(gvol_size), b_voldata(gvol_size), &
282  b_voldata_ref(gvol_size))
283  a_voldata(:) = -1
284  b_voldata_ref(:) = -1
285 
286  qa = 0
287  blk_ofs_accum = 0
288  DO i = 1, nwin
289  ia = int(iveca(i)) + 1
290  gdepth_i = gdepth(ia)
291  src_block_offsets(i) = blk_ofs_accum
292  blk_ofs_accum = blk_ofs_accum + gdepth_i
293  src_block_sizes(i) = gdepth_i
294  p = ig2col_off(ia)
295  DO j = 1, gdepth_i
296  a_voldata(qa + j) = gvoldata(p + j)
297  END DO
298  qa = qa + gdepth_i
299  END DO
300 
301  qb = 0
302  blk_ofs_accum = 0
303  DO i = 1, nwin
304  ib = int(ivecb(i)) + 1
305  gdepth_i = gdepth(ib)
306  dst_block_offsets(i) = blk_ofs_accum
307  blk_ofs_accum = blk_ofs_accum + gdepth_i
308  dst_block_sizes(i) = gdepth_i
309  p = ig2col_off(ib)
310  DO j = 1, gdepth_i
311  b_voldata_ref(qb + j) = gvoldata(p + j)
312  END DO
313  qb = qb + gdepth_i
314  END DO
315 
316  ! redist with blocks
317  block_redist = xt_redist_p2p_blocks_off_custom_new(xmap, &
318  src_block_offsets, src_block_sizes, nwin, &
319  dst_block_offsets, dst_block_sizes, nwin, mpi_integer, config)
320  ! test communicator of redist
321  IF (.NOT. communicators_are_congruent(xt_redist_get_mpi_comm(block_redist), &
322  mpi_comm_world)) &
323  CALL test_abort("error in xt_redist_get_mpi_comm", filename, __line__)
324 
325  CALL check_redist(block_redist, a_voldata, b_voldata, b_voldata_ref)
326 
327  ! redist with blocks but without explicit offsets:
328  block_redist2 = xt_redist_p2p_blocks_custom_new(xmap, &
329  src_block_sizes, nwin, dst_block_sizes, nwin, mpi_integer, config)
330  ! test communicator of redist
331 
332  IF (.NOT. communicators_are_congruent(xt_redist_get_mpi_comm(block_redist2),&
333  mpi_comm_world)) &
334  CALL test_abort("error in xt_redist_get_mpi_comm", filename, __line__)
335 
336  CALL check_redist(block_redist2, a_voldata, b_voldata, b_voldata_ref)
337 
338  ! cleanup
339  CALL xt_redist_delete(block_redist2)
340  CALL xt_redist_delete(block_redist)
341  CALL xt_xmap_delete(xmap)
342  CALL xt_idxlist_delete(idxlist_a)
343  CALL xt_idxlist_delete(idxlist_b)
344  END SUBROUTINE block_redist_test
345 
346 END PROGRAM test_redist_p2p_parallel
347 !
348 ! Local Variables:
349 ! f90-continuation-indent: 5
350 ! coding: utf-8
351 ! indent-tabs-mode: nil
352 ! show-trailing-whitespace: t
353 ! require-trailing-newline: t
354 ! End:
355 !
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_delete(Xt_idxlist idxlist)
Definition: xt_idxlist.c:74
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_p2p_blocks_custom_new(Xt_xmap xmap, const int *src_block_sizes, int src_block_num, const int *dst_block_sizes, int dst_block_num, MPI_Datatype datatype, Xt_config config)
Xt_redist xt_redist_p2p_new(Xt_xmap xmap, MPI_Datatype datatype)
Xt_redist xt_redist_p2p_blocks_off_custom_new(Xt_xmap xmap, const int *src_block_offsets, const int *src_block_sizes, int src_block_num, const int *dst_block_offsets, const int *dst_block_sizes, int dst_block_num, MPI_Datatype datatype, Xt_config config)
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)