Module shr_mpi_mod 1,1
!-------------------------------------------------------------------------------
! PURPOSE: general layer on MPI functions
!-------------------------------------------------------------------------------
use shr_kind_mod
implicit none
private
! PUBLIC: Public interfaces
public :: shr_mpi_chkerr
public :: shr_mpi_send
public :: shr_mpi_recv
public :: shr_mpi_bcast
public :: shr_mpi_sum
public :: shr_mpi_min
public :: shr_mpi_max
public :: shr_mpi_commsize
public :: shr_mpi_commrank
public :: shr_mpi_initialized
public :: shr_mpi_abort
public :: shr_mpi_barrier
public :: shr_mpi_init
public :: shr_mpi_finalize
interface shr_mpi_send ; module procedure &
shr_mpi_sendi0, &
shr_mpi_sendi1, &
shr_mpi_sendr0, &
shr_mpi_sendr1
end interface
interface shr_mpi_recv ; module procedure &
shr_mpi_recvi0, &
shr_mpi_recvi1, &
shr_mpi_recvr0, &
shr_mpi_recvr1
end interface
interface shr_mpi_bcast ; module procedure &
shr_mpi_bcasti0, &
shr_mpi_bcasti1, &
shr_mpi_bcastr0, &
shr_mpi_bcastr1
end interface
interface shr_mpi_sum ; module procedure &
shr_mpi_sumi0, &
shr_mpi_sumi1, &
shr_mpi_sumr0, &
shr_mpi_sumr1, &
shr_mpi_sumr2, &
shr_mpi_sumr3
end interface
interface shr_mpi_min ; module procedure &
shr_mpi_mini0, &
shr_mpi_mini1, &
shr_mpi_minr0, &
shr_mpi_minr1
end interface
interface shr_mpi_max ; module procedure &
shr_mpi_maxi0, &
shr_mpi_maxi1, &
shr_mpi_maxr0, &
shr_mpi_maxr1
end interface
#if (! defined HIDE_MPI)
#include <mpif.h>
! mpi library include file
#endif
CONTAINS
!===============================================================================
SUBROUTINE shr_mpi_chkerr(rcode,string) 64,1
IMPLICIT none
!----- arguments ---
integer(SHR_KIND_IN), intent(in) :: rcode ! input MPI error code
character(*), intent(in) :: string ! message
!----- local ---
character(*),parameter :: subName = '(shr_mpi_chkerr) '
#if (! defined HIDE_MPI)
character(MPI_MAX_ERROR_STRING) :: lstring
#endif
integer(SHR_KIND_IN) :: len
!-------------------------------------------------------------------------------
! PURPOSE: layer on MPI error checking
!-------------------------------------------------------------------------------
#if (! defined HIDE_MPI)
if (rcode /= MPI_SUCCESS) then
! call MPI_ERROR_STRING(rcode,lstring,len)
! write(6,*) trim(subName),":",lstring(1:len)
write(6,*) trim(subName),":",trim(string)
call shr_mpi_abort
(subName)
endif
#endif
END SUBROUTINE shr_mpi_chkerr
!===============================================================================
SUBROUTINE shr_mpi_sendi0(lvec,pid,tag,comm,string),2
IMPLICIT none
!----- arguments ---
integer(SHR_KIND_IN), intent(in) :: lvec ! send value
integer(SHR_KIND_IN), intent(in) :: pid ! pid to send to
integer(SHR_KIND_IN), intent(in) :: tag ! tag
integer(SHR_KIND_IN), intent(in) :: comm ! mpi communicator
character(*),optional,intent(in) :: string ! message
!----- local ---
character(*),parameter :: subName = '(shr_mpi_sendi0) '
integer(SHR_KIND_IN) :: lsize
integer(SHR_KIND_IN) :: ierr
!-------------------------------------------------------------------------------
! PURPOSE: Send a vector of reals
!-------------------------------------------------------------------------------
#if (! defined HIDE_MPI)
lsize = 1
call MPI_SEND(lvec,lsize,MPI_INTEGER,pid,tag,comm,ierr)
if (present(string)) then
call shr_mpi_chkerr
(ierr,subName//trim(string))
else
call shr_mpi_chkerr
(ierr,subName)
endif
#endif
END SUBROUTINE shr_mpi_sendi0
!===============================================================================
!===============================================================================
SUBROUTINE shr_mpi_sendi1(lvec,pid,tag,comm,string),2
IMPLICIT none
!----- arguments ---
integer(SHR_KIND_IN), intent(in) :: lvec(:) ! in/out local values
integer(SHR_KIND_IN), intent(in) :: pid ! pid to send to
integer(SHR_KIND_IN), intent(in) :: tag ! tag
integer(SHR_KIND_IN), intent(in) :: comm ! mpi communicator
character(*),optional,intent(in) :: string ! message
!----- local ---
character(*),parameter :: subName = '(shr_mpi_sendi1) '
integer(SHR_KIND_IN) :: lsize
integer(SHR_KIND_IN) :: ierr
!-------------------------------------------------------------------------------
! PURPOSE: Send a vector of reals
!-------------------------------------------------------------------------------
#if (! defined HIDE_MPI)
lsize = size(lvec)
call MPI_SEND(lvec,lsize,MPI_INTEGER,pid,tag,comm,ierr)
if (present(string)) then
call shr_mpi_chkerr
(ierr,subName//trim(string))
else
call shr_mpi_chkerr
(ierr,subName)
endif
#endif
END SUBROUTINE shr_mpi_sendi1
!===============================================================================
!===============================================================================
SUBROUTINE shr_mpi_sendr0(lvec,pid,tag,comm,string),2
IMPLICIT none
!----- arguments ---
real(SHR_KIND_R8), intent(in) :: lvec ! in/out local values
integer(SHR_KIND_IN), intent(in) :: pid ! pid to send to
integer(SHR_KIND_IN), intent(in) :: tag ! tag
integer(SHR_KIND_IN), intent(in) :: comm ! mpi communicator
character(*),optional,intent(in) :: string ! message
!----- local ---
character(*),parameter :: subName = '(shr_mpi_sendr0) '
integer(SHR_KIND_IN) :: lsize
integer(SHR_KIND_IN) :: ierr
!-------------------------------------------------------------------------------
! PURPOSE: Send a vector of reals
!-------------------------------------------------------------------------------
#if (! defined HIDE_MPI)
lsize = 1
call MPI_SEND(lvec,lsize,MPI_REAL8,pid,tag,comm,ierr)
if (present(string)) then
call shr_mpi_chkerr
(ierr,subName//trim(string))
else
call shr_mpi_chkerr
(ierr,subName)
endif
#endif
END SUBROUTINE shr_mpi_sendr0
!===============================================================================
!===============================================================================
SUBROUTINE shr_mpi_sendr1(lvec,pid,tag,comm,string),2
IMPLICIT none
!----- arguments ---
real(SHR_KIND_R8), intent(in) :: lvec(:) ! in/out local values
integer(SHR_KIND_IN), intent(in) :: pid ! pid to send to
integer(SHR_KIND_IN), intent(in) :: tag ! tag
integer(SHR_KIND_IN), intent(in) :: comm ! mpi communicator
character(*),optional,intent(in) :: string ! message
!----- local ---
character(*),parameter :: subName = '(shr_mpi_sendr1) '
integer(SHR_KIND_IN) :: lsize
integer(SHR_KIND_IN) :: ierr
!-------------------------------------------------------------------------------
! PURPOSE: Send a vector of reals
!-------------------------------------------------------------------------------
#if (! defined HIDE_MPI)
lsize = size(lvec)
call MPI_SEND(lvec,lsize,MPI_REAL8,pid,tag,comm,ierr)
if (present(string)) then
call shr_mpi_chkerr
(ierr,subName//trim(string))
else
call shr_mpi_chkerr
(ierr,subName)
endif
#endif
END SUBROUTINE shr_mpi_sendr1
!===============================================================================
!===============================================================================
SUBROUTINE shr_mpi_recvi0(lvec,pid,tag,comm,string),2
IMPLICIT none
!----- arguments ---
integer(SHR_KIND_IN), intent(out):: lvec ! in/out local values
integer(SHR_KIND_IN), intent(in) :: pid ! pid to recv from
integer(SHR_KIND_IN), intent(in) :: tag ! tag
integer(SHR_KIND_IN), intent(in) :: comm ! mpi communicator
character(*),optional,intent(in) :: string ! message
!----- local ---
character(*),parameter :: subName = '(shr_mpi_recvi0) '
integer(SHR_KIND_IN) :: lsize
#if (! defined HIDE_MPI)
integer(SHR_KIND_IN) :: status(MPI_STATUS_SIZE) ! mpi status info
#endif
integer(SHR_KIND_IN) :: ierr
!-------------------------------------------------------------------------------
! PURPOSE: Recv a vector of reals
!-------------------------------------------------------------------------------
#if (! defined HIDE_MPI)
lsize = 1
call MPI_RECV(lvec,lsize,MPI_INTEGER,pid,tag,comm,status,ierr)
if (present(string)) then
call shr_mpi_chkerr
(ierr,subName//trim(string))
else
call shr_mpi_chkerr
(ierr,subName)
endif
#endif
END SUBROUTINE shr_mpi_recvi0
!===============================================================================
!===============================================================================
SUBROUTINE shr_mpi_recvi1(lvec,pid,tag,comm,string),2
IMPLICIT none
!----- arguments ---
integer(SHR_KIND_IN), intent(out):: lvec(:) ! in/out local values
integer(SHR_KIND_IN), intent(in) :: pid ! pid to recv from
integer(SHR_KIND_IN), intent(in) :: tag ! tag
integer(SHR_KIND_IN), intent(in) :: comm ! mpi communicator
character(*),optional,intent(in) :: string ! message
!----- local ---
character(*),parameter :: subName = '(shr_mpi_recvi1) '
integer(SHR_KIND_IN) :: lsize
#if (! defined HIDE_MPI)
integer(SHR_KIND_IN) :: status(MPI_STATUS_SIZE) ! mpi status info
#endif
integer(SHR_KIND_IN) :: ierr
!-------------------------------------------------------------------------------
! PURPOSE: Recv a vector of reals
!-------------------------------------------------------------------------------
#if (! defined HIDE_MPI)
lsize = size(lvec)
call MPI_RECV(lvec,lsize,MPI_INTEGER,pid,tag,comm,status,ierr)
if (present(string)) then
call shr_mpi_chkerr
(ierr,subName//trim(string))
else
call shr_mpi_chkerr
(ierr,subName)
endif
#endif
END SUBROUTINE shr_mpi_recvi1
!===============================================================================
!===============================================================================
SUBROUTINE shr_mpi_recvr0(lvec,pid,tag,comm,string),2
IMPLICIT none
!----- arguments ---
real(SHR_KIND_R8), intent(out):: lvec ! in/out local values
integer(SHR_KIND_IN), intent(in) :: pid ! pid to recv from
integer(SHR_KIND_IN), intent(in) :: tag ! tag
integer(SHR_KIND_IN), intent(in) :: comm ! mpi communicator
character(*),optional,intent(in) :: string ! message
!----- local ---
character(*),parameter :: subName = '(shr_mpi_recvr0) '
integer(SHR_KIND_IN) :: lsize
#if (! defined HIDE_MPI)
integer(SHR_KIND_IN) :: status(MPI_STATUS_SIZE) ! mpi status info
#endif
integer(SHR_KIND_IN) :: ierr
!-------------------------------------------------------------------------------
! PURPOSE: Recv a vector of reals
!-------------------------------------------------------------------------------
#if (! defined HIDE_MPI)
lsize = 1
call MPI_RECV(lvec,lsize,MPI_REAL8,pid,tag,comm,status,ierr)
if (present(string)) then
call shr_mpi_chkerr
(ierr,subName//trim(string))
else
call shr_mpi_chkerr
(ierr,subName)
endif
#endif
END SUBROUTINE shr_mpi_recvr0
!===============================================================================
!===============================================================================
SUBROUTINE shr_mpi_recvr1(lvec,pid,tag,comm,string),2
IMPLICIT none
!----- arguments ---
real(SHR_KIND_R8), intent(out):: lvec(:) ! in/out local values
integer(SHR_KIND_IN), intent(in) :: pid ! pid to recv from
integer(SHR_KIND_IN), intent(in) :: tag ! tag
integer(SHR_KIND_IN), intent(in) :: comm ! mpi communicator
character(*),optional,intent(in) :: string ! message
!----- local ---
character(*),parameter :: subName = '(shr_mpi_recvr1) '
integer(SHR_KIND_IN) :: lsize
#if (! defined HIDE_MPI)
integer(SHR_KIND_IN) :: status(MPI_STATUS_SIZE) ! mpi status info
#endif
integer(SHR_KIND_IN) :: ierr
!-------------------------------------------------------------------------------
! PURPOSE: Recv a vector of reals
!-------------------------------------------------------------------------------
#if (! defined HIDE_MPI)
lsize = size(lvec)
call MPI_RECV(lvec,lsize,MPI_REAL8,pid,tag,comm,status,ierr)
if (present(string)) then
call shr_mpi_chkerr
(ierr,subName//trim(string))
else
call shr_mpi_chkerr
(ierr,subName)
endif
#endif
END SUBROUTINE shr_mpi_recvr1
!===============================================================================
!===============================================================================
SUBROUTINE shr_mpi_bcasti0(vec,comm,string),2
IMPLICIT none
!----- arguments ---
integer(SHR_KIND_IN), intent(inout):: vec ! vector of 1
integer(SHR_KIND_IN), intent(in) :: comm ! mpi communicator
character(*),optional,intent(in) :: string ! message
!----- local ---
character(*),parameter :: subName = '(shr_mpi_bcasti0) '
integer(SHR_KIND_IN) :: ierr
integer(SHR_KIND_IN) :: lsize
!-------------------------------------------------------------------------------
! PURPOSE: Broadcast an integer
!-------------------------------------------------------------------------------
#if (! defined HIDE_MPI)
lsize = 1
call MPI_BCAST(vec,lsize,MPI_INTEGER,0,comm,ierr)
if (present(string)) then
call shr_mpi_chkerr
(ierr,subName//trim(string))
else
call shr_mpi_chkerr
(ierr,subName)
endif
#endif
END SUBROUTINE shr_mpi_bcasti0
!===============================================================================
!===============================================================================
SUBROUTINE shr_mpi_bcastr0(vec,comm,string),2
IMPLICIT none
!----- arguments ---
real(SHR_KIND_R8), intent(inout):: vec ! vector of 1
integer(SHR_KIND_IN), intent(in) :: comm ! mpi communicator
character(*),optional,intent(in) :: string ! message
!----- local ---
character(*),parameter :: subName = '(shr_mpi_bcastr0) '
integer(SHR_KIND_IN) :: ierr
integer(SHR_KIND_IN) :: lsize
!-------------------------------------------------------------------------------
! PURPOSE: Broadcast a real
!-------------------------------------------------------------------------------
#if (! defined HIDE_MPI)
lsize = 1
call MPI_BCAST(vec,lsize,MPI_REAL8,0,comm,ierr)
if (present(string)) then
call shr_mpi_chkerr
(ierr,subName//trim(string))
else
call shr_mpi_chkerr
(ierr,subName)
endif
#endif
END SUBROUTINE shr_mpi_bcastr0
!===============================================================================
!===============================================================================
SUBROUTINE shr_mpi_bcasti1(vec,comm,string),2
IMPLICIT none
!----- arguments ---
integer(SHR_KIND_IN), intent(inout):: vec(:) ! vector
integer(SHR_KIND_IN), intent(in) :: comm ! mpi communicator
character(*),optional,intent(in) :: string ! message
!----- local ---
character(*),parameter :: subName = '(shr_mpi_bcasti1) '
integer(SHR_KIND_IN) :: ierr
integer(SHR_KIND_IN) :: lsize
!-------------------------------------------------------------------------------
! PURPOSE: Broadcast a vector of integers
!-------------------------------------------------------------------------------
#if (! defined HIDE_MPI)
lsize = size(vec)
call MPI_BCAST(vec,lsize,MPI_INTEGER,0,comm,ierr)
if (present(string)) then
call shr_mpi_chkerr
(ierr,subName//trim(string))
else
call shr_mpi_chkerr
(ierr,subName)
endif
#endif
END SUBROUTINE shr_mpi_bcasti1
!===============================================================================
!===============================================================================
SUBROUTINE shr_mpi_bcastr1(vec,comm,string),2
IMPLICIT none
!----- arguments ---
real(SHR_KIND_R8), intent(inout):: vec(:) ! vector
integer(SHR_KIND_IN), intent(in) :: comm ! mpi communicator
character(*),optional,intent(in) :: string ! message
!----- local ---
character(*),parameter :: subName = '(shr_mpi_bcastr1) '
integer(SHR_KIND_IN) :: ierr
integer(SHR_KIND_IN) :: lsize
!-------------------------------------------------------------------------------
! PURPOSE: Broadcast a vector of reals
!-------------------------------------------------------------------------------
#if (! defined HIDE_MPI)
lsize = size(vec)
call MPI_BCAST(vec,lsize,MPI_REAL8,0,comm,ierr)
if (present(string)) then
call shr_mpi_chkerr
(ierr,subName//trim(string))
else
call shr_mpi_chkerr
(ierr,subName)
endif
#endif
END SUBROUTINE shr_mpi_bcastr1
!===============================================================================
!===============================================================================
SUBROUTINE shr_mpi_sumi0(lvec,gvec,comm,string,all),3
IMPLICIT none
!----- arguments ---
integer(SHR_KIND_IN), intent(in) :: lvec ! in/out local values
integer(SHR_KIND_IN), intent(out):: gvec ! in/out global values
integer(SHR_KIND_IN), intent(in) :: comm ! mpi communicator
character(*),optional,intent(in) :: string ! message
logical, optional,intent(in) :: all ! allreduce if true
!----- local ---
character(*),parameter :: subName = '(shr_mpi_sumi0) '
logical :: lall
character(SHR_KIND_CL) :: lstring
integer(SHR_KIND_IN) :: reduce_type ! mpi reduction type
integer(SHR_KIND_IN) :: lsize
integer(SHR_KIND_IN) :: gsize
integer(SHR_KIND_IN) :: ierr
!-------------------------------------------------------------------------------
! PURPOSE: Finds sum of a distributed vector of values, assume local sum
! already computed
!-------------------------------------------------------------------------------
#if (! defined HIDE_MPI)
reduce_type = MPI_SUM
if (present(all)) then
lall = all
else
lall = .false.
endif
if (present(string)) then
lstring = trim(subName)//":"//trim(string)
else
lstring = trim(subName)
endif
lsize = 1
gsize = 1
if (lsize /= gsize) then
call shr_mpi_abort
(subName//" lsize,gsize incompatable "//trim(string))
endif
if (lall) then
call MPI_ALLREDUCE(lvec,gvec,gsize,MPI_INTEGER,reduce_type,comm,ierr)
call shr_mpi_chkerr
(ierr,trim(lstring)//" MPI_ALLREDUCE")
else
call MPI_REDUCE(lvec,gvec,gsize,MPI_INTEGER,reduce_type,0,comm,ierr)
call shr_mpi_chkerr
(ierr,trim(lstring)//" MPI_REDUCE")
endif
#endif
END SUBROUTINE shr_mpi_sumi0
!===============================================================================
!===============================================================================
SUBROUTINE shr_mpi_sumi1(lvec,gvec,comm,string,all),3
IMPLICIT none
!----- arguments ---
integer(SHR_KIND_IN), intent(in) :: lvec(:) ! in/out local values
integer(SHR_KIND_IN), intent(out):: gvec(:) ! in/out global values
integer(SHR_KIND_IN), intent(in) :: comm ! mpi communicator
character(*),optional,intent(in) :: string ! message
logical, optional,intent(in) :: all ! allreduce if true
!----- local ---
character(*),parameter :: subName = '(shr_mpi_sumi1) '
logical :: lall
character(SHR_KIND_CL) :: lstring
integer(SHR_KIND_IN) :: reduce_type ! mpi reduction type
integer(SHR_KIND_IN) :: lsize
integer(SHR_KIND_IN) :: gsize
integer(SHR_KIND_IN) :: ierr
!-------------------------------------------------------------------------------
! PURPOSE: Finds sum of a distributed vector of values, assume local sum
! already computed
!-------------------------------------------------------------------------------
#if (! defined HIDE_MPI)
reduce_type = MPI_SUM
if (present(all)) then
lall = all
else
lall = .false.
endif
if (present(string)) then
lstring = trim(subName)//":"//trim(string)
else
lstring = trim(subName)
endif
lsize = size(lvec)
gsize = size(gvec)
if (lsize /= gsize) then
call shr_mpi_abort
(subName//" lsize,gsize incompatable "//trim(string))
endif
if (lall) then
call MPI_ALLREDUCE(lvec,gvec,gsize,MPI_INTEGER,reduce_type,comm,ierr)
call shr_mpi_chkerr
(ierr,trim(lstring)//" MPI_ALLREDUCE")
else
call MPI_REDUCE(lvec,gvec,gsize,MPI_INTEGER,reduce_type,0,comm,ierr)
call shr_mpi_chkerr
(ierr,trim(lstring)//" MPI_REDUCE")
endif
#endif
END SUBROUTINE shr_mpi_sumi1
!===============================================================================
!===============================================================================
SUBROUTINE shr_mpi_sumr0(lvec,gvec,comm,string,all),3
IMPLICIT none
!----- arguments ---
real(SHR_KIND_R8), intent(in) :: lvec ! in/out local values
real(SHR_KIND_R8), intent(out):: gvec ! in/out global values
integer(SHR_KIND_IN), intent(in) :: comm ! mpi communicator
character(*),optional,intent(in) :: string ! message
logical, optional,intent(in) :: all ! allreduce if true
!----- local ---
character(*),parameter :: subName = '(shr_mpi_sumr0) '
logical :: lall
character(SHR_KIND_CL) :: lstring
integer(SHR_KIND_IN) :: reduce_type ! mpi reduction type
integer(SHR_KIND_IN) :: lsize
integer(SHR_KIND_IN) :: gsize
integer(SHR_KIND_IN) :: ierr
!-------------------------------------------------------------------------------
! PURPOSE: Finds sum of a distributed vector of values, assume local sum
! already computed
!-------------------------------------------------------------------------------
#if (! defined HIDE_MPI)
reduce_type = MPI_SUM
if (present(all)) then
lall = all
else
lall = .false.
endif
if (present(string)) then
lstring = trim(subName)//":"//trim(string)
else
lstring = trim(subName)
endif
lsize = 1
gsize = 1
if (lsize /= gsize) then
call shr_mpi_abort
(subName//" lsize,gsize incompatable "//trim(string))
endif
if (lall) then
call MPI_ALLREDUCE(lvec,gvec,gsize,MPI_REAL8,reduce_type,comm,ierr)
call shr_mpi_chkerr
(ierr,trim(lstring)//" MPI_ALLREDUCE")
else
call MPI_REDUCE(lvec,gvec,gsize,MPI_REAL8,reduce_type,0,comm,ierr)
call shr_mpi_chkerr
(ierr,trim(lstring)//" MPI_REDUCE")
endif
#endif
END SUBROUTINE shr_mpi_sumr0
!===============================================================================
!===============================================================================
SUBROUTINE shr_mpi_sumr1(lvec,gvec,comm,string,all),3
IMPLICIT none
!----- arguments ---
real(SHR_KIND_R8), intent(in) :: lvec(:) ! in/out local values
real(SHR_KIND_R8), intent(out):: gvec(:) ! in/out global values
integer(SHR_KIND_IN), intent(in) :: comm ! mpi communicator
character(*),optional,intent(in) :: string ! message
logical, optional,intent(in) :: all ! allreduce if true
!----- local ---
character(*),parameter :: subName = '(shr_mpi_sumr1) '
logical :: lall
character(SHR_KIND_CL) :: lstring
integer(SHR_KIND_IN) :: reduce_type ! mpi reduction type
integer(SHR_KIND_IN) :: lsize
integer(SHR_KIND_IN) :: gsize
integer(SHR_KIND_IN) :: ierr
!-------------------------------------------------------------------------------
! PURPOSE: Finds sum of a distributed vector of values, assume local sum
! already computed
!-------------------------------------------------------------------------------
#if (! defined HIDE_MPI)
reduce_type = MPI_SUM
if (present(all)) then
lall = all
else
lall = .false.
endif
if (present(string)) then
lstring = trim(subName)//":"//trim(string)
else
lstring = trim(subName)
endif
lsize = size(lvec)
gsize = size(gvec)
if (lsize /= gsize) then
call shr_mpi_abort
(subName//" lsize,gsize incompatable "//trim(string))
endif
if (lall) then
call MPI_ALLREDUCE(lvec,gvec,gsize,MPI_REAL8,reduce_type,comm,ierr)
call shr_mpi_chkerr
(ierr,trim(lstring)//" MPI_ALLREDUCE")
else
call MPI_REDUCE(lvec,gvec,gsize,MPI_REAL8,reduce_type,0,comm,ierr)
call shr_mpi_chkerr
(ierr,trim(lstring)//" MPI_REDUCE")
endif
#endif
END SUBROUTINE shr_mpi_sumr1
!===============================================================================
!===============================================================================
SUBROUTINE shr_mpi_sumr2(lvec,gvec,comm,string,all),3
IMPLICIT none
!----- arguments ---
real(SHR_KIND_R8), intent(in) :: lvec(:,:)! in/out local values
real(SHR_KIND_R8), intent(out):: gvec(:,:)! in/out global values
integer(SHR_KIND_IN), intent(in) :: comm ! mpi communicator
character(*),optional,intent(in) :: string ! message
logical, optional,intent(in) :: all ! allreduce if true
!----- local ---
character(*),parameter :: subName = '(shr_mpi_sumr2) '
logical :: lall
character(SHR_KIND_CL) :: lstring
integer(SHR_KIND_IN) :: reduce_type ! mpi reduction type
integer(SHR_KIND_IN) :: lsize
integer(SHR_KIND_IN) :: gsize
integer(SHR_KIND_IN) :: ierr
!-------------------------------------------------------------------------------
! PURPOSE: Finds sum of a distributed vector of values, assume local sum
! already computed
!-------------------------------------------------------------------------------
#if (! defined HIDE_MPI)
reduce_type = MPI_SUM
if (present(all)) then
lall = all
else
lall = .false.
endif
if (present(string)) then
lstring = trim(subName)//":"//trim(string)
else
lstring = trim(subName)
endif
lsize = size(lvec)
gsize = size(gvec)
if (lsize /= gsize) then
call shr_mpi_abort
(subName//" lsize,gsize incompatable "//trim(string))
endif
if (lall) then
call MPI_ALLREDUCE(lvec,gvec,gsize,MPI_REAL8,reduce_type,comm,ierr)
call shr_mpi_chkerr
(ierr,trim(lstring)//" MPI_ALLREDUCE")
else
call MPI_REDUCE(lvec,gvec,gsize,MPI_REAL8,reduce_type,0,comm,ierr)
call shr_mpi_chkerr
(ierr,trim(lstring)//" MPI_REDUCE")
endif
#endif
END SUBROUTINE shr_mpi_sumr2
!===============================================================================
!===============================================================================
SUBROUTINE shr_mpi_sumr3(lvec,gvec,comm,string,all),3
IMPLICIT none
!----- arguments ---
real(SHR_KIND_R8), intent(in) :: lvec(:,:,:) ! in/out local values
real(SHR_KIND_R8), intent(out):: gvec(:,:,:) ! in/out global values
integer(SHR_KIND_IN), intent(in) :: comm ! mpi communicator
character(*),optional,intent(in) :: string ! message
logical, optional,intent(in) :: all ! allreduce if true
!----- local ---
character(*),parameter :: subName = '(shr_mpi_sumr3) '
logical :: lall
character(SHR_KIND_CL) :: lstring
integer(SHR_KIND_IN) :: reduce_type ! mpi reduction type
integer(SHR_KIND_IN) :: lsize
integer(SHR_KIND_IN) :: gsize
integer(SHR_KIND_IN) :: ierr
!-------------------------------------------------------------------------------
! PURPOSE: Finds sum of a distributed vector of values, assume local sum
! already computed
!-------------------------------------------------------------------------------
#if (! defined HIDE_MPI)
reduce_type = MPI_SUM
if (present(all)) then
lall = all
else
lall = .false.
endif
if (present(string)) then
lstring = trim(subName)//":"//trim(string)
else
lstring = trim(subName)
endif
lsize = size(lvec)
gsize = size(gvec)
if (lsize /= gsize) then
call shr_mpi_abort
(subName//" lsize,gsize incompatable "//trim(string))
endif
if (lall) then
call MPI_ALLREDUCE(lvec,gvec,gsize,MPI_REAL8,reduce_type,comm,ierr)
call shr_mpi_chkerr
(ierr,trim(lstring)//" MPI_ALLREDUCE")
else
call MPI_REDUCE(lvec,gvec,gsize,MPI_REAL8,reduce_type,0,comm,ierr)
call shr_mpi_chkerr
(ierr,trim(lstring)//" MPI_REDUCE")
endif
#endif
END SUBROUTINE shr_mpi_sumr3
!===============================================================================
!===============================================================================
SUBROUTINE shr_mpi_mini0(lvec,gvec,comm,string,all),3
IMPLICIT none
!----- arguments ---
integer(SHR_KIND_IN), intent(in) :: lvec ! in/out local values
integer(SHR_KIND_IN), intent(out):: gvec ! in/out global values
integer(SHR_KIND_IN), intent(in) :: comm ! mpi communicator
character(*),optional,intent(in) :: string ! message
logical, optional,intent(in) :: all ! allreduce if true
!----- local ---
character(*),parameter :: subName = '(shr_mpi_mini0) '
logical :: lall
character(SHR_KIND_CL) :: lstring
integer(SHR_KIND_IN) :: reduce_type ! mpi reduction type
integer(SHR_KIND_IN) :: lsize
integer(SHR_KIND_IN) :: gsize
integer(SHR_KIND_IN) :: ierr
!-------------------------------------------------------------------------------
! PURPOSE: Finds min of a distributed vector of values, assume local min
! already computed
!-------------------------------------------------------------------------------
#if (! defined HIDE_MPI)
reduce_type = MPI_MIN
if (present(all)) then
lall = all
else
lall = .false.
endif
if (present(string)) then
lstring = trim(subName)//":"//trim(string)
else
lstring = trim(subName)
endif
lsize = 1
gsize = 1
if (lsize /= gsize) then
call shr_mpi_abort
(subName//" lsize,gsize incompatable "//trim(string))
endif
if (lall) then
call MPI_ALLREDUCE(lvec,gvec,gsize,MPI_INTEGER,reduce_type,comm,ierr)
call shr_mpi_chkerr
(ierr,trim(lstring)//" MPI_ALLREDUCE")
else
call MPI_REDUCE(lvec,gvec,gsize,MPI_INTEGER,reduce_type,0,comm,ierr)
call shr_mpi_chkerr
(ierr,trim(lstring)//" MPI_REDUCE")
endif
#endif
END SUBROUTINE shr_mpi_mini0
!===============================================================================
!===============================================================================
SUBROUTINE shr_mpi_mini1(lvec,gvec,comm,string,all),3
IMPLICIT none
!----- arguments ---
integer(SHR_KIND_IN), intent(in) :: lvec(:) ! in/out local values
integer(SHR_KIND_IN), intent(out):: gvec(:) ! in/out global values
integer(SHR_KIND_IN), intent(in) :: comm ! mpi communicator
character(*),optional,intent(in) :: string ! message
logical, optional,intent(in) :: all ! allreduce if true
!----- local ---
character(*),parameter :: subName = '(shr_mpi_mini1) '
logical :: lall
character(SHR_KIND_CL) :: lstring
integer(SHR_KIND_IN) :: reduce_type ! mpi reduction type
integer(SHR_KIND_IN) :: lsize
integer(SHR_KIND_IN) :: gsize
integer(SHR_KIND_IN) :: ierr
!-------------------------------------------------------------------------------
! PURPOSE: Finds min of a distributed vector of values, assume local min
! already computed
!-------------------------------------------------------------------------------
#if (! defined HIDE_MPI)
reduce_type = MPI_MIN
if (present(all)) then
lall = all
else
lall = .false.
endif
if (present(string)) then
lstring = trim(subName)//":"//trim(string)
else
lstring = trim(subName)
endif
lsize = size(lvec)
gsize = size(gvec)
if (lsize /= gsize) then
call shr_mpi_abort
(subName//" lsize,gsize incompatable "//trim(string))
endif
if (lall) then
call MPI_ALLREDUCE(lvec,gvec,gsize,MPI_INTEGER,reduce_type,comm,ierr)
call shr_mpi_chkerr
(ierr,trim(lstring)//" MPI_ALLREDUCE")
else
call MPI_REDUCE(lvec,gvec,gsize,MPI_INTEGER,reduce_type,0,comm,ierr)
call shr_mpi_chkerr
(ierr,trim(lstring)//" MPI_REDUCE")
endif
#endif
END SUBROUTINE shr_mpi_mini1
!===============================================================================
!===============================================================================
SUBROUTINE shr_mpi_minr0(lvec,gvec,comm,string,all),3
IMPLICIT none
!----- arguments ---
real(SHR_KIND_R8), intent(in) :: lvec ! in/out local values
real(SHR_KIND_R8), intent(out):: gvec ! in/out global values
integer(SHR_KIND_IN), intent(in) :: comm ! mpi communicator
character(*),optional,intent(in) :: string ! message
logical, optional,intent(in) :: all ! allreduce if true
!----- local ---
character(*),parameter :: subName = '(shr_mpi_minr0) '
logical :: lall
character(SHR_KIND_CL) :: lstring
integer(SHR_KIND_IN) :: reduce_type ! mpi reduction type
integer(SHR_KIND_IN) :: lsize
integer(SHR_KIND_IN) :: gsize
integer(SHR_KIND_IN) :: ierr
!-------------------------------------------------------------------------------
! PURPOSE: Finds min of a distributed vector of values, assume local min
! already computed
!-------------------------------------------------------------------------------
#if (! defined HIDE_MPI)
reduce_type = MPI_MIN
if (present(all)) then
lall = all
else
lall = .false.
endif
if (present(string)) then
lstring = trim(subName)//":"//trim(string)
else
lstring = trim(subName)
endif
lsize = 1
gsize = 1
if (lsize /= gsize) then
call shr_mpi_abort
(subName//" lsize,gsize incompatable "//trim(string))
endif
if (lall) then
call MPI_ALLREDUCE(lvec,gvec,gsize,MPI_REAL8,reduce_type,comm,ierr)
call shr_mpi_chkerr
(ierr,trim(lstring)//" MPI_ALLREDUCE")
else
call MPI_REDUCE(lvec,gvec,gsize,MPI_REAL8,reduce_type,0,comm,ierr)
call shr_mpi_chkerr
(ierr,trim(lstring)//" MPI_REDUCE")
endif
#endif
END SUBROUTINE shr_mpi_minr0
!===============================================================================
!===============================================================================
SUBROUTINE shr_mpi_minr1(lvec,gvec,comm,string,all),3
IMPLICIT none
!----- arguments ---
real(SHR_KIND_R8), intent(in) :: lvec(:) ! in/out local values
real(SHR_KIND_R8), intent(out):: gvec(:) ! in/out global values
integer(SHR_KIND_IN), intent(in) :: comm ! mpi communicator
character(*),optional,intent(in) :: string ! message
logical, optional,intent(in) :: all ! allreduce if true
!----- local ---
character(*),parameter :: subName = '(shr_mpi_minr1) '
logical :: lall
character(SHR_KIND_CL) :: lstring
integer(SHR_KIND_IN) :: reduce_type ! mpi reduction type
integer(SHR_KIND_IN) :: lsize
integer(SHR_KIND_IN) :: gsize
integer(SHR_KIND_IN) :: ierr
!-------------------------------------------------------------------------------
! PURPOSE: Finds min of a distributed vector of values, assume local min
! already computed
!-------------------------------------------------------------------------------
#if (! defined HIDE_MPI)
reduce_type = MPI_MIN
if (present(all)) then
lall = all
else
lall = .false.
endif
if (present(string)) then
lstring = trim(subName)//":"//trim(string)
else
lstring = trim(subName)
endif
lsize = size(lvec)
gsize = size(gvec)
if (lsize /= gsize) then
call shr_mpi_abort
(subName//" lsize,gsize incompatable "//trim(string))
endif
if (lall) then
call MPI_ALLREDUCE(lvec,gvec,gsize,MPI_REAL8,reduce_type,comm,ierr)
call shr_mpi_chkerr
(ierr,trim(lstring)//" MPI_ALLREDUCE")
else
call MPI_REDUCE(lvec,gvec,gsize,MPI_REAL8,reduce_type,0,comm,ierr)
call shr_mpi_chkerr
(ierr,trim(lstring)//" MPI_REDUCE")
endif
#endif
END SUBROUTINE shr_mpi_minr1
!===============================================================================
!===============================================================================
SUBROUTINE shr_mpi_maxi0(lvec,gvec,comm,string,all),3
IMPLICIT none
!----- arguments ---
integer(SHR_KIND_IN), intent(in) :: lvec ! in/out local values
integer(SHR_KIND_IN), intent(out):: gvec ! in/out global values
integer(SHR_KIND_IN), intent(in) :: comm ! mpi communicator
character(*),optional,intent(in) :: string ! message
logical, optional,intent(in) :: all ! allreduce if true
!----- local ---
character(*),parameter :: subName = '(shr_mpi_maxi0) '
logical :: lall
character(SHR_KIND_CL) :: lstring
integer(SHR_KIND_IN) :: reduce_type ! mpi reduction type
integer(SHR_KIND_IN) :: lsize
integer(SHR_KIND_IN) :: gsize
integer(SHR_KIND_IN) :: ierr
!-------------------------------------------------------------------------------
! PURPOSE: Finds max of a distributed vector of values, assume local max
! already computed
!-------------------------------------------------------------------------------
#if (! defined HIDE_MPI)
reduce_type = MPI_MAX
if (present(all)) then
lall = all
else
lall = .false.
endif
if (present(string)) then
lstring = trim(subName)//":"//trim(string)
else
lstring = trim(subName)
endif
lsize = 1
gsize = 1
if (lsize /= gsize) then
call shr_mpi_abort
(subName//" lsize,gsize incompatable "//trim(string))
endif
if (lall) then
call MPI_ALLREDUCE(lvec,gvec,gsize,MPI_INTEGER,reduce_type,comm,ierr)
call shr_mpi_chkerr
(ierr,trim(lstring)//" MPI_ALLREDUCE")
else
call MPI_REDUCE(lvec,gvec,gsize,MPI_INTEGER,reduce_type,0,comm,ierr)
call shr_mpi_chkerr
(ierr,trim(lstring)//" MPI_REDUCE")
endif
#endif
END SUBROUTINE shr_mpi_maxi0
!===============================================================================
!===============================================================================
SUBROUTINE shr_mpi_maxi1(lvec,gvec,comm,string,all),3
IMPLICIT none
!----- arguments ---
integer(SHR_KIND_IN), intent(in) :: lvec(:) ! in/out local values
integer(SHR_KIND_IN), intent(out):: gvec(:) ! in/out global values
integer(SHR_KIND_IN), intent(in) :: comm ! mpi communicator
character(*),optional,intent(in) :: string ! message
logical, optional,intent(in) :: all ! allreduce if true
!----- local ---
character(*),parameter :: subName = '(shr_mpi_maxi1) '
logical :: lall
character(SHR_KIND_CL) :: lstring
integer(SHR_KIND_IN) :: reduce_type ! mpi reduction type
integer(SHR_KIND_IN) :: lsize
integer(SHR_KIND_IN) :: gsize
integer(SHR_KIND_IN) :: ierr
!-------------------------------------------------------------------------------
! PURPOSE: Finds max of a distributed vector of values, assume local max
! already computed
!-------------------------------------------------------------------------------
#if (! defined HIDE_MPI)
reduce_type = MPI_MAX
if (present(all)) then
lall = all
else
lall = .false.
endif
if (present(string)) then
lstring = trim(subName)//":"//trim(string)
else
lstring = trim(subName)
endif
lsize = size(lvec)
gsize = size(gvec)
if (lsize /= gsize) then
call shr_mpi_abort
(subName//" lsize,gsize incompatable "//trim(string))
endif
if (lall) then
call MPI_ALLREDUCE(lvec,gvec,gsize,MPI_INTEGER,reduce_type,comm,ierr)
call shr_mpi_chkerr
(ierr,trim(lstring)//" MPI_ALLREDUCE")
else
call MPI_REDUCE(lvec,gvec,gsize,MPI_INTEGER,reduce_type,0,comm,ierr)
call shr_mpi_chkerr
(ierr,trim(lstring)//" MPI_REDUCE")
endif
#endif
END SUBROUTINE shr_mpi_maxi1
!===============================================================================
!===============================================================================
SUBROUTINE shr_mpi_maxr0(lvec,gvec,comm,string,all),3
IMPLICIT none
!----- arguments ---
real(SHR_KIND_R8), intent(in) :: lvec ! in/out local values
real(SHR_KIND_R8), intent(out):: gvec ! in/out global values
integer(SHR_KIND_IN), intent(in) :: comm ! mpi communicator
character(*),optional,intent(in) :: string ! message
logical, optional,intent(in) :: all ! allreduce if true
!----- local ---
character(*),parameter :: subName = '(shr_mpi_maxr0) '
logical :: lall
character(SHR_KIND_CL) :: lstring
integer(SHR_KIND_IN) :: reduce_type ! mpi reduction type
integer(SHR_KIND_IN) :: lsize
integer(SHR_KIND_IN) :: gsize
integer(SHR_KIND_IN) :: ierr
!-------------------------------------------------------------------------------
! PURPOSE: Finds max of a distributed vector of values, assume local max
! already computed
!-------------------------------------------------------------------------------
#if (! defined HIDE_MPI)
reduce_type = MPI_MAX
if (present(all)) then
lall = all
else
lall = .false.
endif
if (present(string)) then
lstring = trim(subName)//":"//trim(string)
else
lstring = trim(subName)
endif
lsize = 1
gsize = 1
if (lsize /= gsize) then
call shr_mpi_abort
(subName//" lsize,gsize incompatable "//trim(string))
endif
if (lall) then
call MPI_ALLREDUCE(lvec,gvec,gsize,MPI_REAL8,reduce_type,comm,ierr)
call shr_mpi_chkerr
(ierr,trim(lstring)//" MPI_ALLREDUCE")
else
call MPI_REDUCE(lvec,gvec,gsize,MPI_REAL8,reduce_type,0,comm,ierr)
call shr_mpi_chkerr
(ierr,trim(lstring)//" MPI_REDUCE")
endif
#endif
END SUBROUTINE shr_mpi_maxr0
!===============================================================================
!===============================================================================
SUBROUTINE shr_mpi_maxr1(lvec,gvec,comm,string,all),3
IMPLICIT none
!----- arguments ---
real(SHR_KIND_R8), intent(in) :: lvec(:) ! in/out local values
real(SHR_KIND_R8), intent(out):: gvec(:) ! in/out global values
integer(SHR_KIND_IN), intent(in) :: comm ! mpi communicator
character(*),optional,intent(in) :: string ! message
logical, optional,intent(in) :: all ! allreduce if true
!----- local ---
character(*),parameter :: subName = '(shr_mpi_maxr1) '
logical :: lall
character(SHR_KIND_CL) :: lstring
integer(SHR_KIND_IN) :: reduce_type ! mpi reduction type
integer(SHR_KIND_IN) :: lsize
integer(SHR_KIND_IN) :: gsize
integer(SHR_KIND_IN) :: ierr
!-------------------------------------------------------------------------------
! PURPOSE: Finds max of a distributed vector of values, assume local max
! already computed
!-------------------------------------------------------------------------------
#if (! defined HIDE_MPI)
reduce_type = MPI_MAX
if (present(all)) then
lall = all
else
lall = .false.
endif
if (present(string)) then
lstring = trim(subName)//":"//trim(string)
else
lstring = trim(subName)
endif
lsize = size(lvec)
gsize = size(gvec)
if (lsize /= gsize) then
call shr_mpi_abort
(subName//" lsize,gsize incompatable "//trim(string))
endif
if (lall) then
call MPI_ALLREDUCE(lvec,gvec,gsize,MPI_REAL8,reduce_type,comm,ierr)
call shr_mpi_chkerr
(ierr,trim(lstring)//" MPI_ALLREDUCE")
else
call MPI_REDUCE(lvec,gvec,gsize,MPI_REAL8,reduce_type,0,comm,ierr)
call shr_mpi_chkerr
(ierr,trim(lstring)//" MPI_REDUCE")
endif
#endif
END SUBROUTINE shr_mpi_maxr1
!===============================================================================
!===============================================================================
SUBROUTINE shr_mpi_commsize(comm,size,string),2
IMPLICIT none
!----- arguments ---
integer,intent(in) :: comm
integer,intent(out) :: size
character(*),optional,intent(in) :: string ! message
!----- local ---
character(*),parameter :: subName = '(shr_mpi_commsize) '
integer(SHR_KIND_IN) :: ierr
!-------------------------------------------------------------------------------
! PURPOSE: MPI commsize
!-------------------------------------------------------------------------------
#if (! defined HIDE_MPI)
call MPI_COMM_SIZE(comm,size,ierr)
if (present(string)) then
call shr_mpi_chkerr
(ierr,subName//trim(string))
else
call shr_mpi_chkerr
(ierr,subName)
endif
#endif
END SUBROUTINE shr_mpi_commsize
!===============================================================================
!===============================================================================
SUBROUTINE shr_mpi_commrank(comm,rank,string),2
IMPLICIT none
!----- arguments ---
integer,intent(in) :: comm
integer,intent(out) :: rank
character(*),optional,intent(in) :: string ! message
!----- local ---
character(*),parameter :: subName = '(shr_mpi_commrank) '
integer(SHR_KIND_IN) :: ierr
!-------------------------------------------------------------------------------
! PURPOSE: MPI commrank
!-------------------------------------------------------------------------------
#if (! defined HIDE_MPI)
call MPI_COMM_RANK(comm,rank,ierr)
if (present(string)) then
call shr_mpi_chkerr
(ierr,subName//trim(string))
else
call shr_mpi_chkerr
(ierr,subName)
endif
#endif
END SUBROUTINE shr_mpi_commrank
!===============================================================================
!===============================================================================
SUBROUTINE shr_mpi_initialized(flag,string) 1,2
IMPLICIT none
!----- arguments ---
logical,intent(out) :: flag
character(*),optional,intent(in) :: string ! message
!----- local ---
character(*),parameter :: subName = '(shr_mpi_initialized) '
integer(SHR_KIND_IN) :: ierr
!-------------------------------------------------------------------------------
! PURPOSE: MPI initialized
!-------------------------------------------------------------------------------
#if (! defined HIDE_MPI)
call MPI_INITIALIZED(flag,ierr)
if (present(string)) then
call shr_mpi_chkerr
(ierr,subName//trim(string))
else
call shr_mpi_chkerr
(ierr,subName)
endif
#endif
END SUBROUTINE shr_mpi_initialized
!===============================================================================
!===============================================================================
SUBROUTINE shr_mpi_abort(string,rcode) 19
IMPLICIT none
!----- arguments ---
character(*),optional,intent(in) :: string ! message
integer,optional,intent(in) :: rcode ! optional code
!----- local ---
character(*),parameter :: subName = '(shr_mpi_abort) '
integer(SHR_KIND_IN) :: ierr
!-------------------------------------------------------------------------------
! PURPOSE: MPI abort
!-------------------------------------------------------------------------------
#if (! defined HIDE_MPI)
write(6,*) trim(subName),":",trim(string),rcode
call MPI_ABORT(MPI_COMM_WORLD,rcode,ierr)
#endif
END SUBROUTINE shr_mpi_abort
!===============================================================================
!===============================================================================
SUBROUTINE shr_mpi_barrier(comm,string),2
IMPLICIT none
!----- arguments ---
integer,intent(in) :: comm
character(*),optional,intent(in) :: string ! message
!----- local ---
character(*),parameter :: subName = '(shr_mpi_barrier) '
integer(SHR_KIND_IN) :: ierr
!-------------------------------------------------------------------------------
! PURPOSE: MPI barrier
!-------------------------------------------------------------------------------
#if (! defined HIDE_MPI)
call MPI_BARRIER(comm,ierr)
if (present(string)) then
call shr_mpi_chkerr
(ierr,subName//trim(string))
else
call shr_mpi_chkerr
(ierr,subName)
endif
#endif
END SUBROUTINE shr_mpi_barrier
!===============================================================================
!===============================================================================
SUBROUTINE shr_mpi_init(string),2
IMPLICIT none
!----- arguments ---
character(*),optional,intent(in) :: string ! message
!----- local ---
character(*),parameter :: subName = '(shr_mpi_init) '
integer(SHR_KIND_IN) :: ierr
!-------------------------------------------------------------------------------
! PURPOSE: MPI init
!-------------------------------------------------------------------------------
#if (! defined HIDE_MPI)
call MPI_INIT(ierr)
if (present(string)) then
call shr_mpi_chkerr
(ierr,subName//trim(string))
else
call shr_mpi_chkerr
(ierr,subName)
endif
#endif
END SUBROUTINE shr_mpi_init
!===============================================================================
!===============================================================================
SUBROUTINE shr_mpi_finalize(string),2
IMPLICIT none
!----- arguments ---
character(*),optional,intent(in) :: string ! message
!----- local ---
character(*),parameter :: subName = '(shr_mpi_finalize) '
integer(SHR_KIND_IN) :: ierr
!-------------------------------------------------------------------------------
! PURPOSE: MPI finalize
!-------------------------------------------------------------------------------
#if (! defined HIDE_MPI)
call MPI_FINALIZE(ierr)
if (present(string)) then
call shr_mpi_chkerr
(ierr,subName//trim(string))
else
call shr_mpi_chkerr
(ierr,subName)
endif
#endif
END SUBROUTINE shr_mpi_finalize
!===============================================================================
!===============================================================================
END MODULE shr_mpi_mod