1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47#include "fc_feature_defs.inc"
48PROGRAM test_yaxt
49 USE iso_c_binding, ONLY: c_int
50 USE mpi
58 USE ftest_common, ONLY: test_abort, id_map, factorize, regular_deco, &
59 init_mpi, finish_mpi, cmp_arrays
60
61
62
63#if defined __PGI
67#endif
68 IMPLICIT NONE
69
70 INTEGER, PARAMETER :: g_ie = 8, g_je = 4
71 LOGICAL, PARAMETER :: verbose = .false.
72 INTEGER, PARAMETER :: nlev = 3
73 INTEGER, PARAMETER :: undef_int = -1
74 INTEGER(xt_int_kind), PARAMETER :: undef_index = -1
75 INTEGER, PARAMETER :: nhalo = 1
76 LOGICAL, PARAMETER :: increased_north_halo = .false.
77 LOGICAL, PARAMETER :: with_north_halo = .true.
78
79 INTEGER :: ie, je
80 INTEGER :: p_ioff, p_joff
81 INTEGER :: nprocx, nprocy
82 INTEGER :: nprocs
83 INTEGER :: mype, mypx, mypy
84 LOGICAL :: lroot
85
86 INTEGER(xt_int_kind) :: g_id(g_ie, g_je)
87 CHARACTER(len=*), PARAMETER :: filename = 'test_yaxt.f90'
88
89 INTEGER(xt_int_kind) :: g_tpex(g_ie, g_je)
90 TYPE(xt_xmap) :: xmap_tpex
91 TYPE(xt_redist) :: redist_tpex
92 TYPE(xt_redist) :: redist_surf_tpex
93
94 INTEGER(xt_int_kind), ALLOCATABLE :: loc_id(:,:), loc_tpex(:,:)
95 INTEGER, ALLOCATABLE :: fval(:,:), gval(:,:)
96 INTEGER, ALLOCATABLE :: gval3d(:,:,:)
97 INTEGER, ALLOCATABLE :: id_pos(:,:), pos3d_surf(:,:)
98
99
100 CALL init_all
101
102
103 CALL id_map(g_id)
104
105
106 CALL get_window(g_id, loc_id)
107
108
109 CALL def_exchange(g_id, g_tpex)
110
111
112 CALL get_window(g_tpex, loc_tpex)
113
114
115 CALL general_fsection_test
116
117
118 CALL check_modifiers
119
120
121 CALL gen_template(loc_id, loc_tpex, xmap_tpex)
122
123
124 CALL gen_trans(xmap_tpex, mpi_integer, mpi_integer, redist_tpex)
125
126
127 fval = int(loc_id)
129
130 IF (cmp_arrays(gval, loc_tpex)) &
131 CALL test_abort('array eqivalence test failed', filename, __line__)
132 DEALLOCATE(gval)
133 CALL check_redist_collection_static
134 DEALLOCATE(loc_id)
135
136
137 CALL gen_id_pos(id_pos)
138 CALL gen_id_pos(pos3d_surf)
139 CALL gen_pos3d_surf(pos3d_surf)
140
141
142 CALL gen_off_trans(xmap_tpex, mpi_integer, id_pos(:,:) - 1, &
143 mpi_integer, pos3d_surf, redist_surf_tpex)
144 DEALLOCATE(pos3d_surf, id_pos)
145
146
147 gval3d = -1
149 reshape(fval, (/ ie, 1, je /)), gval3d)
150 DEALLOCATE(fval)
151
152 IF (cmp_arrays(gval3d(:, 1, :), loc_tpex)) &
153 CALL test_abort('surface check failed', filename, __line__)
154
155 IF (any(gval3d(:, 2, :) /= -1)) &
156 CALL test_abort('surface check failed', filename, __line__)
157
158
159 DEALLOCATE(loc_tpex, gval3d)
161
163
165
167 CALL finish_mpi
168
169CONTAINS
170
171#define abort(msg, line) test_abort(msg, filename, line)
172
173 SUBROUTINE check_redist_collection_static
174 INTEGER, PARAMETER :: nr = 2
175 TYPE(xt_redist) :: rvec(nr), rcol
176 INTEGER :: f(ie,je,nr), g(ie,je,nr), ref_g(ie,je,nr)
177 INTEGER(mpi_address_kind) :: f_addr(nr), g_addr(nr)
178 INTEGER(mpi_address_kind) :: f_disp(nr), g_disp(nr)
179
180 INTEGER :: ir, ierror
181 rvec(:) = redist_tpex
182 DO ir = 1, nr
183 CALL mpi_get_address(f(1,1,ir), f_addr(ir), ierror)
184 IF (ierror /= mpi_success) &
185 CALL abort('MPI_GET_ADDRESS failed', __line__)
186 CALL mpi_get_address(g(1,1,ir), g_addr(ir), ierror)
187 IF (ierror /= mpi_success) &
188 CALL abort('MPI_GET_ADDRESS failed', __line__)
189 f_disp(ir) = f_addr(ir) - f_addr(1)
190 g_disp(ir) = g_addr(ir) - g_addr(1)
191 ENDDO
192
194 DO ir = 1, nr
195 f(:,:,ir) = int(loc_id) + (ir-1) * ie*je
196 ENDDO
197
198 ref_g = 0
199 DO ir = 1, nr
201 ENDDO
202
203 g = 0
205 IF (any(g /= ref_g)) CALL abort('(g /= ref_g)', __line__)
207
208 END SUBROUTINE check_redist_collection_static
209
210 SUBROUTINE check_modifiers()
211 TYPE(xt_modifier) :: m_tpex(5)
212 INTEGER :: m_tpex_num
213 TYPE(xt_idxlist) :: loc_id_idxlist
214 INTEGER(xt_int_kind) :: loc_tpex2(ie,je)
215 TYPE(xt_idxlist) :: loc_tpex2_idxlist
216 INTEGER(c_int), ALLOCATABLE :: mstate(:,:)
217
219 ALLOCATE(mstate(ie,je))
220
221
222 CALL def_tpex_mod_via_idxvec(m_tpex, m_tpex_num)
223
224 loc_tpex2_idxlist =
xt_idxmod_new(loc_id_idxlist, m_tpex, m_tpex_num, mstate)
225 loc_tpex2 = -1
227 IF (any(loc_tpex2 /= loc_tpex)) &
228 CALL abort('idx copy does not match', __line__)
230
231
232 loc_tpex2_idxlist =
xt_idxmod_new(loc_id_idxlist, m_tpex, m_tpex_num)
233 loc_tpex2 = -1
235 IF (any(loc_tpex2 /= loc_tpex)) &
236 CALL abort('idx copy does not match', __line__)
238 CALL delete_modifiers(m_tpex(1:m_tpex_num))
239
240
241 CALL def_tpex_mod_via_sections(m_tpex, m_tpex_num)
242 loc_tpex2_idxlist =
xt_idxmod_new(loc_id_idxlist, m_tpex, m_tpex_num, mstate)
243 loc_tpex2 = -1
245
246 IF (any(loc_tpex2 /= loc_tpex)) &
247 CALL abort('idx copy does not match', __line__)
249 CALL delete_modifiers(m_tpex(1:m_tpex_num))
250
251
253 END SUBROUTINE check_modifiers
254
255 SUBROUTINE delete_modifiers(m)
256 TYPE(xt_modifier), INTENT(inout) :: m(:)
257
258 INTEGER :: i
259
260 DO i = 1, SIZE(m)
263 ENDDO
264
265 END SUBROUTINE delete_modifiers
266
267 SUBROUTINE general_fsection_test
268 INTEGER(xt_int_kind), PARAMETER :: gdx = 10_xt_int_kind, gdy=5_xt_int_kind
269 INTEGER, PARAMETER :: ldx = 4, ldy=2
270 INTEGER(xt_int_kind), PARAMETER :: gstart = 1
271 INTEGER(xt_int_kind), PARAMETER :: gsize(2) = (/ gdx, gdy /)
272 TYPE(xt_idxlist) :: global_section, local_section
273 INTEGER(xt_int_kind) :: indices(gdx*gdy), lstart(2)
274 INTEGER :: egis(gdx, gdy)
275 INTEGER :: i, j, idx, p
276
277
278 idx = gstart - 1
279 DO j = 1, gdy
280 DO i = 1, gdx
281 idx = idx + 1
282 egis(i,j) = idx
283 ENDDO
284 ENDDO
285
286 lstart = (/ 1_xt_int_kind, 1_xt_int_kind /)
287
288
289 global_section = xt_idxfsection_new(gstart, gsize, int(gsize), lstart)
290 indices = -1
292 p = 0
293 DO j = 1, gdy
294 DO i = 1, gdx
295 p = p + 1
296 IF (egis(i,j) /= indices(p)) CALL abort('(1) bad indices', __line__)
297 ENDDO
298 ENDDO
300
301
302 local_section = xt_idxfsection_new(gstart, gsize, (/ ldx, ldy /), lstart)
303 indices = -1
305 p = 0
306 DO j = 1, ldy
307 DO i = 1, ldx
308 p = p + 1
309 IF (egis(i,j) /= indices(p)) CALL abort('(2) bad indices', __line__)
310 ENDDO
311 ENDDO
313
314
315 local_section = xt_idxfsection_new(gstart, gsize, &
316 (/ -ldx, ldy /), lstart)
317 indices = -1
319 p = 0
320 DO j = 1, ldy
321 DO i = ldx, 1, -1
322 p = p + 1
323 IF (egis(i,j) /= indices(p)) CALL abort('(3) bad indices', __line__)
324 ENDDO
325 ENDDO
327
328
329 local_section = xt_idxfsection_new(gstart, gsize, (/ ldx, -ldy /), lstart)
330 indices = -1
332 p = 0
333 DO j = ldy, 1, -1
334 DO i = 1, ldx
335 p = p + 1
336 IF (egis(i,j) /= indices(p)) CALL abort('(4) bad indices', __line__)
337 ENDDO
338 ENDDO
340
341
342 local_section = xt_idxfsection_new(gstart, gsize, &
343 (/ -ldx, -ldy /), lstart)
344 indices = -1
346 p = 0
347 DO j = ldy, 1, -1
348 DO i = ldx, 1, -1
349 p = p + 1
350 IF (egis(i,j) /= indices(p)) CALL abort('(5) bad indices', __line__)
351 ENDDO
352 ENDDO
354 END SUBROUTINE general_fsection_test
355
356 SUBROUTINE gen_pos3d_surf(pos)
357 INTEGER, INTENT(inout) :: pos(:,:)
358
359
360
361 INTEGER :: ii,jj, i,j,k, p
362
363 k = 0
364 DO jj=1,je
365 DO ii=1,ie
366 p = pos(ii,jj) - 1
367 j = p/ie
368 i = mod(p,ie)
369 pos(ii,jj) = i + k*ie + j*ie*nlev
370 ENDDO
371 ENDDO
372
373 END SUBROUTINE gen_pos3d_surf
374
375 SUBROUTINE init_all
376 CHARACTER(len=*), PARAMETER :: context = 'init_all: '
377 INTEGER :: ierror
378
379 CALL init_mpi
380
382
383 CALL mpi_comm_size(mpi_comm_world, nprocs, ierror)
384 IF (ierror /= mpi_success) &
385 CALL abort(context//'MPI_COMM_SIZE failed', __line__)
386
387 CALL mpi_comm_rank(mpi_comm_world, mype, ierror)
388 IF (ierror /= mpi_success) &
389 CALL abort(context//'MPI_COMM_RANK failed', __line__)
390 IF (mype==0) THEN
391 lroot = .true.
392 ELSE
393 lroot = .false.
394 ENDIF
395
396 CALL factorize(nprocs, nprocx, nprocy)
397 IF (verbose .AND. lroot) WRITE(0,*) 'nprocx, nprocy=',nprocx, nprocy
398 mypy = mype / nprocx
399 mypx = mod(mype, nprocx)
400
401 CALL deco
402
403 ALLOCATE(fval(ie,je), gval(ie,je))
404 ALLOCATE(loc_id(ie,je), loc_tpex(ie,je))
405 ALLOCATE(id_pos(ie,je), gval3d(ie,nlev,je), pos3d_surf(ie,je))
406
407 fval = undef_int
408 gval = undef_int
409 loc_id = int(undef_int, xt_int_kind)
410 loc_tpex = int(undef_int, xt_int_kind)
411 id_pos = undef_int
412 gval3d = undef_int
413 pos3d_surf = undef_int
414
415 END SUBROUTINE init_all
416
417 SUBROUTINE gen_id_pos(pos)
418 INTEGER, INTENT(out) :: pos(:,:)
419
420 INTEGER :: i,j,p
421
422 p = 0
423 DO j = 1, SIZE(pos,2)
424 DO i = 1, SIZE(pos,1)
425 p = p + 1
426 pos(i,j) = p
427 ENDDO
428 ENDDO
429
430 END SUBROUTINE gen_id_pos
431
432 SUBROUTINE gen_trans(xmap, send_dt, recv_dt, redist)
433 TYPE(xt_xmap), INTENT(in) :: xmap
434 INTEGER,INTENT(in) :: send_dt, recv_dt
435 TYPE(xt_redist),INTENT(out) :: redist
436
437 INTEGER :: dt
438
439 IF (send_dt /= recv_dt) &
440 CALL abort('gen_trans: (send_dt /= recv_dt) unsupported', __line__)
441 dt = send_dt
443
444 END SUBROUTINE gen_trans
445
446 SUBROUTINE gen_off_trans(xmap, send_dt, send_off, recv_dt, recv_off, redist)
447 TYPE(xt_xmap), INTENT(in) :: xmap
448 INTEGER,INTENT(in) :: send_dt, recv_dt
449 INTEGER(c_int),INTENT(in) :: send_off(:,:), recv_off(:,:)
450 TYPE(xt_redist),INTENT(out) :: redist
451
452
453
454
455
456 IF (recv_dt /= send_dt) &
457 CALL abort('(datatype_in /= datatype_out) not supported', __line__)
458
460
461 END SUBROUTINE gen_off_trans
462
463 SUBROUTINE get_window(gval, win)
464 INTEGER(xt_int_kind), INTENT(in) :: gval(:,:)
465 INTEGER(xt_int_kind), INTENT(out) :: win(:,:)
466
467 INTEGER :: i, j, ig, jg
468
469 DO j = 1, je
470 jg = p_joff + j
471 DO i = 1, ie
472 ig = p_ioff + i
473 win(i,j) = gval(ig,jg)
474 ENDDO
475 ENDDO
476
477 END SUBROUTINE get_window
478
479 SUBROUTINE gen_template(local_src_idx, local_dst_idx, xmap)
480 INTEGER(xt_int_kind), INTENT(in) :: local_src_idx(:,:)
481 INTEGER(xt_int_kind), INTENT(in) :: local_dst_idx(:,:)
482 TYPE(xt_xmap), INTENT(out) :: xmap
483
484 TYPE(Xt_idxlist) :: src_idxlist, dst_idxlist
485 INTEGER :: src_num, dst_num
486 INTEGER(xt_int_kind) :: cp_src_idx(g_ie, g_je)
487 INTEGER(xt_int_kind) :: cp_dst_idx(g_ie, g_je)
488
491 IF (src_num /= g_ie*g_je) CALL abort('unexpected src_num', __line__)
493 IF (any(cp_src_idx /= local_src_idx)) &
494 CALL abort('idx copy does not match', __line__)
495
498 IF (dst_num /= g_ie*g_je) CALL abort('unexpected dst_num', __line__)
500 IF (any(cp_dst_idx /= local_dst_idx)) &
501 CALL abort('idx copy does not match', __line__)
502
506
507 END SUBROUTINE gen_template
508
509 SUBROUTINE def_tpex_mod_via_idxvec(mvec, mvec_num)
510 TYPE(xt_modifier), INTENT(out) :: mvec(:)
511 INTEGER, INTENT(out) :: mvec_num
512
513 INTEGER(xt_int_kind) :: g_start_indices(g_ie, g_je)
514 INTEGER(xt_int_kind) :: g_end_indices(g_ie, g_je)
515 TYPE(xt_idxlist) :: g_start_idxlist
516 TYPE(xt_idxlist) :: g_end_idxlist
517
518 IF (SIZE(mvec)<1) &
519 CALL abort('def_tpex_mod_via_idxvec mvec too small', __line__)
520
521 CALL id_map(g_start_indices)
522 g_start_idxlist =
xt_idxvec_new(g_start_indices,
SIZE(g_start_indices))
523
524 CALL def_exchange(g_start_indices, g_end_indices)
525 g_end_idxlist =
xt_idxvec_new(g_end_indices,
SIZE(g_end_indices))
526
527 mvec(1)%extract = g_start_idxlist
528 mvec(1)%subst = g_end_idxlist
529 mvec(1)%mask = 1
530 mvec_num = 1
531
532 END SUBROUTINE def_tpex_mod_via_idxvec
533
534 SUBROUTINE def_tpex_mod_via_sections(mvec, mvec_num)
535 TYPE(xt_modifier), INTENT(out) :: mvec(:)
536 INTEGER, INTENT(out) :: mvec_num
537
538 INTEGER(xt_int_kind), PARAMETER :: gstart_idx = 1_xt_int_kind
539 INTEGER(xt_int_kind), PARAMETER :: gsize(2) &
540 = (/ int(g_ie, xt_int_kind), int(g_je, xt_int_kind) /)
541 INTEGER :: ldx, ldy
542
543 INTEGER(xt_int_kind) :: g_core_is, g_core_ie, g_core_je
544 INTEGER(xt_int_kind) :: north_halo, im
545
546
547 g_core_is = nhalo + 1
548 g_core_ie = g_ie-nhalo
549 g_core_je = g_je-nhalo
550
551 im = 0
552
553
554 IF (with_north_halo) THEN
555
556
557 IF (increased_north_halo) THEN
558 north_halo = nhalo+1
559 ELSE
560 north_halo = nhalo
561 ENDIF
562
563 IF (2*north_halo > g_core_je) &
564 CALL test_abort('def_tpex_mod_via_sections: grid too small '//&
565 '(or halo too large) for tripolar north exchange',filename, __line__)
566
567 im = im + 1_xt_int_kind
568 IF (SIZE(mvec)<im) CALL abort('(SIZE(mvec)<im)', __line__)
569
570 ldx = int(g_core_ie - g_core_is + 1)
571 ldy = int(north_halo)
572 mvec(im)%extract = xt_idxfsection_new(gstart_idx, gsize, &
573 (/ ldx, ldy /), (/g_core_is, 1_xt_int_kind/))
574 mvec(im)%subst = xt_idxfsection_new(gstart_idx, gsize, &
575 (/ -ldx, -ldy /), (/g_core_is, north_halo+1_xt_int_kind/))
576 mvec(im)%mask = 1
577
578
579 im = im + 1_xt_int_kind
580 IF (SIZE(mvec)<im) CALL abort('(SIZE(mvec)<im)', __line__)
581 ldx = 1
582 ldy = int(north_halo)
583 mvec(im)%extract = xt_idxfsection_new(gstart_idx, gsize, &
584 (/ ldx, ldy /), (/1_xt_int_kind, 1_xt_int_kind/))
585 mvec(im)%subst = xt_idxfsection_new(gstart_idx, gsize, &
586 (/ -ldx, -ldy /), (/2_xt_int_kind, north_halo+1_xt_int_kind/))
587 mvec(im)%mask = 1
588
589 im = im + 1_xt_int_kind
590 IF (SIZE(mvec)<im) CALL abort('(SIZE(mvec)<im)', __line__)
591 ldx = 1
592 ldy = int(north_halo)
593 mvec(im)%extract = xt_idxfsection_new(gstart_idx, gsize, &
594 (/ ldx, ldy /), (/int(g_ie, xt_int_kind), 1_xt_int_kind/))
595 mvec(im)%subst = xt_idxfsection_new(gstart_idx, gsize, &
596 (/ -ldx, -ldy /), &
597 (/int(g_ie - 1, xt_int_kind), north_halo+1_xt_int_kind/))
598 mvec(im)%mask = 1
599
600 ELSE
601
602
603
604 ENDIF
605
606
607 ldx = nhalo
608 ldy = int(int(g_je, xt_int_kind) - north_halo)
609 im = im + 1_xt_int_kind
610 IF (SIZE(mvec)<im) CALL abort('(SIZE(mvec)<im)', __line__)
611 mvec(im)%extract = xt_idxfsection_new(gstart_idx, gsize, &
612 (/ ldx, ldy /), (/1_xt_int_kind, north_halo+1_xt_int_kind/))
613 mvec(im)%subst = xt_idxfsection_new(gstart_idx, gsize, &
614 (/ ldx, ldy /), &
615 (/ g_core_ie - int(ldx, xt_int_kind) + 1_xt_int_kind, &
616 north_halo + 1_xt_int_kind/))
617 mvec(im)%mask = 1
618
619 im = im + 1_xt_int_kind
620 IF (SIZE(mvec)<im) CALL abort('(SIZE(mvec)<im)', __line__)
621 mvec(im)%extract = xt_idxfsection_new(gstart_idx, gsize, &
622 (/ ldx, ldy /), (/g_core_ie+1_xt_int_kind, north_halo+1_xt_int_kind/))
623 mvec(im)%subst = xt_idxfsection_new(gstart_idx, gsize, &
624 (/ ldx, ldy /), (/ int(ldx, xt_int_kind) + 1_xt_int_kind, &
625 & north_halo+1_xt_int_kind/))
626 mvec(im)%mask = 1
627
628 mvec_num = int(im)
629
630 END SUBROUTINE def_tpex_mod_via_sections
631
632 SUBROUTINE def_exchange(id_in, id_out)
633 INTEGER(xt_int_kind), INTENT(in) :: id_in(:,:)
634 INTEGER(xt_int_kind), INTENT(out) :: id_out(:,:)
635
636 INTEGER :: i, j
637 INTEGER :: g_core_is, g_core_ie, g_core_js, g_core_je
638 INTEGER :: north_halo
639
640
641 g_core_is = nhalo + 1
642 g_core_ie = g_ie-nhalo
643 g_core_js = nhalo + 1
644 g_core_je = g_je-nhalo
645
646
647 id_out = undef_index
648 id_out(g_core_is:g_core_ie, g_core_js:g_core_je) &
649 = id_in(g_core_is:g_core_ie, g_core_js:g_core_je)
650
651 IF (with_north_halo) THEN
652
653
654 IF (increased_north_halo) THEN
655 north_halo = nhalo+1
656 ELSE
657 north_halo = nhalo
658 ENDIF
659
660 IF (2*north_halo > g_core_je) &
661 CALL test_abort('def_exchange: grid too small (or halo too large)'//&
662 'for tripolar north exchange', filename, __line__)
663 DO j = 1, north_halo
664 DO i = g_core_is, g_core_ie
665 id_out(i,j) = id_out(g_core_ie + (g_core_is-i), 2*north_halo + (1-j))
666 ENDDO
667 ENDDO
668
669 ELSE
670
671 DO j = 1, nhalo
672 DO i = nhalo+1, g_ie-nhalo
673 id_out(i,j) = id_in(i,j)
674 ENDDO
675 ENDDO
676
677 ENDIF
678
679
680 DO j = g_core_je+1, g_je
681 DO i = nhalo+1, g_ie-nhalo
682 id_out(i,j) = id_in(i,j)
683 ENDDO
684 ENDDO
685
686
687 DO j = 1, g_je
688 DO i = 1, nhalo
689 id_out(g_core_is-i,j) = id_out(g_core_ie+(1-i),j)
690 ENDDO
691 DO i = 1, nhalo
692 id_out(g_core_ie+i,j) = id_out(nhalo+i,j)
693 ENDDO
694 ENDDO
695
696 IF (any(id_out == undef_index)) &
697 CALL abort('found undefined indices', __line__)
698
699 END SUBROUTINE def_exchange
700
701 SUBROUTINE deco
702 INTEGER :: cx0(0:nprocx-1), cxn(0:nprocx-1)
703 INTEGER :: cy0(0:nprocy-1), cyn(0:nprocy-1)
704
705 CALL regular_deco(g_ie-2*nhalo, cx0, cxn)
706 CALL regular_deco(g_je-2*nhalo, cy0, cyn)
707
708
709 ie = cxn(mypx) + 2*nhalo
710 je = cyn(mypy) + 2*nhalo
711 p_ioff = cx0(mypx)
712 p_joff = cy0(mypy)
713
714 END SUBROUTINE deco
715
716END PROGRAM test_yaxt
717
718
719
720
721
722
723
724
725
void xt_initialize(MPI_Comm default_comm)
void xt_idxlist_get_indices(Xt_idxlist idxlist, Xt_int *indices)
void xt_idxlist_delete(Xt_idxlist idxlist)
#define xt_idxlist_get_num_indices(idxlist)
Xt_idxlist xt_idxmod_new(Xt_idxlist patch_idxlist, struct Xt_modifier *modifier, int modifier_num, int *mstate)
generates a new index list based on an index list and a sequence of modifiers
Xt_idxlist xt_idxvec_new(const Xt_int *idxlist, int num_indices)
void xt_redist_delete(Xt_redist redist)
void xt_redist_s_exchange(Xt_redist redist, int num_arrays, const void *const src_data[], void *const dst_data[])
Xt_redist xt_redist_collection_static_new(Xt_redist *redists, int num_redists, const MPI_Aint src_displacements[num_redists], const MPI_Aint dst_displacements[num_redists], MPI_Comm comm)
Xt_redist xt_redist_p2p_off_new(Xt_xmap xmap, const int *src_offsets, const int *dst_offsets, MPI_Datatype datatype)
Xt_redist xt_redist_p2p_new(Xt_xmap xmap, MPI_Datatype datatype)
void xt_xmap_delete(Xt_xmap xmap)
Xt_xmap xt_xmap_all2all_new(Xt_idxlist src_idxlist, Xt_idxlist dst_idxlist, MPI_Comm comm)