FMS  2025.04
Flexible Modeling System
mpp_pset.F90
1 !***********************************************************************
2 !* Apache License 2.0
3 !*
4 !* This file is part of the GFDL Flexible Modeling System (FMS).
5 !*
6 !* Licensed under the Apache License, Version 2.0 (the "License");
7 !* you may not use this file except in compliance with the License.
8 !* You may obtain a copy of the License at
9 !*
10 !* http://www.apache.org/licenses/LICENSE-2.0
11 !*
12 !* FMS is distributed in the hope that it will be useful, but WITHOUT
13 !* WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied;
14 !* without even the implied warranty of MERCHANTABILITY or FITNESS FOR A
15 !* PARTICULAR PURPOSE. See the License for the specific language
16 !* governing permissions and limitations under the License.
17 !***********************************************************************
18 #ifdef test_mpp_pset
19 !PSET_DEBUG is always turned on in the test program
20 #define PSET_DEBUG
21 #endif
22 
23 !> @defgroup mpp_pset_mod mpp_pset_mod
24 !> @ingroup mpp
25 !> @brief Handles PSETs(Persistent Shared-memory Execution Threads) for mpp modules
26 !!
27 !! @author V. Balaji (v.balaji@noaa.gov)
28 !! @date 2006-01-15
29 
30 !> @addtogroup mpp_pset_mod mpp_pset_mod
31 !> @{
32 module mpp_pset_mod
33 #include <fms_platform.h>
34  use mpp_mod, only: mpp_pe, mpp_npes, mpp_root_pe, mpp_send, mpp_recv, &
35  mpp_sync, mpp_error, fatal, warning, stdout, stderr, mpp_chksum, &
36  mpp_declare_pelist, mpp_get_current_pelist, mpp_set_current_pelist, &
37  mpp_init, comm_tag_1, comm_tag_2, comm_tag_3, mpp_sync_self
38  implicit none
39  private
40 
41 !private variables
42  integer :: pe
43  integer :: commID !MPI communicator, copy here from pset
44  logical :: verbose=.false.
45  logical :: module_is_initialized=.false.
46  character(len=256) :: text
47 !generic interfaces
49  module procedure mpp_pset_broadcast_ptr_scalar
50  module procedure mpp_pset_broadcast_ptr_array
51  end interface
52  interface mpp_send_ptr
53  module procedure mpp_send_ptr_scalar
54  module procedure mpp_send_ptr_array
55  end interface
56  interface mpp_recv_ptr
57  module procedure mpp_recv_ptr_scalar
58  module procedure mpp_recv_ptr_array
59  end interface
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
65  end interface
66 !public type
67  type :: mpp_pset_type
68  private
69  integer :: npset !number of PSETs
70  integer :: next_in_pset, prev_in_pset !next and prev PE in PSET (cyclic)
71  integer :: root_in_pset !PE designated to be the root within PSET
72  logical :: root !true if you are the root PSET
73  integer :: pos !position of current PE within pset
74 !stack is allocated by root
75 !it is then mapped to mpp_pset_stack by mpp_pset_broadcast_ptr
76  real, allocatable :: stack(:)
77  integer, allocatable :: pelist(:) !base PElist
78  integer, allocatable :: root_pelist(:) !a PElist of all the roots
79  integer, allocatable :: pset(:) !PSET IDs
80  integer(POINTER_KIND) :: p_stack
81  integer :: lstack, maxstack, hiWM !current stack length, max, hiWM
82  integer :: commID
83  character(len=32) :: name
84  logical :: initialized=.false.
85  end type mpp_pset_type
86 !public types
87  public :: mpp_pset_type
88 !public variables
89 !public member functions
90  public :: mpp_pset_create, mpp_pset_sync, mpp_pset_broadcast, &
91  mpp_pset_broadcast_ptr, mpp_pset_check_ptr, mpp_pset_segment_array, &
92  mpp_pset_stack_push, mpp_pset_stack_reset, mpp_pset_print_chksum, &
93  mpp_pset_delete, mpp_pset_root, mpp_pset_numroots, mpp_pset_init, &
94  mpp_pset_get_root_pelist, mpp_pset_print_stack_chksum
95 
96 
97 contains
98  subroutine mpp_pset_init
99  module_is_initialized = .true.
100  end subroutine mpp_pset_init
101 
102  subroutine mpp_pset_create(npset,pset,stacksize,pelist, commID)
103 !create PSETs
104 ! called by all PEs in parent pelist
105 ! mpset must be exact divisor of npes
106  integer, intent(in) :: npset !number of PSETs per set
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
111 
112  integer :: npes, my_commID
113  integer :: i, j, k, out_unit, errunit
114  integer, allocatable :: my_pelist(:), root_pelist(:)
115 
116  call mpp_init()
117  call mpp_pset_init()
118 
119 #ifdef PSET_DEBUG
120  verbose=.true.
121 #endif
122  out_unit = stdout()
123  errunit = stderr()
124  pe = mpp_pe()
125  if(present(pelist)) then
126  npes = size(pelist(:))
127  else
128  npes = mpp_npes()
129  endif
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
134  call mpp_error( fatal, text )
135  end if
136 
137  !configure out root_pelist
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')
143  my_pelist = pelist
144  my_commid = commid
145  else
146  call mpp_get_current_pelist(my_pelist, commid = my_commid)
147  endif
148  do i = 0,npes/npset-1
149  root_pelist(i) = my_pelist(npset*i)
150  enddo
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!' )
155  pset%npset = npset
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
160 !create the root PElist
161  pset%root_pelist = root_pelist
162  allocate( pset%pset(0:npset-1) )
163  do i = 0,npes/npset-1
164  k = npset*i
165 !designate the root PE, next PE, prev PE
166  do j = 0,npset-1
167  if( pe.EQ.pset%pelist(k+j) )then
168  pset%pset(:) = pset%pelist(k:k+npset-1)
169  pset%pos = j
170  pset%root_in_pset = pset%root_pelist(i)
171  if( j.EQ.0 )then
172  pset%prev_in_pset = pset%pelist(k+npset-1)
173  else
174  pset%prev_in_pset = pset%pelist(k+j-1)
175  end if
176  if( j.EQ.npset-1 )then
177  pset%next_in_pset = pset%pelist(k)
178  else
179  pset%next_in_pset = pset%pelist(k+j+1)
180  end if
181  end if
182  end do
183  end do
184 
185  pset%root = pe.EQ.pset%root_in_pset
186 
187 !stack
188  pset%hiWM = 0 !initialize hi-water-mark
189  pset%maxstack = 1000000 !default
190  if( PRESENT(stacksize) )pset%maxstack = stacksize
191  write( out_unit,'(a,i8)' ) &
192  'MPP_PSET_CREATE: setting stacksize=', pset%maxstack
193  if( pset%root )then
194  allocate( pset%stack(pset%maxstack) )
195 #ifdef use_CRI_pointers
196  pset%p_stack = loc(pset%stack)
197 #endif
198  end if
199  pset%initialized = .true. !must be called before using pset
200  call mpp_pset_broadcast_ptr(pset,pset%p_stack)
201  endif
202 
203  call mpp_declare_pelist(root_pelist)
204 
205  if( verbose )then
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(:)
210  end if
211  end subroutine mpp_pset_create
212 
213  subroutine mpp_pset_delete(pset)
214  type(mpp_pset_type), intent(inout) :: pset
215  integer :: out_unit
216 
217  out_unit = stdout()
218  if( .NOT.pset%initialized )call mpp_error( fatal, &
219  'MPP_PSET_DELETE: called with uninitialized PSET.' )
220 !deallocate arrays...
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
227 !... and set status flag
228  pset%initialized = .false.
229  end subroutine mpp_pset_delete
230 
231  subroutine mpp_send_ptr_scalar( ptr, pe )
232  integer(POINTER_KIND), intent(in) :: ptr
233  integer, intent(in) :: pe
234 
235 !currently only wraps mpp_send
236 !on some architectures, mangling might occur
237  call mpp_send( ptr, pe, tag=comm_tag_1 )
238  end subroutine mpp_send_ptr_scalar
239 
240  subroutine mpp_send_ptr_array( ptr, pe )
241  integer(POINTER_KIND), intent(in) :: ptr(:)
242  integer, intent(in) :: pe
243 
244 !currently only wraps mpp_send
245 !on some architectures, mangling might occur
246  call mpp_send( ptr, size(ptr), pe, tag=comm_tag_2 )
247  end subroutine mpp_send_ptr_array
248 
249  subroutine mpp_recv_ptr_scalar( ptr, pe )
250  integer(POINTER_KIND), intent(inout) :: ptr
251  integer, intent(in) :: pe
252 
253  call mpp_recv( ptr, pe, tag=comm_tag_1 )
254  call mpp_translate_remote_ptr( ptr, pe )
255  return
256  end subroutine mpp_recv_ptr_scalar
257 
258  subroutine mpp_recv_ptr_array( ptr, pe )
259  integer(POINTER_KIND), intent(inout) :: ptr(:)
260  integer, intent(in) :: pe
261  integer :: i
262 
263  call mpp_recv( ptr, size(ptr), pe, tag=comm_tag_2 )
264  do i = 1, size(ptr)
265  call mpp_translate_remote_ptr( ptr(i), pe )
266  end do
267  return
268  end subroutine mpp_recv_ptr_array
269 
270  subroutine mpp_translate_remote_ptr( ptr, pe )
271 !modifies the received pointer to correct numerical address
272  integer(POINTER_KIND), intent(inout) :: ptr
273  integer, intent(in) :: pe
274  return
275  end subroutine mpp_translate_remote_ptr
276 
277  subroutine mpp_pset_sync(pset)
278 !this is a replacement for mpp_sync, doing syncs across
279 !shared arrays without calling mpp_sync
280  type(mpp_pset_type), intent(in) :: pset
281 
282  if( .NOT.pset%initialized )call mpp_error( fatal, &
283  'MPP_PSET_SYNC: called with uninitialized PSET.' )
284 !currently does mpp_sync!!! slow!!!
285 !try and make a lightweight pset sync
286  call mpp_sync
287  end subroutine mpp_pset_sync
288 
289  subroutine mpp_pset_broadcast(pset,a)
290 !broadcast value on the root to its sub-threads
291  type(mpp_pset_type), intent(in) :: pset
292  real, intent(inout) :: a
293  integer :: i
294 
295  if( .NOT.pset%initialized )call mpp_error( fatal, &
296  'MPP_PSET_BROADCAST: called with uninitialized PSET.' )
297  if( pset%root )then
298  do i = 1,pset%npset-1
299  call mpp_send( a, pset%pset(i), tag=comm_tag_3 )
300  end do
301  else
302  call mpp_recv( a, pset%root_in_pset, tag=comm_tag_3 )
303  end if
304  call mpp_pset_sync(pset)
305  end subroutine mpp_pset_broadcast
306 
307  subroutine mpp_pset_broadcast_ptr_scalar(pset,ptr)
308 !create a shared array by broadcasting pointer
309 !root allocates memory and passes pointer in
310 !on return all other PSETs will have the pointer to a shared object
311  type(mpp_pset_type), intent(in) :: pset
312  integer(POINTER_KIND), intent(inout) :: ptr
313  integer :: i
314 
315  if( .NOT.pset%initialized )call mpp_error( fatal, &
316  'MPP_PSET_BROADCAST_PTR: called with uninitialized PSET.' )
317  commid = pset%commID !pass to mpp_translate_remote_ptr
318  if( pset%root )then
319  do i = 1,pset%npset-1
320  call mpp_send_ptr( ptr, pset%pset(i) )
321  end do
322  else
323  call mpp_recv_ptr( ptr, pset%root_in_pset )
324  end if
325  call mpp_sync_self()
326  end subroutine mpp_pset_broadcast_ptr_scalar
327 
328  subroutine mpp_pset_broadcast_ptr_array(pset,ptr)
329 !create a shared array by broadcasting pointer
330 !root allocates memory and passes pointer in
331 !on return all other PSETs will have the pointer to a shared object
332  type(mpp_pset_type), intent(in) :: pset
333  integer(POINTER_KIND), intent(inout) :: ptr(:)
334  integer :: i
335 
336  if( .NOT.pset%initialized )call mpp_error( fatal, &
337  'MPP_PSET_BROADCAST_PTR: called with uninitialized PSET.' )
338  commid = pset%commID !pass to mpp_translate_remote_ptr
339  if( pset%root )then
340  do i = 1,pset%npset-1
341  call mpp_send_ptr( ptr, pset%pset(i) )
342  end do
343  else
344  call mpp_recv_ptr( ptr, pset%root_in_pset )
345  end if
346  call mpp_sync_self()
347 
348  end subroutine mpp_pset_broadcast_ptr_array
349 
350  subroutine mpp_pset_check_ptr(pset,ptr)
351 !checks if the supplied pointer is indeed shared
352  type(mpp_pset_type), intent(in) :: pset
353 #ifdef use_CRI_pointers
354  real :: dummy
355  pointer( ptr, dummy )
356 #else
357  integer(POINTER_KIND), intent(in) :: ptr
358 #endif
359 #ifdef PSET_DEBUG
360  integer(POINTER_KIND) :: p
361  integer :: i
362  if( .NOT.pset%initialized )call mpp_error( fatal, &
363  'MPP_PSET_CHECK_PTR: called with uninitialized PSET.' )
364  commid = pset%commID !pass to mpp_translate_remote_ptr
365 !check if this is a shared pointer
366  p = ptr
367  if( pset%root )then
368  do i = 1,pset%npset-1
369  call mpp_send_ptr( p, pset%pset(i) )
370  end do
371  else
372  call mpp_recv_ptr( p, pset%root_in_pset )
373  end if
374  call mpp_pset_sync(pset)
375  if( p.NE.ptr )call mpp_error( fatal, &
376  'MPP_PSET_CHECK_PTR: pointers do not match!' )
377 #else
378 !do nothing if the debug CPP flag isn't on
379 #endif
380  end subroutine mpp_pset_check_ptr
381 
382  subroutine mpp_pset_segment_array( pset, ls, le, lsp, lep )
383 !given input indices ls, le, returns indices lsp, lep
384 !so that segments span the range ls:le with no overlaps.
385 !attempts load balance: also some PSETs might get lsp>lep
386 !so that do-loops will be null
387  type(mpp_pset_type), intent(in) :: pset
388  integer, intent(in) :: ls, le
389  integer, intent(out) :: lsp, lep
390  integer :: i
391 
392  if( .NOT.pset%initialized )call mpp_error( fatal, &
393  'MPP_PSET_SEGMENT_ARRAY: called with uninitialized PSET.' )
394 #ifdef PSET_DEBUG
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
399  call mpp_error( warning, text )
400  end if
401 #endif
402  lep = ls-1 !initialize so that lsp is correct on first pass
403  do i = 0,pset%pos
404  lsp = lep + 1
405  lep = lsp + ceiling( real(le-lsp+1)/(pset%npset-i) ) - 1
406  end do
407  end subroutine mpp_pset_segment_array
408 
409  subroutine mpp_pset_stack_push( pset, ptr, len )
410 !mpp_malloc specialized for shared arrays
411 !len is the length of the required array
412 !lstack is the stack already in play
413 !user should zero lstack (call mpp_pset_stack_reset) when the stack is to be cleared
414  type(mpp_pset_type), intent(inout) :: pset
415  integer, intent(in) :: len
416 #ifdef use_CRI_pointers
417  real :: dummy
418  pointer( ptr, dummy )
419  real :: stack(pset%maxstack)
420  pointer( p, stack )
421 
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
429  call mpp_error( fatal, text )
430  end if
431  p = pset%p_stack !point stack to shared stack pointer
432  ptr = loc( stack(pset%lstack+1) )
433  call mpp_pset_check_ptr(pset,ptr) !make sure ptr is the same across PSETs
434  pset%lstack = pset%lstack + len
435  pset%hiWM = max( pset%hiWM, pset%lstack )
436 #else
437  integer(POINTER_KIND), intent(out) :: ptr
438  call mpp_error( fatal, &
439  'MPP_PSET_STACK_PUSH only works with Cray pointers.' )
440 #endif
441  end subroutine mpp_pset_stack_push
442 
443  subroutine mpp_pset_stack_reset(pset)
444  type(mpp_pset_type), intent(inout) :: pset
445 !reset stack... will reuse any temporary arrays! USE WITH CARE
446 !next few lines are to zero stack contents...
447 !but it's better noone tries to use uninitialized stack variables!
448 ! integer :: l1, l2
449 ! real :: mpp_pset_stack(maxstack)
450 ! pointer( p_mpp_pset_stack, mpp_pset_stack )
451 ! p_mpp_pset_stack = ptr_mpp_pset_stack
452 ! call mpp_pset_array_segment( 1, lstack, l1, l2 )
453 ! mpp_pset_stack(l1:l2) = 0.
454  if( .NOT.pset%initialized )call mpp_error( fatal, &
455  'MPP_PSET_STACK_RESET: called with uninitialized PSET.' )
456  pset%lstack = 0
457  end subroutine mpp_pset_stack_reset
458 
459  subroutine mpp_pset_print_chksum_1d(pset, caller, array)
460 !print a checksum of an array
461 !pass the whole domain seen by root PSET
462 !add lines to check on shared array?
463  type(mpp_pset_type), intent(in) :: pset
464  character(len=*), intent(in) :: caller
465  real, intent(in) :: array(:)
466 
467 #ifdef PSET_DEBUG
468  integer :: errunit
469  logical :: do_print
470  integer(LONG_KIND) :: chksum
471 
472  if( .NOT.pset%initialized )call mpp_error( fatal, &
473  'MPP_PSET_PRINT_CHKSUM: called with uninitialized PSET.' )
474  errunit = stderr()
475 
476  if( pset%root )then
477  do_print = pe.EQ.mpp_root_pe() !set to T to print from all PEs
478  call mpp_set_current_pelist(pset%root_pelist)
479  chksum = mpp_chksum( array )
480  if( do_print ) &
481  write( errunit, '(a,z18)' )trim(caller)//' chksum=', chksum
482  end if
483  call mpp_set_current_pelist(pset%pelist)
484 #endif
485  return
486  end subroutine mpp_pset_print_chksum_1d
487 
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 )
495  p = loc(array)
496 #else
497  array1d = transfer( array, array1d )
498 #endif
499  call mpp_pset_print_chksum(pset, caller, array1d)
500  end subroutine mpp_pset_print_chksum_2d
501 
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(:,:,:) !overload for other ranks
506  real :: array1D( size(array) )
507 #ifdef use_CRI_pointers
508  pointer( p, array1d )
509  p = loc(array)
510 #else
511  array1d = transfer( array, array1d )
512 #endif
513  call mpp_pset_print_chksum(pset, caller, array1d)
514  end subroutine mpp_pset_print_chksum_3d
515 
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 )
523  p = loc(array)
524 #else
525  array1d = transfer( array, array1d )
526 #endif
527  call mpp_pset_print_chksum(pset, caller, array1d)
528  end subroutine mpp_pset_print_chksum_4d
529 
530  subroutine mpp_pset_print_stack_chksum( pset, caller )
531  type(mpp_pset_type), intent(in) :: pset
532  character(len=*), intent(in) :: caller
533 
534  if( .NOT.pset%initialized )call mpp_error( fatal, &
535  'MPP_PSET_PRINT_STACK_CHKSUM: called with uninitialized PSET.' )
536  call mpp_pset_print_chksum( pset, trim(caller)//' stack', &
537  pset%stack(1:pset%lstack) )
538  end subroutine mpp_pset_print_stack_chksum
539 
540 !accessor functions
541  function mpp_pset_root(pset)
542  logical :: mpp_pset_root
543  type(mpp_pset_type), intent(in) :: pset
544 
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
549 
550  function mpp_pset_numroots(pset)
551 !necessary to export root_pelist: caller needs to pre-allocate
552  integer :: mpp_pset_numroots
553  type(mpp_pset_type), intent(in) :: pset
554 
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
559 
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
564 
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 )
572  end if
573  pelist(:) = pset%root_pelist(:)
574  if( PRESENT(commid) )then
575 #ifdef use_libMPI
576  commid = pset%commID
577 #else
578  call mpp_error( warning, &
579  'MPP_PSET_GET_ROOT_PELIST: commID is only defined under -Duse_libMPI.' )
580 #endif
581  end if
582  end subroutine mpp_pset_get_root_pelist
583 
584 end module mpp_pset_mod
585 !> @}
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.
Definition: mpp_util.inc:42
subroutine mpp_set_current_pelist(pelist, no_sync)
Set context pelist.
Definition: mpp_util.inc:498
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.
Definition: mpp_util.inc:50
subroutine mpp_declare_pelist(pelist, name, commID)
Declare a pelist.
Definition: mpp_util.inc:469
integer function mpp_npes()
Returns processor count for current pelist.
Definition: mpp_util.inc:420
integer function mpp_pe()
Returns processor ID.
Definition: mpp_util.inc:406
subroutine mpp_sync(pelist, do_self)
Synchronize PEs in list.
Calculate parallel checksums.
Definition: mpp.F90:1213
Error handler.
Definition: mpp.F90:381
Receive data from another PE.
Definition: mpp.F90:950
Send data to a receiving PE.
Definition: mpp.F90:1017