Yet Another eXchange Tool 0.11.4
Loading...
Searching...
No Matches
test_yaxt.f90
1
12
13!
14! Keywords:
15! Maintainer: Jörg Behrens <behrens@dkrz.de>
16! Moritz Hanke <hanke@dkrz.de>
17! Thomas Jahns <jahns@dkrz.de>
18! URL: https://dkrz-sw.gitlab-pages.dkrz.de/yaxt/
19!
20! Redistribution and use in source and binary forms, with or without
21! modification, are permitted provided that the following conditions are
22! met:
23!
24! Redistributions of source code must retain the above copyright notice,
25! this list of conditions and the following disclaimer.
26!
27! Redistributions in binary form must reproduce the above copyright
28! notice, this list of conditions and the following disclaimer in the
29! documentation and/or other materials provided with the distribution.
30!
31! Neither the name of the DKRZ GmbH nor the names of its contributors
32! may be used to endorse or promote products derived from this software
33! without specific prior written permission.
34!
35! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS
36! IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
37! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A
38! PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER
39! OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
40! EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
41! PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
42! PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
43! LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
44! NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
45! SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
46!
47#include "fc_feature_defs.inc"
48PROGRAM test_yaxt
49 USE iso_c_binding, ONLY: c_int
50 USE mpi
51 USE yaxt, ONLY: xt_idxlist, xt_idxlist_delete, xt_idxvec_new, &
57 xt_idxfsection_new, xt_redist_collection_static_new
58 USE ftest_common, ONLY: test_abort, id_map, factorize, regular_deco, &
59 init_mpi, finish_mpi, cmp_arrays
60
61 ! PGI compilers up to at least version 15 do not handle generic
62 ! interfaces correctly
63#if defined __PGI
64 USE xt_redist_int_i2, ONLY: xt_redist_s_exchange
65 USE xt_redist_int_i4, ONLY: xt_redist_s_exchange
66 USE xt_redist_int_i8, ONLY: xt_redist_s_exchange
67#endif
68 IMPLICIT NONE
69
70 INTEGER, PARAMETER :: g_ie = 8, g_je = 4! global extents including halos
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 ! 1dim. halo border size
76 LOGICAL, PARAMETER :: increased_north_halo = .false.
77 LOGICAL, PARAMETER :: with_north_halo = .true.
78
79 INTEGER :: ie, je ! local extents, including halos
80 INTEGER :: p_ioff, p_joff ! offsets within global domain
81 INTEGER :: nprocx, nprocy ! process space extents
82 INTEGER :: nprocs ! == nprocx*nprocy
83 INTEGER :: mype, mypx, mypy ! process rank, process coords within (0:, 0:) process space
84 LOGICAL :: lroot ! true only for proc 0
85
86 INTEGER(xt_int_kind) :: g_id(g_ie, g_je) ! global id
87 CHARACTER(len=*), PARAMETER :: filename = 'test_yaxt.f90'
88 ! global "tripolar-like" toy bounds exchange
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 ! mpi & decomposition & allocate mem:
100 CALL init_all
101
102 ! full global index space:
103 CALL id_map(g_id)
104
105 ! local window of global index space:
106 CALL get_window(g_id, loc_id)
107
108 ! define bounds exchange for full global index space
109 CALL def_exchange(g_id, g_tpex)
110
111 ! local window of global bounds exchange:
112 CALL get_window(g_tpex, loc_tpex)
113
114 ! check interface to idxsection:
115 CALL general_fsection_test
116
117 ! compare current index construction with modifier results
118 CALL check_modifiers
119
120 ! template: loc_id -> loc_tpex
121 CALL gen_template(loc_id, loc_tpex, xmap_tpex) ! todo rename template to xmap
122
123 ! transposition: loc_id:data -> loc_tpex:data
124 CALL gen_trans(xmap_tpex, mpi_integer, mpi_integer, redist_tpex)
125
126 ! test 2d-to-2d transposition:
127 fval = int(loc_id)
128 CALL xt_redist_s_exchange(redist_tpex, fval, gval)
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 ! define positions of surface elements within (i,k,j) array
137 CALL gen_id_pos(id_pos)
138 CALL gen_id_pos(pos3d_surf)
139 CALL gen_pos3d_surf(pos3d_surf)
140
141 ! generate surface transposition:
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 ! 2d to surface boundsexchange:
147 gval3d = -1
148 CALL xt_redist_s_exchange(redist_surf_tpex, &
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 ! check sub surface:
155 IF (any(gval3d(:, 2, :) /= -1)) &
156 CALL test_abort('surface check failed', filename, __line__)
157
158 ! cleanup:
159 DEALLOCATE(loc_tpex, gval3d)
160 CALL xt_xmap_delete(xmap_tpex)
161
162 CALL xt_redist_delete(redist_tpex)
163
164 CALL xt_redist_delete(redist_surf_tpex)
165
166 CALL xt_finalize()
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
193 rcol = xt_redist_collection_static_new(rvec, nr, f_disp, g_disp, mpi_comm_world)
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
200 CALL xt_redist_s_exchange(rvec(ir), f(:,:,ir), ref_g(:,:,ir))
201 ENDDO
202
203 g = 0
204 CALL xt_redist_s_exchange(rcol, f, g)
205 IF (any(g /= ref_g)) CALL abort('(g /= ref_g)', __line__)
206 CALL xt_redist_delete(rcol)
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
218 loc_id_idxlist = xt_idxvec_new(loc_id, SIZE(loc_id))
219 ALLOCATE(mstate(ie,je))
220
221 ! use one simple modifier to define index transfer:
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
226 CALL xt_idxlist_get_indices(loc_tpex2_idxlist, loc_tpex2)
227 IF (any(loc_tpex2 /= loc_tpex)) &
228 CALL abort('idx copy does not match', __line__)
229 CALL xt_idxlist_delete(loc_tpex2_idxlist)
230
231 ! test call without mstate
232 loc_tpex2_idxlist = xt_idxmod_new(loc_id_idxlist, m_tpex, m_tpex_num)
233 loc_tpex2 = -1
234 CALL xt_idxlist_get_indices(loc_tpex2_idxlist, loc_tpex2)
235 IF (any(loc_tpex2 /= loc_tpex)) &
236 CALL abort('idx copy does not match', __line__)
237 CALL xt_idxlist_delete(loc_tpex2_idxlist)
238 CALL delete_modifiers(m_tpex(1:m_tpex_num))
239
240 ! use compact modifiers to define index transfer:
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
244 CALL xt_idxlist_get_indices(loc_tpex2_idxlist, loc_tpex2)
245
246 IF (any(loc_tpex2 /= loc_tpex)) &
247 CALL abort('idx copy does not match', __line__)
248 CALL xt_idxlist_delete(loc_tpex2_idxlist)
249 CALL delete_modifiers(m_tpex(1:m_tpex_num))
250
251 ! cleanup:
252 CALL xt_idxlist_delete(loc_id_idxlist)
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)
261 CALL xt_idxlist_delete(m(i)%extract)
262 CALL xt_idxlist_delete(m(i)%subst)
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 ! prepare explicit global index space
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 ! check case: local section == global section
289 global_section = xt_idxfsection_new(gstart, gsize, int(gsize), lstart)
290 indices = -1
291 CALL xt_idxlist_get_indices(global_section, indices)
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
299 CALL xt_idxlist_delete(global_section)
300
301 ! check case: simple subsection
302 local_section = xt_idxfsection_new(gstart, gsize, (/ ldx, ldy /), lstart)
303 indices = -1
304 CALL xt_idxlist_get_indices(local_section, indices)
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
312 CALL xt_idxlist_delete(local_section)
313
314 ! check case: i-reverse subsection
315 local_section = xt_idxfsection_new(gstart, gsize, &
316 (/ -ldx, ldy /), lstart)
317 indices = -1
318 CALL xt_idxlist_get_indices(local_section, indices)
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
326 CALL xt_idxlist_delete(local_section)
327
328 ! check case: j-reverse subsection
329 local_section = xt_idxfsection_new(gstart, gsize, (/ ldx, -ldy /), lstart)
330 indices = -1
331 CALL xt_idxlist_get_indices(local_section, indices)
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
339 CALL xt_idxlist_delete(local_section)
340
341 ! check case: ij-reverse subsection
342 local_section = xt_idxfsection_new(gstart, gsize, &
343 (/ -ldx, -ldy /), lstart)
344 indices = -1
345 CALL xt_idxlist_get_indices(local_section, indices)
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
353 CALL xt_idxlist_delete(local_section)
354 END SUBROUTINE general_fsection_test
355
356 SUBROUTINE gen_pos3d_surf(pos)
357 INTEGER, INTENT(inout) :: pos(:,:)
358 ! positions for zero based arrays (ECHAM grid point dim order)
359 ! old pos = i + j*ie
360 ! new pos = i + k*ie + j*ie*nlev
361 INTEGER :: ii,jj, i,j,k, p
362
363 k = 0 ! surface
364 DO jj=1,je
365 DO ii=1,ie
366 p = pos(ii,jj) - 1 ! shift to 0-based index
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
381 CALL xt_initialize(mpi_comm_world)
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
442 redist = xt_redist_p2p_new(xmap, 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 !INTEGER :: send_offsets(SIZE(send_off)), recv_offsets(SIZE(recv_off))
453
454 !send_offsets = RESHAPE(send_off, (/SIZE(send_off)/) )
455 !recv_offsets = RESHAPE(recv_off, (/SIZE(recv_off)/) )
456 IF (recv_dt /= send_dt) &
457 CALL abort('(datatype_in /= datatype_out) not supported', __line__)
458
459 redist = xt_redist_p2p_off_new(xmap, send_off, recv_off, send_dt);
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
489 src_idxlist = xt_idxvec_new(local_src_idx, g_ie * g_je)
490 src_num = int(xt_idxlist_get_num_indices(src_idxlist))
491 IF (src_num /= g_ie*g_je) CALL abort('unexpected src_num', __line__)
492 CALL xt_idxlist_get_indices(src_idxlist, cp_src_idx)
493 IF (any(cp_src_idx /= local_src_idx)) &
494 CALL abort('idx copy does not match', __line__)
495
496 dst_idxlist = xt_idxvec_new(local_dst_idx, g_ie * g_je)
497 dst_num = int(xt_idxlist_get_num_indices(dst_idxlist))
498 IF (dst_num /= g_ie*g_je) CALL abort('unexpected dst_num', __line__)
499 CALL xt_idxlist_get_indices(dst_idxlist, cp_dst_idx)
500 IF (any(cp_dst_idx /= local_dst_idx)) &
501 CALL abort('idx copy does not match', __line__)
502
503 xmap = xt_xmap_all2all_new(src_idxlist, dst_idxlist, mpi_comm_world)
504 CALL xt_idxlist_delete(src_idxlist)
505 CALL xt_idxlist_delete(dst_idxlist)
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 ! global core domain:
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 ! global tripolar boundsexchange:
554 IF (with_north_halo) THEN
555
556 ! north inversion, (maybe with increased north halo)
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 ! north border exchange without ew-halos
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 ! 1. north edge:
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 ! 2. north edge:
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 ! nothing to do at the north border
603
604 ENDIF
605
606 ! PBC below north border
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 ! global core domain:
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 ! global tripolar boundsexchange:
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 ! north inversion, (maybe with increased north halo)
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 ! south: no change
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 ! PBC
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 ! process local deco variables:
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! Local Variables:
719! f90-continuation-indent: 5
720! coding: utf-8
721! indent-tabs-mode: nil
722! show-trailing-whitespace: t
723! require-trailing-newline: t
724! End:
725!
void xt_initialize(MPI_Comm default_comm)
Definition xt_init.c:70
void xt_finalize(void)
Definition xt_init.c:92
void xt_idxlist_get_indices(Xt_idxlist idxlist, Xt_int *indices)
Definition xt_idxlist.c:113
void xt_idxlist_delete(Xt_idxlist idxlist)
Definition xt_idxlist.c:75
#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
Definition xt_idxmod.c:62
Xt_idxlist xt_idxvec_new(const Xt_int *idxlist, int num_indices)
Definition xt_idxvec.c:213
void xt_redist_delete(Xt_redist redist)
Definition xt_redist.c:74
void xt_redist_s_exchange(Xt_redist redist, int num_arrays, const void *const src_data[], void *const dst_data[])
Definition xt_redist.c:79
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)
Definition xt_xmap.c:86
Xt_xmap xt_xmap_all2all_new(Xt_idxlist src_idxlist, Xt_idxlist dst_idxlist, MPI_Comm comm)