49 USE iso_c_binding,
ONLY: c_int
58 USE ftest_common,
ONLY: test_abort, id_map, factorize, regular_deco, &
59 finish_mpi, cmp_arrays
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.
80 INTEGER :: p_ioff, p_joff
81 INTEGER :: nprocx, nprocy
83 INTEGER :: mype, mypx, mypy
86 INTEGER(xt_int_kind) :: g_id(g_ie, g_je)
87 CHARACTER(len=*),
PARAMETER :: filename =
'test_yaxt.f90'
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
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(:,:)
106 CALL get_window(g_id, loc_id)
109 CALL def_exchange(g_id, g_tpex)
112 CALL get_window(g_tpex, loc_tpex)
115 CALL general_fsection_test
121 CALL gen_template(loc_id, loc_tpex, xmap_tpex)
124 CALL gen_trans(xmap_tpex, mpi_integer, mpi_integer, redist_tpex)
130 IF (cmp_arrays(gval, loc_tpex)) &
131 CALL test_abort(
'array eqivalence test failed', filename, __line__)
133 CALL check_redist_collection_static
137 CALL gen_id_pos(id_pos)
138 CALL gen_id_pos(pos3d_surf)
139 CALL gen_pos3d_surf(pos3d_surf)
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)
149 reshape(fval, (/ ie, 1, je /)), gval3d)
152 IF (cmp_arrays(gval3d(:, 1, :), loc_tpex)) &
153 CALL test_abort(
'surface check failed', filename, __line__)
155 IF (any(gval3d(:, 2, :) /= -1)) &
156 CALL test_abort(
'surface check failed', filename, __line__)
159 DEALLOCATE(loc_tpex, gval3d)
171 #define abort(msg, line) test_abort(msg, filename, line)
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)
180 INTEGER :: ir, ierror
181 rvec(:) = redist_tpex
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)
195 f(:,:,ir) = int(loc_id) + (ir-1) * ie*je
205 IF (any(g /= ref_g))
CALL abort(
'(g /= ref_g)', __line__)
208 END SUBROUTINE check_redist_collection_static
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(:,:)
219 ALLOCATE(mstate(ie,je))
222 CALL def_tpex_mod_via_idxvec(m_tpex, m_tpex_num)
224 loc_tpex2_idxlist =
xt_idxmod_new(loc_id_idxlist, m_tpex, m_tpex_num, mstate)
227 IF (any(loc_tpex2 /= loc_tpex)) &
228 CALL abort(
'idx copy does not match', __line__)
232 loc_tpex2_idxlist =
xt_idxmod_new(loc_id_idxlist, m_tpex, m_tpex_num)
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))
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)
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))
253 END SUBROUTINE check_modifiers
255 SUBROUTINE delete_modifiers(m)
256 TYPE(xt_modifier),
INTENT(inout) :: m(:)
265 END SUBROUTINE delete_modifiers
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
286 lstart = (/ 1_xt_int_kind, 1_xt_int_kind /)
289 global_section = xt_idxfsection_new(gstart, gsize, int(gsize), lstart)
296 IF (egis(i,j) /= indices(p))
CALL abort(
'(1) bad indices', __line__)
302 local_section = xt_idxfsection_new(gstart, gsize, (/ ldx, ldy /), lstart)
309 IF (egis(i,j) /= indices(p))
CALL abort(
'(2) bad indices', __line__)
315 local_section = xt_idxfsection_new(gstart, gsize, &
316 (/ -ldx, ldy /), lstart)
323 IF (egis(i,j) /= indices(p))
CALL abort(
'(3) bad indices', __line__)
329 local_section = xt_idxfsection_new(gstart, gsize, (/ ldx, -ldy /), lstart)
336 IF (egis(i,j) /= indices(p))
CALL abort(
'(4) bad indices', __line__)
342 local_section = xt_idxfsection_new(gstart, gsize, &
343 (/ -ldx, -ldy /), lstart)
350 IF (egis(i,j) /= indices(p))
CALL abort(
'(5) bad indices', __line__)
354 END SUBROUTINE general_fsection_test
356 SUBROUTINE gen_pos3d_surf(pos)
357 INTEGER,
INTENT(inout) :: pos(:,:)
361 INTEGER :: ii,jj, i,j,k, p
369 pos(ii,jj) = i + k*ie + j*ie*nlev
373 END SUBROUTINE gen_pos3d_surf
376 CHARACTER(len=*),
PARAMETER :: context =
'init_all: '
379 CALL mpi_init(ierror)
380 IF (ierror /= mpi_success) &
381 CALL abort(context//
'MPI_INIT failed', __line__)
385 CALL mpi_comm_size(mpi_comm_world, nprocs, ierror)
386 IF (ierror /= mpi_success) &
387 CALL abort(context//
'MPI_COMM_SIZE failed', __line__)
389 CALL mpi_comm_rank(mpi_comm_world, mype, ierror)
390 IF (ierror /= mpi_success) &
391 CALL abort(context//
'MPI_COMM_RANK failed', __line__)
398 CALL factorize(nprocs, nprocx, nprocy)
399 IF (verbose .AND. lroot)
WRITE(0,*)
'nprocx, nprocy=',nprocx, nprocy
401 mypx = mod(mype, nprocx)
408 ALLOCATE(fval(ie,je), gval(ie,je))
409 ALLOCATE(loc_id(ie,je), loc_tpex(ie,je))
410 ALLOCATE(id_pos(ie,je), gval3d(ie,nlev,je), pos3d_surf(ie,je))
414 loc_id = int(undef_int, xt_int_kind)
415 loc_tpex = int(undef_int, xt_int_kind)
418 pos3d_surf = undef_int
420 END SUBROUTINE init_all
422 SUBROUTINE gen_id_pos(pos)
423 INTEGER,
INTENT(out) :: pos(:,:)
428 DO j = 1,
SIZE(pos,2)
429 DO i = 1,
SIZE(pos,1)
435 END SUBROUTINE gen_id_pos
437 SUBROUTINE gen_trans(xmap, send_dt, recv_dt, redist)
438 TYPE(xt_xmap),
INTENT(in) :: xmap
439 INTEGER,
INTENT(in) :: send_dt, recv_dt
440 TYPE(xt_redist),
INTENT(out) :: redist
444 IF (send_dt /= recv_dt) &
445 CALL abort(
'gen_trans: (send_dt /= recv_dt) unsupported', __line__)
450 END SUBROUTINE gen_trans
452 SUBROUTINE gen_off_trans(xmap, send_dt, send_off, recv_dt, recv_off, redist)
453 TYPE(xt_xmap),
INTENT(in) :: xmap
454 INTEGER,
INTENT(in) :: send_dt, recv_dt
455 INTEGER(c_int),
INTENT(in) :: send_off(:,:), recv_off(:,:)
456 TYPE(xt_redist),
INTENT(out) :: redist
462 IF (recv_dt /= send_dt) &
463 CALL abort(
'(datatype_in /= datatype_out) not supported', __line__)
468 END SUBROUTINE gen_off_trans
470 SUBROUTINE get_window(gval, win)
471 INTEGER(xt_int_kind),
INTENT(in) :: gval(:,:)
472 INTEGER(xt_int_kind),
INTENT(out) :: win(:,:)
474 INTEGER :: i, j, ig, jg
480 win(i,j) = gval(ig,jg)
484 END SUBROUTINE get_window
486 SUBROUTINE gen_template(local_src_idx, local_dst_idx, xmap)
487 INTEGER(xt_int_kind),
INTENT(in) :: local_src_idx(:,:)
488 INTEGER(xt_int_kind),
INTENT(in) :: local_dst_idx(:,:)
489 TYPE(xt_xmap),
INTENT(out) :: xmap
491 TYPE(Xt_idxlist) :: src_idxlist, dst_idxlist
492 INTEGER :: src_num, dst_num
493 INTEGER(xt_int_kind) :: cp_src_idx(g_ie, g_je)
494 INTEGER(xt_int_kind) :: cp_dst_idx(g_ie, g_je)
498 IF (src_num /= g_ie*g_je)
CALL abort(
'unexpected src_num', __line__)
500 IF (any(cp_src_idx /= local_src_idx)) &
501 CALL abort(
'idx copy does not match', __line__)
505 IF (dst_num /= g_ie*g_je)
CALL abort(
'unexpected dst_num', __line__)
507 IF (any(cp_dst_idx /= local_dst_idx)) &
508 CALL abort(
'idx copy does not match', __line__)
514 END SUBROUTINE gen_template
516 SUBROUTINE def_tpex_mod_via_idxvec(mvec, mvec_num)
517 TYPE(xt_modifier),
INTENT(out) :: mvec(:)
518 INTEGER,
INTENT(out) :: mvec_num
520 INTEGER(xt_int_kind) :: g_start_indices(g_ie, g_je)
521 INTEGER(xt_int_kind) :: g_end_indices(g_ie, g_je)
522 TYPE(xt_idxlist) :: g_start_idxlist
523 TYPE(xt_idxlist) :: g_end_idxlist
526 CALL abort(
'def_tpex_mod_via_idxvec mvec too small', __line__)
528 CALL id_map(g_start_indices)
529 g_start_idxlist =
xt_idxvec_new(g_start_indices,
SIZE(g_start_indices))
531 CALL def_exchange(g_start_indices, g_end_indices)
532 g_end_idxlist =
xt_idxvec_new(g_end_indices,
SIZE(g_end_indices))
534 mvec(1)%extract = g_start_idxlist
535 mvec(1)%subst = g_end_idxlist
539 END SUBROUTINE def_tpex_mod_via_idxvec
541 SUBROUTINE def_tpex_mod_via_sections(mvec, mvec_num)
542 TYPE(xt_modifier),
INTENT(out) :: mvec(:)
543 INTEGER,
INTENT(out) :: mvec_num
545 INTEGER(xt_int_kind),
PARAMETER :: gstart_idx = 1_xt_int_kind
546 INTEGER(xt_int_kind),
PARAMETER :: gsize(2) &
547 = (/ int(g_ie, xt_int_kind), int(g_je, xt_int_kind) /)
550 INTEGER(xt_int_kind) :: g_core_is, g_core_ie, g_core_je
551 INTEGER(xt_int_kind) :: north_halo, im
554 g_core_is = nhalo + 1
555 g_core_ie = g_ie-nhalo
556 g_core_je = g_je-nhalo
561 IF (with_north_halo)
THEN
564 IF (increased_north_halo)
THEN
570 IF (2*north_halo > g_core_je) &
571 CALL test_abort(
'def_tpex_mod_via_sections: grid too small '//&
572 '(or halo too large) for tripolar north exchange',filename, __line__)
574 im = im + 1_xt_int_kind
575 IF (
SIZE(mvec)<im)
CALL abort(
'(SIZE(mvec)<im)', __line__)
577 ldx = int(g_core_ie - g_core_is + 1)
578 ldy = int(north_halo)
579 mvec(im)%extract = xt_idxfsection_new(gstart_idx, gsize, &
580 (/ ldx, ldy /), (/g_core_is, 1_xt_int_kind/))
581 mvec(im)%subst = xt_idxfsection_new(gstart_idx, gsize, &
582 (/ -ldx, -ldy /), (/g_core_is, north_halo+1_xt_int_kind/))
586 im = im + 1_xt_int_kind
587 IF (
SIZE(mvec)<im)
CALL abort(
'(SIZE(mvec)<im)', __line__)
589 ldy = int(north_halo)
590 mvec(im)%extract = xt_idxfsection_new(gstart_idx, gsize, &
591 (/ ldx, ldy /), (/1_xt_int_kind, 1_xt_int_kind/))
592 mvec(im)%subst = xt_idxfsection_new(gstart_idx, gsize, &
593 (/ -ldx, -ldy /), (/2_xt_int_kind, north_halo+1_xt_int_kind/))
596 im = im + 1_xt_int_kind
597 IF (
SIZE(mvec)<im)
CALL abort(
'(SIZE(mvec)<im)', __line__)
599 ldy = int(north_halo)
600 mvec(im)%extract = xt_idxfsection_new(gstart_idx, gsize, &
601 (/ ldx, ldy /), (/int(g_ie, xt_int_kind), 1_xt_int_kind/))
602 mvec(im)%subst = xt_idxfsection_new(gstart_idx, gsize, &
604 (/int(g_ie - 1, xt_int_kind), north_halo+1_xt_int_kind/))
615 ldy = int(int(g_je, xt_int_kind) - north_halo)
616 im = im + 1_xt_int_kind
617 IF (
SIZE(mvec)<im)
CALL abort(
'(SIZE(mvec)<im)', __line__)
618 mvec(im)%extract = xt_idxfsection_new(gstart_idx, gsize, &
619 (/ ldx, ldy /), (/1_xt_int_kind, north_halo+1_xt_int_kind/))
620 mvec(im)%subst = xt_idxfsection_new(gstart_idx, gsize, &
622 (/ g_core_ie - int(ldx, xt_int_kind) + 1_xt_int_kind, &
623 north_halo + 1_xt_int_kind/))
626 im = im + 1_xt_int_kind
627 IF (
SIZE(mvec)<im)
CALL abort(
'(SIZE(mvec)<im)', __line__)
628 mvec(im)%extract = xt_idxfsection_new(gstart_idx, gsize, &
629 (/ ldx, ldy /), (/g_core_ie+1_xt_int_kind, north_halo+1_xt_int_kind/))
630 mvec(im)%subst = xt_idxfsection_new(gstart_idx, gsize, &
631 (/ ldx, ldy /), (/ int(ldx, xt_int_kind) + 1_xt_int_kind, &
632 & north_halo+1_xt_int_kind/))
637 END SUBROUTINE def_tpex_mod_via_sections
639 SUBROUTINE def_exchange(id_in, id_out)
640 INTEGER(xt_int_kind),
INTENT(in) :: id_in(:,:)
641 INTEGER(xt_int_kind),
INTENT(out) :: id_out(:,:)
644 INTEGER :: g_core_is, g_core_ie, g_core_js, g_core_je
645 INTEGER :: north_halo
648 g_core_is = nhalo + 1
649 g_core_ie = g_ie-nhalo
650 g_core_js = nhalo + 1
651 g_core_je = g_je-nhalo
655 id_out(g_core_is:g_core_ie, g_core_js:g_core_je) &
656 = id_in(g_core_is:g_core_ie, g_core_js:g_core_je)
658 IF (with_north_halo)
THEN
661 IF (increased_north_halo)
THEN
667 IF (2*north_halo > g_core_je) &
668 CALL test_abort(
'def_exchange: grid too small (or halo too large)'//&
669 'for tripolar north exchange', filename, __line__)
671 DO i = g_core_is, g_core_ie
672 id_out(i,j) = id_out(g_core_ie + (g_core_is-i), 2*north_halo + (1-j))
679 DO i = nhalo+1, g_ie-nhalo
680 id_out(i,j) = id_in(i,j)
687 DO j = g_core_je+1, g_je
688 DO i = nhalo+1, g_ie-nhalo
689 id_out(i,j) = id_in(i,j)
696 id_out(g_core_is-i,j) = id_out(g_core_ie+(1-i),j)
699 id_out(g_core_ie+i,j) = id_out(nhalo+i,j)
703 IF (any(id_out == undef_index)) &
704 CALL abort(
'found undefined indices', __line__)
706 END SUBROUTINE def_exchange
709 INTEGER :: cx0(0:nprocx-1), cxn(0:nprocx-1)
710 INTEGER :: cy0(0:nprocy-1), cyn(0:nprocy-1)
712 CALL regular_deco(g_ie-2*nhalo, cx0, cxn)
713 CALL regular_deco(g_je-2*nhalo, cy0, cyn)
716 ie = cxn(mypx) + 2*nhalo
717 je = cyn(mypy) + 2*nhalo
723 END PROGRAM test_yaxt
void xt_initialize(MPI_Comm default_comm)
int xt_idxlist_get_num_indices(Xt_idxlist idxlist)
void xt_idxlist_get_indices(Xt_idxlist idxlist, Xt_int *indices)
void xt_idxlist_delete(Xt_idxlist 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 **src_data, void **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)