diff --git a/CHANGELOG.md b/CHANGELOG.md index c7aafc6484..c2d0e0a4d0 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -6,6 +6,12 @@ and this project uses `yyyy.rr[.pp]`, where `yyyy` is the year a patch is releas `rr` is a sequential release number (starting from `01`), and an optional two-digit sequential patch number (starting from `01`). +## [2021.02.01] - 2021-07-09 +### Fixed +- MPP: Fixed uninitialized variables for data domains in mpp domain broadcast routines +### Added +- MPP: Restored deleted pset functionality needed by GFDL SCM by reinstating mpp_peset.F90 + ## [2021.02] - 2021-05-20 ### Added - FMS2_IO: Added fms2_io support for boundary condition restarts. `register_restart_region_2d` and `register_restart_region_3d` were added to fms2_io’s `register_restart_field` interface and `read_restart_bc` and `write_restart_bc` subroutines were added to read and write boundary conditions restarts. See [test_fms/fms2_io/test_bc_restart.F90](https://github.com/NOAA-GFDL/FMS/blob/9d55115a331685e4c6e01f2dfb3b770a9f80fa37/test_fms/fms2_io/test_bc_restart.F90) for sample usage. diff --git a/configure.ac b/configure.ac index 6e0b63d1e8..7963297efe 100644 --- a/configure.ac +++ b/configure.ac @@ -25,7 +25,7 @@ AC_PREREQ([2.69]) # Initialize with name, version, and support email address. AC_INIT([GFDL FMS Library], - [2021.02.0], + [2021.02.01], [gfdl.climate.model.info@noaa.gov], [FMS], [https://www.gfdl.noaa.gov/fms]) diff --git a/mpp/Makefile.am b/mpp/Makefile.am index 7f1987b516..150127c315 100644 --- a/mpp/Makefile.am +++ b/mpp/Makefile.am @@ -40,6 +40,7 @@ libmpp_la_SOURCES = \ mpp_efp.F90 \ mpp_memutils.F90 \ mpp_memuse.c \ + mpp_pset.F90 \ include/group_update_pack.inc \ include/group_update_unpack.inc \ include/mpp_alltoall_mpi.h \ @@ -171,7 +172,7 @@ mpp_utilities_mod.$(FC_MODEXT): \ mpp_domains_mod.$(FC_MODEXT): \ mpp_data_mod.$(FC_MODEXT) \ mpp_parameter_mod.$(FC_MODEXT) \ - mpp_mod.$(FC_MODEXT) \ + mpp_mod.$(FC_MODEXT) \ mpp_memutils_mod.$(FC_MODEXT) \ mpp_efp_mod.$(FC_MODEXT) \ include/mpp_define_nest_domains.inc \ @@ -233,6 +234,7 @@ mpp_io_mod.$(FC_MODEXT): \ include/mpp_io_unstructured_read.inc mpp_efp_mod.$(FC_MODEXT): mpp_parameter_mod.$(FC_MODEXT) mpp_mod.$(FC_MODEXT) mpp_memutils_mod.$(FC_MODEXT): mpp_mod.$(FC_MODEXT) +mpp_pset_mod.$(FC_MODEXT): mpp_mod.$(FC_MODEXT) # Mod files are built and then installed as headers. BUILT_SOURCES = \ @@ -243,7 +245,8 @@ BUILT_SOURCES = \ mpp_memutils_mod.$(FC_MODEXT) \ mpp_efp_mod.$(FC_MODEXT) \ mpp_domains_mod.$(FC_MODEXT) \ - mpp_io_mod.$(FC_MODEXT) -nodist_include_HEADERS = $(BUILD_SOURCES) + mpp_io_mod.$(FC_MODEXT) \ + mpp_pset_mod.$(FC_MODEXT) +nodist_include_HEADERS = $(BUILT_SOURCES) include $(top_srcdir)/mkmods.mk diff --git a/mpp/include/mpp_domains_misc.inc b/mpp/include/mpp_domains_misc.inc index 69899a3b29..331896384d 100644 --- a/mpp/include/mpp_domains_misc.inc +++ b/mpp/include/mpp_domains_misc.inc @@ -536,6 +536,10 @@ end subroutine init_nonblock_type domain%x%compute%end = -1 domain%y%compute%begin = 1 domain%y%compute%end = -1 + domain%x%data %begin = -1 + domain%x%data %end = -1 + domain%y%data %begin = -1 + domain%y%data %end = -1 domain%x%global %begin = -1 domain%x%global %end = -1 domain%y%global %begin = -1 @@ -577,6 +581,10 @@ end subroutine init_nonblock_type domain%list(listpos)%y%compute%end = msg(5) domain%list(listpos)%tile_id(1) = msg(6) if(domain%x(1)%global%begin < 0) then + domain%x(1)%data %begin = msg(2) + domain%x(1)%data %end = msg(3) + domain%y(1)%data %begin = msg(4) + domain%y(1)%data %end = msg(5) domain%x(1)%global%begin = msg(2) domain%x(1)%global%end = msg(3) domain%y(1)%global%begin = msg(4) @@ -592,6 +600,10 @@ end subroutine init_nonblock_type endif domain%ntiles = msg(12) else + domain%x(1)%data %begin = msg(2) - msg(7) + domain%x(1)%data %end = msg(3) + msg(8) + domain%y(1)%data %begin = msg(4) - msg(9) + domain%y(1)%data %end = msg(5) + msg(10) domain%x(1)%global%begin = min(domain%x(1)%global%begin, msg(2)) domain%x(1)%global%end = max(domain%x(1)%global%end, msg(3)) domain%y(1)%global%begin = min(domain%y(1)%global%begin, msg(4)) @@ -651,6 +663,10 @@ end subroutine init_nonblock_type domain_out%x%compute%end = -1 domain_out%y%compute%begin = 1 domain_out%y%compute%end = -1 + domain_out%x%data %begin = -1 + domain_out%x%data %end = -1 + domain_out%y%data %begin = -1 + domain_out%y%data %end = -1 domain_out%x%global %begin = -1 domain_out%x%global %end = -1 domain_out%y%global %begin = -1 @@ -698,6 +714,10 @@ end subroutine init_nonblock_type domain_out%list(listpos)%y%compute%end = msg(5) domain_out%list(listpos)%tile_id(1) = msg(6) if(domain_out%x(1)%global%begin < 0) then + domain_out%x(1)%data %begin = msg(2) + domain_out%x(1)%data %end = msg(3) + domain_out%y(1)%data %begin = msg(4) + domain_out%y(1)%data %end = msg(5) domain_out%x(1)%global%begin = msg(2) domain_out%x(1)%global%end = msg(3) domain_out%y(1)%global%begin = msg(4) @@ -713,6 +733,10 @@ end subroutine init_nonblock_type endif domain_out%ntiles = msg(12) else + domain_out%x(1)%data %begin = msg(2) - msg(7) + domain_out%x(1)%data %end = msg(3) + msg(8) + domain_out%y(1)%data %begin = msg(4) - msg(9) + domain_out%y(1)%data %end = msg(5) + msg(10) domain_out%x(1)%global%begin = min(domain_out%x(1)%global%begin, msg(2)) domain_out%x(1)%global%end = max(domain_out%x(1)%global%end, msg(3)) domain_out%y(1)%global%begin = min(domain_out%y(1)%global%begin, msg(4)) @@ -776,6 +800,10 @@ end subroutine init_nonblock_type domain%x%compute%end = -1 domain%y%compute%begin = 0 domain%y%compute%end = -1 + domain%x%data %begin = 0 + domain%x%data %end = -1 + domain%y%data %begin = 0 + domain%y%data %end = -1 domain%x%global %begin = 0 domain%x%global %end = -1 domain%y%global %begin = 0 @@ -827,10 +855,10 @@ end subroutine init_nonblock_type domain%list(listpos)%y%compute%begin = msg(4) domain%list(listpos)%y%compute%end = msg(5) domain%list(listpos)%tile_id(1) = msg(6) - domain%list(listpos)%x%global%begin = msg(12) - domain%list(listpos)%x%global%end = msg(13) - domain%list(listpos)%y%global%begin = msg(14) - domain%list(listpos)%y%global%end = msg(15) + domain%list(listpos)%x%global %begin = msg(12) + domain%list(listpos)%x%global %end = msg(13) + domain%list(listpos)%y%global %begin = msg(14) + domain%list(listpos)%y%global %end = msg(15) listpos = listpos + 1 if( debug )write( errunit,* )'PE ', mpp_pe(), 'received domain from PE ', msg(1), 'is,ie,js,je=', msg(2:5) end if @@ -889,6 +917,10 @@ end subroutine init_nonblock_type domain%x%compute%end = -1 domain%y%compute%begin = 0 domain%y%compute%end = -1 + domain%x%data %begin = 0 + domain%x%data %end = -1 + domain%y%data %begin = 0 + domain%y%data %end = -1 domain%x%global %begin = 0 domain%x%global %end = -1 domain%y%global %begin = 0 @@ -927,6 +959,10 @@ end subroutine init_nonblock_type if( .NOT.native .AND. msg(1).NE.NULL_PE .AND. tile_coarse==msg(16) )then domain%list(listpos)%pe = msg(1) if(domain%x(1)%compute%begin == 0) then + domain%x(1)%data %begin = msg(2) - msg(7) + domain%x(1)%data %end = msg(3) + msg(8) + domain%y(1)%data %begin = msg(4) - msg(9) + domain%y(1)%data %end = msg(5) + msg(10) domain%x(1)%global%begin = msg(12) domain%x(1)%global%end = msg(13) domain%y(1)%global%begin = msg(14) @@ -947,10 +983,10 @@ end subroutine init_nonblock_type domain%list(listpos)%y%compute%begin = msg(4) domain%list(listpos)%y%compute%end = msg(5) domain%list(listpos)%tile_id(1) = msg(6) - domain%list(listpos)%x%global%begin = msg(12) - domain%list(listpos)%x%global%end = msg(13) - domain%list(listpos)%y%global%begin = msg(14) - domain%list(listpos)%y%global%end = msg(15) + domain%list(listpos)%x%global %begin = msg(12) + domain%list(listpos)%x%global %end = msg(13) + domain%list(listpos)%y%global %begin = msg(14) + domain%list(listpos)%y%global %end = msg(15) listpos = listpos + 1 if( debug )write( errunit,* )'PE ', mpp_pe(), 'received domain from PE ', msg(1), 'is,ie,js,je=', msg(2:5) end if diff --git a/mpp/mpp_pset.F90 b/mpp/mpp_pset.F90 new file mode 100644 index 0000000000..a23d75d19c --- /dev/null +++ b/mpp/mpp_pset.F90 @@ -0,0 +1,581 @@ +!*********************************************************************** +!* GNU Lesser General Public License +!* +!* This file is part of the GFDL Flexible Modeling System (FMS). +!* +!* FMS is free software: you can redistribute it and/or modify it under +!* the terms of the GNU Lesser General Public License as published by +!* the Free Software Foundation, either version 3 of the License, or (at +!* your option) any later version. +!* +!* FMS is distributed in the hope that it will be useful, but WITHOUT +!* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or +!* FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License +!* for more details. +!* +!* You should have received a copy of the GNU Lesser General Public +!* License along with FMS. If not, see . +!*********************************************************************** +! module within MPP for handling PSETs: +! PSET: Persistent Shared-memory Execution Thread +! +! AUTHOR: V. Balaji (v.balaji@noaa.gov) +! DATE: 2006-01-15 +#ifdef test_mpp_pset +!PSET_DEBUG is always turned on in the test program +#define PSET_DEBUG +#endif + +module mpp_pset_mod +#include + use mpp_mod, only: mpp_pe, mpp_npes, mpp_root_pe, mpp_send, mpp_recv, & + mpp_sync, mpp_error, FATAL, WARNING, stdout, stderr, mpp_chksum, & + mpp_declare_pelist, mpp_get_current_pelist, mpp_set_current_pelist, & + mpp_init, COMM_TAG_1, COMM_TAG_2, COMM_TAG_3, mpp_sync_self + implicit none + private + +!private variables + integer :: pe + integer :: commID !MPI communicator, copy here from pset + logical :: verbose=.FALSE. + logical :: module_is_initialized=.FALSE. + character(len=256) :: text +!generic interfaces + interface mpp_pset_broadcast_ptr + module procedure mpp_pset_broadcast_ptr_scalar + module procedure mpp_pset_broadcast_ptr_array + end interface + interface mpp_send_ptr + module procedure mpp_send_ptr_scalar + module procedure mpp_send_ptr_array + end interface + interface mpp_recv_ptr + module procedure mpp_recv_ptr_scalar + module procedure mpp_recv_ptr_array + end interface + interface mpp_pset_print_chksum + module procedure mpp_pset_print_chksum_1D + module procedure mpp_pset_print_chksum_2D + module procedure mpp_pset_print_chksum_3D + module procedure mpp_pset_print_chksum_4D + end interface +!public type + type :: mpp_pset_type + private + integer :: npset !number of PSETs + integer :: next_in_pset, prev_in_pset !next and prev PE in PSET (cyclic) + integer :: root_in_pset !PE designated to be the root within PSET + logical :: root !true if you are the root PSET + integer :: pos !position of current PE within pset +!stack is allocated by root +!it is then mapped to mpp_pset_stack by mpp_pset_broadcast_ptr + real, allocatable :: stack(:) + integer, allocatable :: pelist(:) !base PElist + integer, allocatable :: root_pelist(:) !a PElist of all the roots + integer, allocatable :: pset(:) !PSET IDs + integer(POINTER_KIND) :: p_stack + integer :: lstack, maxstack, hiWM !current stack length, max, hiWM + integer :: commID + character(len=32) :: name + logical :: initialized=.FALSE. + end type mpp_pset_type +!public types + public :: mpp_pset_type +!public variables +!public member functions + public :: mpp_pset_create, mpp_pset_sync, mpp_pset_broadcast, & + mpp_pset_broadcast_ptr, mpp_pset_check_ptr, mpp_pset_segment_array, & + mpp_pset_stack_push, mpp_pset_stack_reset, mpp_pset_print_chksum, & + mpp_pset_delete, mpp_pset_root, mpp_pset_numroots, mpp_pset_init, & + mpp_pset_get_root_pelist, mpp_pset_print_stack_chksum + + +contains + subroutine mpp_pset_init + module_is_initialized = .TRUE. + end subroutine mpp_pset_init + + subroutine mpp_pset_create(npset,pset,stacksize,pelist, commID) +!create PSETs +! called by all PEs in parent pelist +! mpset must be exact divisor of npes + integer, intent(in) :: npset !number of PSETs per set + type(mpp_pset_type), intent(inout) :: pset + integer, intent(in), optional :: stacksize + integer, intent(in), optional :: pelist(:) + integer, intent(in), optional :: commID + + integer :: npes, my_commID + integer :: i, j, k, out_unit, errunit + integer, allocatable :: my_pelist(:), root_pelist(:) + + call mpp_init() + call mpp_pset_init() + +#ifdef PSET_DEBUG + verbose=.TRUE. +#endif + out_unit = stdout() + errunit = stderr() + pe = mpp_pe() + if(present(pelist)) then + npes = size(pelist(:)) + else + npes = mpp_npes() + endif + if( mod(npes,npset).NE.0 )then + write( text,'(a,2i6)' ) & + 'MPP_PSET_CREATE: PSET size (npset) must divide npes exactly:'// & + ' npset, npes=', npset, npes + call mpp_error( FATAL, text ) + end if + + !configure out root_pelist + allocate(my_pelist(0:npes-1) ) + allocate(root_pelist(0:npes/npset-1) ) + if(present(pelist)) then + if(.not. present(commID)) call mpp_error(FATAL, & + 'MPP_PSET_CREATE: when pelist is present, commID should also be present') + my_pelist = pelist + my_commID = commID + else + call mpp_get_current_pelist(my_pelist, commID = my_commID) + endif + do i = 0,npes/npset-1 + root_pelist(i) = my_pelist(npset*i) + enddo + write( out_unit,'(a,i6)' )'MPP_PSET_CREATE creating PSETs... npset=', npset + if(ANY(my_pelist == pe) ) then + if( pset%initialized )call mpp_error( FATAL, & + 'MPP_PSET_CREATE: PSET already initialized!' ) + pset%npset = npset + allocate( pset%pelist(0:npes-1) ) + allocate( pset%root_pelist(0:npes/npset-1) ) + pset%commID = my_commID + pset%pelist = my_pelist +!create the root PElist + pset%root_pelist = root_pelist + allocate( pset%pset(0:npset-1) ) + do i = 0,npes/npset-1 + k = npset*i +!designate the root PE, next PE, prev PE + do j = 0,npset-1 + if( pe.EQ.pset%pelist(k+j) )then + pset%pset(:) = pset%pelist(k:k+npset-1) + pset%pos = j + pset%root_in_pset = pset%root_pelist(i) + if( j.EQ.0 )then + pset%prev_in_pset = pset%pelist(k+npset-1) + else + pset%prev_in_pset = pset%pelist(k+j-1) + end if + if( j.EQ.npset-1 )then + pset%next_in_pset = pset%pelist(k) + else + pset%next_in_pset = pset%pelist(k+j+1) + end if + end if + end do + end do + + pset%root = pe.EQ.pset%root_in_pset + +!stack + pset%hiWM = 0 !initialize hi-water-mark + pset%maxstack = 1000000 !default + if( PRESENT(stacksize) )pset%maxstack = stacksize + write( out_unit,'(a,i8)' ) & + 'MPP_PSET_CREATE: setting stacksize=', pset%maxstack + if( pset%root )then + allocate( pset%stack(pset%maxstack) ) +#ifdef use_CRI_pointers + pset%p_stack = LOC(pset%stack) +#endif + end if + pset%initialized = .TRUE. !must be called before using pset + call mpp_pset_broadcast_ptr(pset,pset%p_stack) + endif + + call mpp_declare_pelist(root_pelist) + + if( verbose )then + write( errunit,'(a,4i6)' )'MPP_PSET_CREATE: pe, root, next, prev=', & + pe, pset%root_in_pset, pset%next_in_pset, pset%prev_in_pset + write( errunit,* )'PE ', pe, ' pset=', pset%pset(:) + write( out_unit,* )'root pelist=', pset%root_pelist(:) + end if + end subroutine mpp_pset_create + + subroutine mpp_pset_delete(pset) + type(mpp_pset_type), intent(inout) :: pset + integer :: out_unit + + out_unit = stdout() + if( .NOT.pset%initialized )call mpp_error( FATAL, & + 'MPP_PSET_DELETE: called with uninitialized PSET.' ) +!deallocate arrays... + deallocate( pset%pelist ) + deallocate( pset%root_pelist ) + deallocate( pset%pset ) + if( pset%root )deallocate( pset%stack ) + write( out_unit, '(a,i10)' ) & + 'Deleting PSETs... stack high-water-mark=', pset%hiWM +!... and set status flag + pset%initialized = .FALSE. + end subroutine mpp_pset_delete + + subroutine mpp_send_ptr_scalar( ptr, pe ) + integer(POINTER_KIND), intent(in) :: ptr + integer, intent(in) :: pe + +!currently only wraps mpp_send +!on some architectures, mangling might occur + call mpp_send( ptr, pe, tag=COMM_TAG_1 ) + end subroutine mpp_send_ptr_scalar + + subroutine mpp_send_ptr_array( ptr, pe ) + integer(POINTER_KIND), intent(in) :: ptr(:) + integer, intent(in) :: pe + +!currently only wraps mpp_send +!on some architectures, mangling might occur + call mpp_send( ptr, size(ptr), pe, tag=COMM_TAG_2 ) + end subroutine mpp_send_ptr_array + + subroutine mpp_recv_ptr_scalar( ptr, pe ) + integer(POINTER_KIND), intent(inout) :: ptr + integer, intent(in) :: pe + + call mpp_recv( ptr, pe, tag=COMM_TAG_1 ) + call mpp_translate_remote_ptr( ptr, pe ) + return + end subroutine mpp_recv_ptr_scalar + + subroutine mpp_recv_ptr_array( ptr, pe ) + integer(POINTER_KIND), intent(inout) :: ptr(:) + integer, intent(in) :: pe + integer :: i + + call mpp_recv( ptr, size(ptr), pe, tag=COMM_TAG_2 ) + do i = 1, size(ptr) + call mpp_translate_remote_ptr( ptr(i), pe ) + end do + return + end subroutine mpp_recv_ptr_array + + subroutine mpp_translate_remote_ptr( ptr, pe ) +!modifies the received pointer to correct numerical address + integer(POINTER_KIND), intent(inout) :: ptr + integer, intent(in) :: pe + return + end subroutine mpp_translate_remote_ptr + + subroutine mpp_pset_sync(pset) +!this is a replacement for mpp_sync, doing syncs across +!shared arrays without calling mpp_sync + type(mpp_pset_type), intent(in) :: pset + + if( .NOT.pset%initialized )call mpp_error( FATAL, & + 'MPP_PSET_SYNC: called with uninitialized PSET.' ) +!currently does mpp_sync!!! slow!!! +!try and make a lightweight pset sync + call mpp_sync + end subroutine mpp_pset_sync + + subroutine mpp_pset_broadcast(pset,a) +!broadcast value on the root to its sub-threads + type(mpp_pset_type), intent(in) :: pset + real, intent(inout) :: a + integer :: i + + if( .NOT.pset%initialized )call mpp_error( FATAL, & + 'MPP_PSET_BROADCAST: called with uninitialized PSET.' ) + if( pset%root )then + do i = 1,pset%npset-1 + call mpp_send( a, pset%pset(i), tag=COMM_TAG_3 ) + end do + else + call mpp_recv( a, pset%root_in_pset, tag=COMM_TAG_3 ) + end if + call mpp_pset_sync(pset) + end subroutine mpp_pset_broadcast + + subroutine mpp_pset_broadcast_ptr_scalar(pset,ptr) +!create a shared array by broadcasting pointer +!root allocates memory and passes pointer in +!on return all other PSETs will have the pointer to a shared object + type(mpp_pset_type), intent(in) :: pset + integer(POINTER_KIND), intent(inout) :: ptr + integer :: i + + if( .NOT.pset%initialized )call mpp_error( FATAL, & + 'MPP_PSET_BROADCAST_PTR: called with uninitialized PSET.' ) + commID = pset%commID !pass to mpp_translate_remote_ptr + if( pset%root )then + do i = 1,pset%npset-1 + call mpp_send_ptr( ptr, pset%pset(i) ) + end do + else + call mpp_recv_ptr( ptr, pset%root_in_pset ) + end if + call mpp_sync_self() + end subroutine mpp_pset_broadcast_ptr_scalar + + subroutine mpp_pset_broadcast_ptr_array(pset,ptr) +!create a shared array by broadcasting pointer +!root allocates memory and passes pointer in +!on return all other PSETs will have the pointer to a shared object + type(mpp_pset_type), intent(in) :: pset + integer(POINTER_KIND), intent(inout) :: ptr(:) + integer :: i + + if( .NOT.pset%initialized )call mpp_error( FATAL, & + 'MPP_PSET_BROADCAST_PTR: called with uninitialized PSET.' ) + commID = pset%commID !pass to mpp_translate_remote_ptr + if( pset%root )then + do i = 1,pset%npset-1 + call mpp_send_ptr( ptr, pset%pset(i) ) + end do + else + call mpp_recv_ptr( ptr, pset%root_in_pset ) + end if + call mpp_sync_self() + + end subroutine mpp_pset_broadcast_ptr_array + + subroutine mpp_pset_check_ptr(pset,ptr) +!checks if the supplied pointer is indeed shared + type(mpp_pset_type), intent(in) :: pset +#ifdef use_CRI_pointers + real :: dummy + pointer( ptr, dummy ) +#else + integer(POINTER_KIND), intent(in) :: ptr +#endif +#ifdef PSET_DEBUG + integer(POINTER_KIND) :: p + integer :: i + if( .NOT.pset%initialized )call mpp_error( FATAL, & + 'MPP_PSET_CHECK_PTR: called with uninitialized PSET.' ) + commID = pset%commID !pass to mpp_translate_remote_ptr +!check if this is a shared pointer + p = ptr + if( pset%root )then + do i = 1,pset%npset-1 + call mpp_send_ptr( p, pset%pset(i) ) + end do + else + call mpp_recv_ptr( p, pset%root_in_pset ) + end if + call mpp_pset_sync(pset) + if( p.NE.ptr )call mpp_error( FATAL, & + 'MPP_PSET_CHECK_PTR: pointers do not match!' ) +#else +!do nothing if the debug CPP flag isn't on +#endif + end subroutine mpp_pset_check_ptr + + subroutine mpp_pset_segment_array( pset, ls, le, lsp, lep ) +!given input indices ls, le, returns indices lsp, lep +!so that segments span the range ls:le with no overlaps. +!attempts load balance: also some PSETs might get lsp>lep +!so that do-loops will be null + type(mpp_pset_type), intent(in) :: pset + integer, intent(in) :: ls, le + integer, intent(out) :: lsp, lep + integer :: i + + if( .NOT.pset%initialized )call mpp_error( FATAL, & + 'MPP_PSET_SEGMENT_ARRAY: called with uninitialized PSET.' ) +#ifdef PSET_DEBUG + if( le-ls+1.LT.pset%npset )then + write( text,'(3(a,i6))' ) & + 'MPP_PSET_ARRAY_SEGMENT: parallel range (', ls, ',', le, & + ') is smaller than the number of threads:', pset%npset + call mpp_error( WARNING, text ) + end if +#endif + lep = ls-1 !initialize so that lsp is correct on first pass + do i = 0,pset%pos + lsp = lep + 1 + lep = lsp + CEILING( REAL(le-lsp+1)/(pset%npset-i) ) - 1 + end do + end subroutine mpp_pset_segment_array + + subroutine mpp_pset_stack_push( pset, ptr, len ) +!mpp_malloc specialized for shared arrays +!len is the length of the required array +!lstack is the stack already in play +!user should zero lstack (call mpp_pset_stack_reset) when the stack is to be cleared + type(mpp_pset_type), intent(inout) :: pset + integer, intent(in) :: len +#ifdef use_CRI_pointers + real :: dummy + pointer( ptr, dummy ) + real :: stack(pset%maxstack) + pointer( p, stack ) + + if( .NOT.pset%initialized )call mpp_error( FATAL, & + 'MPP_PSET_STACK_PUSH: called with uninitialized PSET.' ) + if( pset%lstack+len.GT.pset%maxstack )then + write( text, '(a,3i12)' ) & + 'MPP_PSET_STACK_PUSH: mpp_pset_stack overflow: '// & + 'len+lstack.GT.maxstack. len, lstack, maxstack=', & + len, pset%lstack, pset%maxstack + call mpp_error( FATAL, text ) + end if + p = pset%p_stack !point stack to shared stack pointer + ptr = LOC( stack(pset%lstack+1) ) + call mpp_pset_check_ptr(pset,ptr) !make sure ptr is the same across PSETs + pset%lstack = pset%lstack + len + pset%hiWM = max( pset%hiWM, pset%lstack ) +#else + integer(POINTER_KIND), intent(out) :: ptr + call mpp_error( FATAL, & + 'MPP_PSET_STACK_PUSH only works with Cray pointers.' ) +#endif + end subroutine mpp_pset_stack_push + + subroutine mpp_pset_stack_reset(pset) + type(mpp_pset_type), intent(inout) :: pset +!reset stack... will reuse any temporary arrays! USE WITH CARE +!next few lines are to zero stack contents... +!but it's better noone tries to use uninitialized stack variables! +! integer :: l1, l2 +! real :: mpp_pset_stack(maxstack) +! pointer( p_mpp_pset_stack, mpp_pset_stack ) +! p_mpp_pset_stack = ptr_mpp_pset_stack +! call mpp_pset_array_segment( 1, lstack, l1, l2 ) +! mpp_pset_stack(l1:l2) = 0. + if( .NOT.pset%initialized )call mpp_error( FATAL, & + 'MPP_PSET_STACK_RESET: called with uninitialized PSET.' ) + pset%lstack = 0 + end subroutine mpp_pset_stack_reset + + subroutine mpp_pset_print_chksum_1D(pset, caller, array) +!print a checksum of an array +!pass the whole domain seen by root PSET +!add lines to check on shared array? + type(mpp_pset_type), intent(in) :: pset + character(len=*), intent(in) :: caller + real, intent(in) :: array(:) + integer :: errunit + +#ifdef PSET_DEBUG + logical :: do_print + integer(LONG_KIND) :: chksum + + if( .NOT.pset%initialized )call mpp_error( FATAL, & + 'MPP_PSET_PRINT_CHKSUM: called with uninitialized PSET.' ) + errunit = stderr() + + if( pset%root )then + do_print = pe.EQ.mpp_root_pe() !set to T to print from all PEs + call mpp_set_current_pelist(pset%root_pelist) + chksum = mpp_chksum( array ) + if( do_print ) & + write( errunit, '(a,z18)' )trim(caller)//' chksum=', chksum + end if + call mpp_set_current_pelist(pset%pelist) +#endif + return + end subroutine mpp_pset_print_chksum_1D + + subroutine mpp_pset_print_chksum_2D(pset, caller, array) + type(mpp_pset_type), intent(in) :: pset + character(len=*), intent(in) :: caller + real, intent(in) :: array(:,:) + real :: array1D( size(array) ) +#ifdef use_CRI_pointers + pointer( p, array1D ) + p = LOC(array) +#else + array1D = TRANSFER( array, array1D ) +#endif + call mpp_pset_print_chksum(pset, caller, array1D) + end subroutine mpp_pset_print_chksum_2D + + subroutine mpp_pset_print_chksum_3D(pset, caller, array) + type(mpp_pset_type), intent(in) :: pset + character(len=*), intent(in) :: caller + real, intent(in) :: array(:,:,:) !overload for other ranks + real :: array1D( size(array) ) +#ifdef use_CRI_pointers + pointer( p, array1D ) + p = LOC(array) +#else + array1D = TRANSFER( array, array1D ) +#endif + call mpp_pset_print_chksum(pset, caller, array1D) + end subroutine mpp_pset_print_chksum_3D + + subroutine mpp_pset_print_chksum_4D(pset, caller, array) + type(mpp_pset_type), intent(in) :: pset + character(len=*), intent(in) :: caller + real, intent(in) :: array(:,:,:,:) + real :: array1D( size(array) ) +#ifdef use_CRI_pointers + pointer( p, array1D ) + p = LOC(array) +#else + array1D = TRANSFER( array, array1D ) +#endif + call mpp_pset_print_chksum(pset, caller, array1D) + end subroutine mpp_pset_print_chksum_4D + + subroutine mpp_pset_print_stack_chksum( pset, caller ) + type(mpp_pset_type), intent(in) :: pset + character(len=*), intent(in) :: caller + + if( .NOT.pset%initialized )call mpp_error( FATAL, & + 'MPP_PSET_PRINT_STACK_CHKSUM: called with uninitialized PSET.' ) + call mpp_pset_print_chksum( pset, trim(caller)//' stack', & + pset%stack(1:pset%lstack) ) + end subroutine mpp_pset_print_stack_chksum + +!accessor functions + function mpp_pset_root(pset) + logical :: mpp_pset_root + type(mpp_pset_type), intent(in) :: pset + + if( .NOT.pset%initialized )call mpp_error( FATAL, & + 'MPP_PSET_ROOT: called with uninitialized PSET.' ) + mpp_pset_root = pset%root + end function mpp_pset_root + + function mpp_pset_numroots(pset) +!necessary to export root_pelist: caller needs to pre-allocate + integer :: mpp_pset_numroots + type(mpp_pset_type), intent(in) :: pset + + if( .NOT.pset%initialized )call mpp_error( FATAL, & + 'MPP_PSET_NUMROOTS: called with uninitialized PSET.' ) + mpp_pset_numroots = size(pset%root_pelist) + end function mpp_pset_numroots + + subroutine mpp_pset_get_root_pelist(pset,pelist,commID) + type(mpp_pset_type), intent(in) :: pset + integer, intent(out) :: pelist(:) + integer, intent(out), optional :: commID + + if( .NOT.pset%initialized )call mpp_error( FATAL, & + 'MPP_PSET_GET_ROOT_PELIST: called with uninitialized PSET.' ) + if( size(pelist).NE.size(pset%root_pelist) )then + write( text,'(a,2i6)' ) & + 'pelist argument has wrong size: requested, actual=', & + size(pelist), size(pset%root_pelist) + call mpp_error( FATAL, 'MPP_PSET_GET_ROOT_PELIST: '//text ) + end if + pelist(:) = pset%root_pelist(:) + if( PRESENT(commID) )then +#ifdef use_libMPI + commID = pset%commID +#else + call mpp_error( WARNING, & + 'MPP_PSET_GET_ROOT_PELIST: commID is only defined under -Duse_libMPI.' ) +#endif + end if + end subroutine mpp_pset_get_root_pelist + +end module mpp_pset_mod