33 #include <fms_platform.h>
44 logical :: verbose=.false.
45 logical :: module_is_initialized=.false.
46 character(len=256) :: text
49 module procedure mpp_pset_broadcast_ptr_scalar
50 module procedure mpp_pset_broadcast_ptr_array
53 module procedure mpp_send_ptr_scalar
54 module procedure mpp_send_ptr_array
57 module procedure mpp_recv_ptr_scalar
58 module procedure mpp_recv_ptr_array
61 module procedure mpp_pset_print_chksum_1d
62 module procedure mpp_pset_print_chksum_2d
63 module procedure mpp_pset_print_chksum_3d
64 module procedure mpp_pset_print_chksum_4d
70 integer :: next_in_pset, prev_in_pset
71 integer :: root_in_pset
76 real,
allocatable :: stack(:)
77 integer,
allocatable :: pelist(:)
78 integer,
allocatable :: root_pelist(:)
79 integer,
allocatable :: pset(:)
80 integer(POINTER_KIND) :: p_stack
81 integer :: lstack, maxstack, hiWM
83 character(len=32) :: name
84 logical :: initialized=.false.
90 public :: mpp_pset_create, mpp_pset_sync, mpp_pset_broadcast, &
93 mpp_pset_delete, mpp_pset_root, mpp_pset_numroots, mpp_pset_init, &
94 mpp_pset_get_root_pelist, mpp_pset_print_stack_chksum
98 subroutine mpp_pset_init
99 module_is_initialized = .true.
100 end subroutine mpp_pset_init
102 subroutine mpp_pset_create(npset,pset,stacksize,pelist, commID)
106 integer,
intent(in) :: npset
107 type(mpp_pset_type),
intent(inout) :: pset
108 integer,
intent(in),
optional :: stacksize
109 integer,
intent(in),
optional :: pelist(:)
110 integer,
intent(in),
optional :: commID
112 integer :: npes, my_commID
113 integer :: i, j, k, out_unit, errunit
114 integer,
allocatable :: my_pelist(:), root_pelist(:)
125 if(
present(pelist))
then
126 npes =
size(pelist(:))
130 if( mod(npes,npset).NE.0 )
then
131 write( text,
'(a,2i6)' ) &
132 'MPP_PSET_CREATE: PSET size (npset) must divide npes exactly:'// &
133 ' npset, npes=', npset, npes
138 allocate(my_pelist(0:npes-1) )
139 allocate(root_pelist(0:npes/npset-1) )
140 if(
present(pelist))
then
141 if(.not.
present(commid))
call mpp_error(fatal, &
142 'MPP_PSET_CREATE: when pelist is present, commID should also be present')
146 call mpp_get_current_pelist(my_pelist, commid = my_commid)
148 do i = 0,npes/npset-1
149 root_pelist(i) = my_pelist(npset*i)
151 write( out_unit,
'(a,i6)' )
'MPP_PSET_CREATE creating PSETs... npset=', npset
152 if(any(my_pelist == pe) )
then
153 if( pset%initialized )
call mpp_error( fatal, &
154 'MPP_PSET_CREATE: PSET already initialized!' )
156 allocate( pset%pelist(0:npes-1) )
157 allocate( pset%root_pelist(0:npes/npset-1) )
158 pset%commID = my_commid
159 pset%pelist = my_pelist
161 pset%root_pelist = root_pelist
162 allocate( pset%pset(0:npset-1) )
163 do i = 0,npes/npset-1
167 if( pe.EQ.pset%pelist(k+j) )
then
168 pset%pset(:) = pset%pelist(k:k+npset-1)
170 pset%root_in_pset = pset%root_pelist(i)
172 pset%prev_in_pset = pset%pelist(k+npset-1)
174 pset%prev_in_pset = pset%pelist(k+j-1)
176 if( j.EQ.npset-1 )
then
177 pset%next_in_pset = pset%pelist(k)
179 pset%next_in_pset = pset%pelist(k+j+1)
185 pset%root = pe.EQ.pset%root_in_pset
189 pset%maxstack = 1000000
190 if(
PRESENT(stacksize) )pset%maxstack = stacksize
191 write( out_unit,
'(a,i8)' ) &
192 'MPP_PSET_CREATE: setting stacksize=', pset%maxstack
194 allocate( pset%stack(pset%maxstack) )
195 #ifdef use_CRI_pointers
196 pset%p_stack = loc(pset%stack)
199 pset%initialized = .true.
206 write( errunit,
'(a,4i6)' )
'MPP_PSET_CREATE: pe, root, next, prev=', &
207 pe, pset%root_in_pset, pset%next_in_pset, pset%prev_in_pset
208 write( errunit,* )
'PE ', pe,
' pset=', pset%pset(:)
209 write( out_unit,* )
'root pelist=', pset%root_pelist(:)
211 end subroutine mpp_pset_create
213 subroutine mpp_pset_delete(pset)
214 type(mpp_pset_type),
intent(inout) :: pset
218 if( .NOT.pset%initialized )
call mpp_error( fatal, &
219 'MPP_PSET_DELETE: called with uninitialized PSET.' )
221 deallocate( pset%pelist )
222 deallocate( pset%root_pelist )
223 deallocate( pset%pset )
224 if( pset%root )
deallocate( pset%stack )
225 write( out_unit,
'(a,i10)' ) &
226 'Deleting PSETs... stack high-water-mark=', pset%hiWM
228 pset%initialized = .false.
229 end subroutine mpp_pset_delete
231 subroutine mpp_send_ptr_scalar( ptr, pe )
232 integer(POINTER_KIND),
intent(in) :: ptr
233 integer,
intent(in) :: pe
237 call mpp_send( ptr, pe, tag=comm_tag_1 )
238 end subroutine mpp_send_ptr_scalar
240 subroutine mpp_send_ptr_array( ptr, pe )
241 integer(POINTER_KIND),
intent(in) :: ptr(:)
242 integer,
intent(in) :: pe
246 call mpp_send( ptr,
size(ptr), pe, tag=comm_tag_2 )
247 end subroutine mpp_send_ptr_array
249 subroutine mpp_recv_ptr_scalar( ptr, pe )
250 integer(POINTER_KIND),
intent(inout) :: ptr
251 integer,
intent(in) :: pe
253 call mpp_recv( ptr, pe, tag=comm_tag_1 )
254 call mpp_translate_remote_ptr( ptr, pe )
256 end subroutine mpp_recv_ptr_scalar
258 subroutine mpp_recv_ptr_array( ptr, pe )
259 integer(POINTER_KIND),
intent(inout) :: ptr(:)
260 integer,
intent(in) :: pe
263 call mpp_recv( ptr,
size(ptr), pe, tag=comm_tag_2 )
265 call mpp_translate_remote_ptr( ptr(i), pe )
268 end subroutine mpp_recv_ptr_array
270 subroutine mpp_translate_remote_ptr( ptr, pe )
272 integer(POINTER_KIND),
intent(inout) :: ptr
273 integer,
intent(in) :: pe
275 end subroutine mpp_translate_remote_ptr
277 subroutine mpp_pset_sync(pset)
280 type(mpp_pset_type),
intent(in) :: pset
282 if( .NOT.pset%initialized )
call mpp_error( fatal, &
283 'MPP_PSET_SYNC: called with uninitialized PSET.' )
287 end subroutine mpp_pset_sync
289 subroutine mpp_pset_broadcast(pset,a)
291 type(mpp_pset_type),
intent(in) :: pset
292 real,
intent(inout) :: a
295 if( .NOT.pset%initialized )
call mpp_error( fatal, &
296 'MPP_PSET_BROADCAST: called with uninitialized PSET.' )
298 do i = 1,pset%npset-1
299 call mpp_send( a, pset%pset(i), tag=comm_tag_3 )
302 call mpp_recv( a, pset%root_in_pset, tag=comm_tag_3 )
304 call mpp_pset_sync(pset)
305 end subroutine mpp_pset_broadcast
307 subroutine mpp_pset_broadcast_ptr_scalar(pset,ptr)
311 type(mpp_pset_type),
intent(in) :: pset
312 integer(POINTER_KIND),
intent(inout) :: ptr
315 if( .NOT.pset%initialized )
call mpp_error( fatal, &
316 'MPP_PSET_BROADCAST_PTR: called with uninitialized PSET.' )
319 do i = 1,pset%npset-1
326 end subroutine mpp_pset_broadcast_ptr_scalar
328 subroutine mpp_pset_broadcast_ptr_array(pset,ptr)
332 type(mpp_pset_type),
intent(in) :: pset
333 integer(POINTER_KIND),
intent(inout) :: ptr(:)
336 if( .NOT.pset%initialized )
call mpp_error( fatal, &
337 'MPP_PSET_BROADCAST_PTR: called with uninitialized PSET.' )
340 do i = 1,pset%npset-1
348 end subroutine mpp_pset_broadcast_ptr_array
350 subroutine mpp_pset_check_ptr(pset,ptr)
352 type(mpp_pset_type),
intent(in) :: pset
353 #ifdef use_CRI_pointers
355 pointer( ptr, dummy )
357 integer(POINTER_KIND),
intent(in) :: ptr
360 integer(POINTER_KIND) :: p
362 if( .NOT.pset%initialized )
call mpp_error( fatal, &
363 'MPP_PSET_CHECK_PTR: called with uninitialized PSET.' )
368 do i = 1,pset%npset-1
374 call mpp_pset_sync(pset)
376 'MPP_PSET_CHECK_PTR: pointers do not match!' )
380 end subroutine mpp_pset_check_ptr
382 subroutine mpp_pset_segment_array( pset, ls, le, lsp, lep )
387 type(mpp_pset_type),
intent(in) :: pset
388 integer,
intent(in) :: ls, le
389 integer,
intent(out) :: lsp, lep
392 if( .NOT.pset%initialized )
call mpp_error( fatal, &
393 'MPP_PSET_SEGMENT_ARRAY: called with uninitialized PSET.' )
395 if( le-ls+1.LT.pset%npset )
then
396 write( text,
'(3(a,i6))' ) &
397 'MPP_PSET_ARRAY_SEGMENT: parallel range (', ls,
',', le, &
398 ') is smaller than the number of threads:', pset%npset
405 lep = lsp + ceiling( real(le-lsp+1)/(pset%npset-i) ) - 1
407 end subroutine mpp_pset_segment_array
409 subroutine mpp_pset_stack_push( pset, ptr, len )
414 type(mpp_pset_type),
intent(inout) :: pset
415 integer,
intent(in) :: len
416 #ifdef use_CRI_pointers
418 pointer( ptr, dummy )
419 real :: stack(pset%maxstack)
422 if( .NOT.pset%initialized )
call mpp_error( fatal, &
423 'MPP_PSET_STACK_PUSH: called with uninitialized PSET.' )
424 if( pset%lstack+len.GT.pset%maxstack )
then
425 write( text,
'(a,3i12)' ) &
426 'MPP_PSET_STACK_PUSH: mpp_pset_stack overflow: '// &
427 .GT.
'len+lstackmaxstack. len, lstack, maxstack=', &
428 len, pset%lstack, pset%maxstack
432 ptr = loc( stack(pset%lstack+1) )
433 call mpp_pset_check_ptr(pset,ptr)
434 pset%lstack = pset%lstack + len
435 pset%hiWM = max( pset%hiWM, pset%lstack )
437 integer(POINTER_KIND),
intent(out) :: ptr
439 'MPP_PSET_STACK_PUSH only works with Cray pointers.' )
441 end subroutine mpp_pset_stack_push
443 subroutine mpp_pset_stack_reset(pset)
444 type(mpp_pset_type),
intent(inout) :: pset
454 if( .NOT.pset%initialized )
call mpp_error( fatal, &
455 'MPP_PSET_STACK_RESET: called with uninitialized PSET.' )
457 end subroutine mpp_pset_stack_reset
459 subroutine mpp_pset_print_chksum_1d(pset, caller, array)
463 type(mpp_pset_type),
intent(in) :: pset
464 character(len=*),
intent(in) :: caller
465 real,
intent(in) :: array(:)
470 integer(LONG_KIND) :: chksum
472 if( .NOT.pset%initialized )
call mpp_error( fatal, &
473 'MPP_PSET_PRINT_CHKSUM: called with uninitialized PSET.' )
477 do_print = pe.EQ.mpp_root_pe()
481 write( errunit,
'(a,z18)' )trim(caller)//
' chksum=', chksum
486 end subroutine mpp_pset_print_chksum_1d
488 subroutine mpp_pset_print_chksum_2d(pset, caller, array)
489 type(mpp_pset_type),
intent(in) :: pset
490 character(len=*),
intent(in) :: caller
491 real,
intent(in) :: array(:,:)
492 real :: array1D( size(array) )
493 #ifdef use_CRI_pointers
494 pointer( p, array1d )
497 array1d = transfer( array, array1d )
500 end subroutine mpp_pset_print_chksum_2d
502 subroutine mpp_pset_print_chksum_3d(pset, caller, array)
503 type(mpp_pset_type),
intent(in) :: pset
504 character(len=*),
intent(in) :: caller
505 real,
intent(in) :: array(:,:,:)
506 real :: array1D( size(array) )
507 #ifdef use_CRI_pointers
508 pointer( p, array1d )
511 array1d = transfer( array, array1d )
514 end subroutine mpp_pset_print_chksum_3d
516 subroutine mpp_pset_print_chksum_4d(pset, caller, array)
517 type(mpp_pset_type),
intent(in) :: pset
518 character(len=*),
intent(in) :: caller
519 real,
intent(in) :: array(:,:,:,:)
520 real :: array1D( size(array) )
521 #ifdef use_CRI_pointers
522 pointer( p, array1d )
525 array1d = transfer( array, array1d )
528 end subroutine mpp_pset_print_chksum_4d
530 subroutine mpp_pset_print_stack_chksum( pset, caller )
531 type(mpp_pset_type),
intent(in) :: pset
532 character(len=*),
intent(in) :: caller
534 if( .NOT.pset%initialized )
call mpp_error( fatal, &
535 'MPP_PSET_PRINT_STACK_CHKSUM: called with uninitialized PSET.' )
537 pset%stack(1:pset%lstack) )
538 end subroutine mpp_pset_print_stack_chksum
541 function mpp_pset_root(pset)
542 logical :: mpp_pset_root
543 type(mpp_pset_type),
intent(in) :: pset
545 if( .NOT.pset%initialized )
call mpp_error( fatal, &
546 'MPP_PSET_ROOT: called with uninitialized PSET.' )
547 mpp_pset_root = pset%root
548 end function mpp_pset_root
550 function mpp_pset_numroots(pset)
552 integer :: mpp_pset_numroots
553 type(mpp_pset_type),
intent(in) :: pset
555 if( .NOT.pset%initialized )
call mpp_error( fatal, &
556 'MPP_PSET_NUMROOTS: called with uninitialized PSET.' )
557 mpp_pset_numroots =
size(pset%root_pelist)
558 end function mpp_pset_numroots
560 subroutine mpp_pset_get_root_pelist(pset,pelist,commID)
561 type(mpp_pset_type),
intent(in) :: pset
562 integer,
intent(out) :: pelist(:)
563 integer,
intent(out),
optional :: commID
565 if( .NOT.pset%initialized )
call mpp_error( fatal, &
566 'MPP_PSET_GET_ROOT_PELIST: called with uninitialized PSET.' )
567 if(
size(pelist).NE.
size(pset%root_pelist) )
then
568 write( text,
'(a,2i6)' ) &
569 'pelist argument has wrong size: requested, actual=', &
570 size(pelist),
size(pset%root_pelist)
571 call mpp_error( fatal,
'MPP_PSET_GET_ROOT_PELIST: '//text )
573 pelist(:) = pset%root_pelist(:)
574 if(
PRESENT(commid) )
then
579 'MPP_PSET_GET_ROOT_PELIST: commID is only defined under -Duse_libMPI.' )
582 end subroutine mpp_pset_get_root_pelist
584 end module mpp_pset_mod
subroutine mpp_sync_self(pelist, check, request, msg_size, msg_type)
This is to check if current PE's outstanding puts are complete but we can't use shmem_fence because w...
integer function stdout()
This function returns the current standard fortran unit numbers for output.
subroutine mpp_set_current_pelist(pelist, no_sync)
Set context pelist.
subroutine mpp_init(flags, localcomm, test_level, alt_input_nml_path)
Initialize the mpp_mod module. Must be called before any usage.
integer function stderr()
This function returns the current standard fortran unit numbers for error messages.
subroutine mpp_declare_pelist(pelist, name, commID)
Declare a pelist.
integer function mpp_npes()
Returns processor count for current pelist.
integer function mpp_pe()
Returns processor ID.
subroutine mpp_sync(pelist, do_self)
Synchronize PEs in list.
Calculate parallel checksums.
Receive data from another PE.
Send data to a receiving PE.