! #define TRACEBACK write (gol,'("in ",a," (",a,", line",i5,")")') rname, __FILE__, __LINE__; call goErr #define IF_NOTOK_RETURN(action) if (status/=0) then; TRACEBACK; action; return; end if #define IF_NOTOK_MPI(action) if (ierr/=MPI_SUCCESS) then; TRACEBACK; action; return; end if #define IF_ERROR_RETURN(action) if (status> 0) then; TRACEBACK; action; return; end if ! #include "tm5.inc" ! !----------------------------------------------------------------------------- ! TM5 ! !----------------------------------------------------------------------------- !BOP ! ! !MODULE: PARTOOLS ! ! !DESCRIPTION: MPI general tools. This module: ! ! - holds MPI constants, specifically defining some if we are not running with MPI ! ! - initializes, finalizes and aborts MPI (wrappers, so it works if not runnning MPI): ! ! TM5_MPI_Init, TM5_MPI_Init2 (init row/col communicators) ! TM5_MPI_Done ! TM5_MPI_Abort ! ! - provides wrappers around some of MPI broadcast/barrier/reduce calls : ! ! PAR_BARRIER : ! PAR_STOPMPI : ! PAR_BROADCAST_STATUS : ! PAR_BROADCAST : ! PAR_REDUCE : reduction is done across the array. Result size = 1 ! PAR_REDUCE_ELEMENT : reduction is done for each element. Result size = size of original array ! ! These are general tools, that deal with MPI, regardless of the domain decomposition used. ! ! - Some of these call can be limited to the 'row' or 'col' communicators, but ! this is partially implemented (although trivial to extent if needed) ! !\\ !\\ ! !INTERFACE: ! MODULE PARTOOLS ! ! !USES: ! use GO, only : gol, goPr, goErr use dims, only : nregions_all #ifdef MPI use mpi #endif IMPLICIT NONE ! ! !PUBLIC DATA MEMBERS: ! integer, parameter :: PAR_OPER_SUM = 100 ! mpi reduce flag integer :: localComm ! global communicator (equal to MPI_COMM_WORLD if not coupled model) integer :: row_comm ! row communicator integer :: col_comm ! column communicator integer :: my_real ! platform dependent reference to real values for MPI integer :: myid ! rank in localComm integer :: npes ! total number of PE's integer :: npe_lat ! #pes in lat direction integer :: npe_lon ! #pes in long direction integer :: ierr ! return status of MPI routine calls integer :: root ! myid of root in localComm integer :: RowRootId ! myid of row_root of my row (in localComm) #ifndef MPI integer, parameter :: MPI_SUCCESS = 0 integer, parameter :: MPI_INFO_NULL = 0 integer, parameter :: MPI_CHARACTER = 0 integer, parameter :: MPI_INTEGER = 0 integer, parameter :: MPI_PROC_NULL = -999 #endif logical :: isRoot, isRowRoot ! is process global root, row root? character(len=6) :: procname ! character keys for each processor ! ! !PRIVATE DATA MEMBERS: ! integer, private :: irow, jcol ! rank in column and row communicators character(len=*), parameter, private :: mname='ParTools' ! ! !INTERFACE: ! INTERFACE Par_Broadcast MODULE PROCEDURE Par_Broadcast_i0 MODULE PROCEDURE Par_Broadcast_i1 MODULE PROCEDURE Par_Broadcast_s MODULE PROCEDURE Par_Broadcast_l0 MODULE PROCEDURE Par_Broadcast_l1 MODULE PROCEDURE Par_Broadcast_r0 MODULE PROCEDURE Par_Broadcast_r1 MODULE PROCEDURE Par_Broadcast_r2 MODULE PROCEDURE Par_Broadcast_r3 END INTERFACE INTERFACE Par_Reduce MODULE PROCEDURE Par_Reduce_i0 MODULE PROCEDURE Par_Reduce_r0 MODULE PROCEDURE Par_Reduce_r1 ! r1, r2, ... are useful only to get results that do not depend ! MODULE PROCEDURE Par_Reduce_r2 ! on the number of processors (when summing). It is cheaper to MODULE PROCEDURE Par_Reduce_r3 ! do the operation on each proc, and then reduce a single value. END INTERFACE INTERFACE Par_Reduce_Element MODULE PROCEDURE Par_Reduce_element_r1 MODULE PROCEDURE Par_Reduce_element_r2 MODULE PROCEDURE Par_Reduce_element_r3 END INTERFACE ! ! !REVISION HISTORY: ! 18 Jan 2012 - P. Le Sager - revamped for lon-lat MPI domain decomposition ! ! !REMARKS: ! !EOP !------------------------------------------------------------------------ CONTAINS !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: TM5_MPI_Init ! ! !DESCRIPTION: initializes MPI !\\ !\\ ! !INTERFACE: ! SUBROUTINE TM5_MPI_Init( status, comm ) ! ! !OUTPUT PARAMETERS: ! integer, intent(out) :: status ! ! !INPUT PARAMETERS: ! integer, intent(in), optional :: comm ! ! !REVISION HISTORY: ! !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/TM5_MPI_Init' #ifdef MPI ! communicator provided, for example by prism coupler ? if ( present(comm) ) then localComm = comm else ! init mpi here to set MPI_COMM_WORLD etc call MPI_INIT( ierr ) IF_NOTOK_MPI(status=1) localComm = MPI_COMM_WORLD end if ! obtain number of proceses: call MPI_COMM_SIZE( localComm, npes, ierr ) IF_NOTOK_MPI(status=1) ! obtain proces number: call MPI_COMM_RANK( localComm, myid, ierr ) IF_NOTOK_MPI(status=1) ! set root in localComm to PE 0; 'real' default to double precision root = 0 my_real = MPI_DOUBLE_PRECISION #else localComm = 0 ! dummy communicator npes = 1 ! single processor root = 0 ! id for root processor myid = root ! only one processor, so this is always root #endif ! set processor names: pe0000, pe0001, ... write (procname,'("pe",i4.4)') myid isRoot = ( myid == root ) ! ok status = 0 END SUBROUTINE TM5_MPI_Init !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: TM5_MPI_INIT2 ! ! !DESCRIPTION: Step #2 of initialization, which requires input from rcfile. ! It defines the extra row and column communicators. !\\ !\\ ! !INTERFACE: ! SUBROUTINE TM5_MPI_Init2( nlon, nlat, status ) ! ! !INPUT PARAMETERS: ! integer, intent(in) :: nlat, nlon ! # PE on each direction ! ! !OUTPUT PARAMETERS: ! integer, intent(out) :: status ! ! !REVISION HISTORY: ! 21 Feb 2012 - P. Le Sager - ! !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/TM5_MPI_Init2' integer :: irow, jcol npe_lat=nlat npe_lon=nlon #ifdef MPI irow = myid/npe_lon jcol = mod(myid,npe_lon) ! NOTE: could use myid as key for rank designation if that facilitates some comm call MPI_COMM_SPLIT(localComm, irow, jcol, row_comm, ierr) IF_NOTOK_MPI(status=1) isRowRoot = jcol==0 RowRootId = (myid/npe_lon)*npe_lon ! define column communicator [not needed yet] call MPI_COMM_SPLIT(localComm, jcol, irow, col_comm, ierr) IF_NOTOK_MPI(status=1) #else row_comm = 0 ! dummy communicator col_comm = 0 ! dummy communicator #endif status = 0 END SUBROUTINE TM5_MPI_Init2 !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: TM5_MPI_Done ! ! !DESCRIPTION: finalizes MPI !\\ !\\ ! !INTERFACE: ! SUBROUTINE TM5_MPI_Done( status, comm ) ! ! !OUTPUT PARAMETERS: ! integer, intent(out) :: status ! ! !INPUT PARAMETERS: ! integer, intent(in), optional :: comm ! ! !REVISION HISTORY: ! !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/TM5_MPI_Done' #ifdef MPI ! shut-down communication, only if communicator is not provided if (.not. present(comm) ) then !write (*,'("call MPI_Finalize from proces ",i6)') myid call MPI_Finalize( ierr ) IF_NOTOK_MPI(status=1) end if #endif ! ok status = 0 END SUBROUTINE TM5_MPI_Done !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: TM5_MPI_Abort ! ! !DESCRIPTION: !\\ !\\ ! !INTERFACE: ! subroutine TM5_MPI_Abort( errorcode, status ) ! ! !USES: ! use GO, only : goExit ! ! !INPUT PARAMETERS: ! integer, intent(in) :: errorcode ! ! !OUTPUT PARAMETERS: ! integer, intent(out) :: status ! ! !REVISION HISTORY: ! 18 Jan 2012 - P. Le Sager - ! ! !REMARKS: ! ! (pls, 8-4-2011) Sometimes the code does not return from MPI_Abort, for ! example when a problem reading restart files occurs. From the doc: ! ------------------------------------------------- ! "Before the error value is returned, the current MPI error handler is ! called. By default, this error handler aborts the MPI job, except for ! I/O function errors." ! ------------------------------------------------- ! so, the only way to nicely abort is to close files when an i/o error ! occurs. Done with a new macro for problems with reading/writing restart ! (see tm5_restart.F90). Check if your module/routine is prone to ! i/o error, and apply a similar patch. ! !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/TM5_MPI_Abort' #ifdef MPI ! emergency break ... call MPI_Abort( localComm, errorcode, ierr ) IF_NOTOK_MPI(status=1) #else ! system exit: call goExit( errorcode ) #endif ! ok status = 0 end subroutine TM5_MPI_Abort !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: Par_Barrier ! ! !DESCRIPTION: !\\ !\\ ! !INTERFACE: ! SUBROUTINE PAR_BARRIER( ROW, COL) ! ! !INPUT PARAMETERS: ! logical, intent(in), optional :: row, col ! to limit call to sub-communicator ! ! !REVISION HISTORY: ! !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/PAR_BARRIER' #ifdef MPI integer :: l_comm l_comm=localComm if(present(row)) then if(row) l_comm=row_comm end if if(present(col)) then if(col) l_comm=col_comm end if call mpi_barrier( l_comm, ierr ) IF_NOTOK_MPI(write (gol,*)"error MPI_BARRIER in PAR_BARRIER";call goErr) #endif END SUBROUTINE PAR_BARRIER !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: Par_StopMPI ! ! !DESCRIPTION: !\\ !\\ ! !INTERFACE: ! subroutine Par_StopMPI ! ! !REVISION HISTORY: ! !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/Par_StopMPI' write (*,'("WARNING - (par_)stopmpi should be avoided; please trace back to main program ...")') #ifdef MPI ! shut down mpi communication: call mpi_finalize(ierr) IF_NOTOK_MPI(write (gol,*)"error MPI_FINALIZE in Par_StopMPI";call goErr) #endif ! fortran stop .... stop 'Fortran STOP in Par_StopMPI ...' end subroutine Par_StopMPI !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: Par_Broadcast_Status ! ! !DESCRIPTION: Broadcast integer istat from id PE to all. Same as ! Par_Broadcast_i(istat, id, istat) !\\ !\\ ! !INTERFACE: ! subroutine Par_Broadcast_Status( istat, id ) ! ! !INPUT/OUTPUT PARAMETERS: ! integer, intent(inout) :: istat ! ! !INPUT PARAMETERS: ! integer, intent(in) :: id ! ! !REVISION HISTORY: ! ! !REMARKS: ! !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/Par_Broadcast_Status' integer :: status #ifdef MPI ! send the input status to all other processes: call Par_Broadcast( istat, status, id=id ) IF_NOTOK_RETURN(istat=1) #endif end subroutine Par_Broadcast_Status !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: Par_Broadcast_i0 ! ! !DESCRIPTION: broadcast one integer (scalar case) !\\ !\\ ! !INTERFACE: ! subroutine Par_Broadcast_i0( i, status, ID, ROW, COL ) ! ! !INPUT/OUTPUT PARAMETERS: ! integer, intent(inout) :: i integer, intent(out) :: status integer, intent(in), optional :: id ! broadcaster ID (default is communicator root) logical, intent(in), optional :: row ! limit to PE on row_comm (default is global localComm) logical, intent(in), optional :: col ! limit to PE on col_comm (default is global localComm) ! ! !REVISION HISTORY: ! 26 Mar 2012 - P. Le Sager - added ROW, COL options. Made ID optional. ! !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/Par_Broadcast_i' logical :: l_row, l_col integer :: l_root, l_comm, l_id #ifdef MPI l_row=.false. ; if(present(row)) l_row=row l_col=.false. ; if(present(col)) l_col=col if(l_row)then l_root = 0 ! by our own design l_comm = row_comm else if(l_col) then l_root = 0 ! by our own design l_comm = col_comm else l_root = root l_comm = localComm end if l_id = l_root if(present(id)) l_id=id call MPI_BCast( i, 1, MPI_INTEGER, l_id, l_comm, ierr ) IF_NOTOK_MPI(status=1) #endif status = 0 end subroutine Par_Broadcast_i0 !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: Par_Broadcast_i1 ! ! !DESCRIPTION: Broadcast 1-D array of integers. !\\ !\\ ! !INTERFACE: ! subroutine Par_Broadcast_i1( i, status, ID, ROW, COL ) ! ! !INPUT/OUTPUT PARAMETERS: ! integer, intent(inout) :: i(:) integer, intent(out) :: status integer, intent(in), optional :: id ! broadcaster ID (default is communicator root) logical, intent(in), optional :: row ! limit to PE on row_comm (default is global localComm) logical, intent(in), optional :: col ! limit to PE on col_comm (default is global localComm) ! !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/Par_Broadcast_i' logical :: l_row, l_col integer :: l_root, l_comm, l_id #ifdef MPI l_row=.false. ; if(present(row)) l_row=row l_col=.false. ; if(present(col)) l_col=col if(l_row)then l_root = 0 ! by our own design l_comm = row_comm else if(l_col) then l_root = 0 ! by our own design l_comm = col_comm else l_root = root l_comm = localComm end if l_id = l_root if(present(id)) l_id=id call MPI_BCast( i, size(i), MPI_INTEGER, l_id, l_comm, ierr ) IF_NOTOK_MPI(status=1) #endif status = 0 end subroutine Par_Broadcast_i1 !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: Par_Broadcast_s ! ! !DESCRIPTION: !\\ !\\ ! !INTERFACE: ! subroutine Par_Broadcast_s( s, status, ID, ROW, COL ) ! ! !INPUT/OUTPUT PARAMETERS: ! character(len=*), intent(inout) :: s integer, intent(out) :: status integer, intent(in), optional :: id ! broadcaster ID (default is communicator root) logical, intent(in), optional :: row ! limit to PE on row_comm (default is global localComm) logical, intent(in), optional :: col ! limit to PE on col_comm (default is global localComm) ! ! !REVISION HISTORY: ! 26 Mar 2012 - P. Le Sager - added ROW, COL options. Made ID optional. ! !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/Par_Broadcast_s' logical :: l_row, l_col integer :: l_root, l_comm, l_id #ifdef MPI l_row=.false. ; if(present(row)) l_row=row l_col=.false. ; if(present(col)) l_col=col if(l_row)then l_root = 0 ! by our own design l_comm = row_comm else if(l_col) then l_root = 0 ! by our own design l_comm = col_comm else l_root = root l_comm = localComm end if l_id = l_root if(present(id)) l_id=id call MPI_BCast( s, len(s), MPI_CHARACTER, l_id, l_comm, ierr ) IF_NOTOK_MPI(status=1) #endif status = 0 end subroutine Par_Broadcast_s !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: Par_Broadcast_l0 ! ! !DESCRIPTION: !\\ !\\ ! !INTERFACE: ! subroutine Par_Broadcast_l0( x, status, ID, ROW, COL ) ! ! !INPUT/OUTPUT PARAMETERS: ! logical, intent(inout) :: x integer, intent(out) :: status integer, intent(in), optional :: id ! broadcaster ID (default is communicator root) logical, intent(in), optional :: row ! limit to PE on row_comm (default is global localComm) logical, intent(in), optional :: col ! limit to PE on col_comm (default is global localComm) ! ! !REVISION HISTORY: ! 26 Mar 2012 - P. Le Sager - added ROW, COL options. Made ID optional. ! !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/Par_Broadcast_l0' logical :: l_row, l_col integer :: l_root, l_comm, l_id #ifdef MPI l_row=.false. ; if(present(row)) l_row=row l_col=.false. ; if(present(col)) l_col=col if(l_row)then l_root = 0 ! by our own design l_comm = row_comm else if(l_col) then l_root = 0 ! by our own design l_comm = col_comm else l_root = root l_comm = localComm end if l_id = l_root if(present(id)) l_id=id call mpi_bcast( x, 1, MPI_LOGICAL, l_id, l_comm, ierr ) IF_NOTOK_MPI(status=1) #endif status = 0 end subroutine Par_Broadcast_l0 !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: Par_Broadcast_l1 ! ! !DESCRIPTION: !\\ !\\ ! !INTERFACE: ! subroutine Par_Broadcast_l1( x, status, ID, ROW, COL ) ! ! !INPUT/OUTPUT PARAMETERS: ! logical, intent(inout) :: x(:) integer, intent(out) :: status integer, intent(in), optional :: id ! broadcaster ID (default is communicator root) logical, intent(in), optional :: row ! limit to PE on row_comm (default is global localComm) logical, intent(in), optional :: col ! limit to PE on col_comm (default is global localComm) ! ! !REVISION HISTORY: ! 26 Mar 2012 - P. Le Sager - added ROW, COL options. Made ID optional. ! !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/Par_Broadcast_l1' logical :: l_row, l_col integer :: l_root, l_comm, l_id #ifdef MPI l_row=.false. ; if(present(row)) l_row=row l_col=.false. ; if(present(col)) l_col=col if(l_row)then l_root = 0 ! by our own design l_comm = row_comm else if(l_col) then l_root = 0 ! by our own design l_comm = col_comm else l_root = root l_comm = localComm end if l_id = l_root if(present(id)) l_id=id call mpi_bcast( x, size(x), MPI_LOGICAL, l_id, l_comm, ierr ) IF_NOTOK_MPI(status=1) #endif status = 0 end subroutine Par_Broadcast_l1 !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: Par_Broadcast_r0 ! ! !DESCRIPTION: !\\ !\\ ! !INTERFACE: ! subroutine Par_Broadcast_r0( x, status, ID, ROW, COL ) ! ! !INPUT/OUTPUT PARAMETERS: ! real, intent(inout) :: x integer, intent(out) :: status integer, intent(in), optional :: id ! broadcaster ID (default is communicator root) logical, intent(in), optional :: row ! limit to PE on row_comm (default is global localComm) logical, intent(in), optional :: col ! limit to PE on col_comm (default is global localComm) ! ! !REVISION HISTORY: ! 26 Mar 2012 - P. Le Sager - added ROW, COL options. Made ID optional. ! !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/Par_Broadcast_r0' logical :: l_row, l_col integer :: l_root, l_comm, l_id #ifdef MPI l_row=.false. ; if(present(row)) l_row=row l_col=.false. ; if(present(col)) l_col=col if(l_row)then l_root = 0 ! by our own design l_comm = row_comm else if(l_col) then l_root = 0 ! by our own design l_comm = col_comm else l_root = root l_comm = localComm end if l_id = l_root if(present(id)) l_id=id call mpi_bcast( x, 1, my_real, l_id, l_comm, ierr ) IF_NOTOK_MPI(status=1) #endif status = 0 end subroutine Par_Broadcast_r0 !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: Par_Broadcast_r1 ! ! !DESCRIPTION: !\\ !\\ ! !INTERFACE: ! subroutine Par_Broadcast_r1( x, status, ID, ROW, COL ) ! ! !INPUT/OUTPUT PARAMETERS: ! real, intent(inout) :: x(:) integer, intent(out) :: status integer, intent(in), optional :: id ! broadcaster ID (default is communicator root) logical, intent(in), optional :: row ! limit to PE on row_comm (default is global localComm) logical, intent(in), optional :: col ! limit to PE on col_comm (default is global localComm) ! ! !REVISION HISTORY: ! 26 Mar 2012 - P. Le Sager - added ROW, COL options. Made ID optional. ! !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/Par_Broadcast_r1' logical :: l_row, l_col integer :: l_root, l_comm, l_id #ifdef MPI l_row=.false. ; if(present(row)) l_row=row l_col=.false. ; if(present(col)) l_col=col if(l_row)then l_root = 0 ! by our own design l_comm = row_comm else if(l_col) then l_root = 0 ! by our own design l_comm = col_comm else l_root = root l_comm = localComm end if l_id = l_root if(present(id)) l_id=id call mpi_bcast( x, size(x), my_real, l_id, l_comm, ierr ) IF_NOTOK_MPI(status=1) #endif status = 0 end subroutine Par_Broadcast_r1 !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: Par_Broadcast_r2 ! ! !DESCRIPTION: !\\ !\\ ! !INTERFACE: ! subroutine Par_Broadcast_r2( x, status, ID, ROW, COL ) ! ! !INPUT/OUTPUT PARAMETERS: ! real, intent(inout) :: x(:,:) integer, intent(out) :: status integer, intent(in), optional :: id ! broadcaster ID (default is communicator root) logical, intent(in), optional :: row ! limit to PE on row_comm (default is global localComm) logical, intent(in), optional :: col ! limit to PE on col_comm (default is global localComm) ! ! !REVISION HISTORY: ! 26 Mar 2012 - P. Le Sager - added ROW, COL options. Made ID optional. ! !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/Par_Broadcast_r2' logical :: l_row, l_col integer :: l_root, l_comm, l_id #ifdef MPI l_row=.false. ; if(present(row)) l_row=row l_col=.false. ; if(present(col)) l_col=col if(l_row)then l_root = 0 ! by our own design l_comm = row_comm else if(l_col) then l_root = 0 ! by our own design l_comm = col_comm else l_root = root l_comm = localComm end if l_id = l_root if(present(id)) l_id=id call mpi_bcast( x, size(x), my_real, l_id, l_comm, ierr ) IF_NOTOK_MPI(status=1) #endif status = 0 end subroutine Par_Broadcast_r2 !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: Par_Broadcast_r3 ! ! !DESCRIPTION: !\\ !\\ ! !INTERFACE: ! subroutine Par_Broadcast_r3( x, status, ID, ROW, COL ) ! ! !INPUT/OUTPUT PARAMETERS: ! real, intent(inout) :: x(:,:,:) integer, intent(out) :: status integer, intent(in), optional :: id ! broadcaster ID (default is communicator root) logical, intent(in), optional :: row ! limit to PE on row_comm (default is global localComm) logical, intent(in), optional :: col ! limit to PE on col_comm (default is global localComm) ! ! !REVISION HISTORY: ! 26 Mar 2012 - P. Le Sager - added ROW, COL options. Made ID optional. ! !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/Par_Broadcast_r3' logical :: l_row, l_col integer :: l_root, l_comm, l_id #ifdef MPI l_row=.false. ; if(present(row)) l_row=row l_col=.false. ; if(present(col)) l_col=col if(l_row)then l_root = 0 ! by our own design l_comm = row_comm else if(l_col) then l_root = 0 ! by our own design l_comm = col_comm else l_root = root l_comm = localComm end if l_id = l_root if(present(id)) l_id=id call mpi_bcast( x, size(x), my_real, l_id, l_comm, ierr ) IF_NOTOK_MPI(status=1) #endif status = 0 end subroutine Par_Broadcast_r3 !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: PAR_REDUCE_I0 ! ! !DESCRIPTION: Wrapper around MPI_REDUCE or MPI_ALLREDUCE. ! ! Apply to a SINGLE INTEGER. !\\ !\\ ! !INTERFACE: ! SUBROUTINE PAR_REDUCE_I0( DATA, OP, RESULTAT, STATUS, ALL, ROW ) ! ! !INPUT PARAMETERS: ! integer, intent(in) :: data character(len=3), intent(in) :: op ! 'MAX', 'MIN' or 'SUM' logical, intent(in), optional :: all ! use MPI_ALLREDUCE instead of MPI_REDUCE logical, intent(in), optional :: row ! limit to PE on row_comm ! ! !OUTPUT PARAMETERS: ! integer, intent(out) :: resultat integer, intent(out) :: status ! ! !REVISION HISTORY: ! 19 Jun 2013 - Ph. Le Sager - v0 (copy from par_reduce_r0) ! !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/PAR_REDUCE_I0' logical :: l_row, l_all integer :: l_root, l_comm, l_id, l_op #ifdef MPI l_row=.false. ; if(present(row)) l_row=row l_all=.false. ; if(present(all)) l_all=all if(l_row)then l_root = 0 ! by our own design l_comm = row_comm l_id = jcol else l_root = root l_comm = localComm l_id = myid end if ! degenerate cases first if(l_row.and.(npe_lon==1)) then resultat = data else SELECT CASE( OP ) case('sum', 'Sum', 'SUM') l_op = MPI_SUM case('max', 'Max', 'MAX') l_op = MPI_MAX case('min', 'Min', 'MIN') l_op = MPI_MIN case default write(gol,*) 'UNSUPPORTED OPERATION :', OP; status=1 IF_NOTOK_RETURN(status=1) END SELECT if (l_all) then call MPI_ALLREDUCE(data, resultat, 1, MPI_INTEGER, l_op, l_comm, ierr) IF_NOTOK_MPI(status=1) else call MPI_REDUCE(data, resultat, 1, MPI_INTEGER, l_op, l_root, l_comm, ierr) IF_NOTOK_MPI(status=1) end if end if #else resultat = data #endif status=0 END SUBROUTINE PAR_REDUCE_I0 !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: PAR_REDUCE_R0 ! ! !DESCRIPTION: Wrapper around MPI_REDUCE or MPI_ALLREDUCE. ! ! Apply to a SINGLE REAL. !\\ !\\ ! !INTERFACE: ! SUBROUTINE PAR_REDUCE_R0( DATA, OP, RESULTAT, STATUS, ALL, ROW, COL ) ! ! !INPUT PARAMETERS: ! real, intent(in) :: data character(len=3), intent(in) :: op ! 'MAX', 'MIN' or 'SUM' logical, intent(in), optional :: all ! use MPI_ALLREDUCE instead of MPI_REDUCE logical, intent(in), optional :: row ! limit to PE on row_comm logical, intent(in), optional :: col ! limit to PE on col_comm ! ! !OUTPUT PARAMETERS: ! real, intent(out) :: resultat integer, intent(out) :: status ! ! !REVISION HISTORY: ! 01 Nov 2011 - P. Le Sager - v0 ! !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/PAR_REDUCE_R0' logical :: l_row, l_all, l_col integer :: l_root, l_comm, l_id, l_op #ifdef MPI l_col=.false. ; if(present(col)) l_col=col l_row=.false. ; if(present(row)) l_row=row l_all=.false. ; if(present(all)) l_all=all if(l_row)then l_root = 0 ! by our own design l_comm = row_comm l_id = jcol else if(l_col) then l_root = 0 ! by our own design l_comm = col_comm l_id = irow else l_root = root l_comm = localComm l_id = myid end if ! degenerate cases first if(l_row.and.(npe_lon==1)) then resultat = data else SELECT CASE( OP ) case('sum', 'Sum', 'SUM') l_op = MPI_SUM case('max', 'Max', 'MAX') l_op = MPI_MAX case('min', 'Min', 'MIN') l_op = MPI_MIN case default write(gol,*) 'UNSUPPORTED OPERATION :', OP; status=1 IF_NOTOK_RETURN(status=1) END SELECT if (l_all) then call MPI_ALLREDUCE(data, resultat, 1, my_real, l_op, l_comm, ierr) IF_NOTOK_MPI(status=1) else call MPI_REDUCE(data, resultat, 1, my_real, l_op, l_root, l_comm, ierr) IF_NOTOK_MPI(status=1) end if end if #else resultat = data #endif status=0 END SUBROUTINE PAR_REDUCE_R0 !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: PAR_REDUCE_R1 ! ! !DESCRIPTION: Global reduction. Data are gathered on root, where OP is ! then done. Apply to a 1D REAL. !\\ !\\ ! !INTERFACE: ! SUBROUTINE PAR_REDUCE_R1( DATA, OP, RESULTAT, STATUS, ALL, ROW ) ! ! !INPUT PARAMETERS: ! real, intent(in) :: data(:) character(len=3), intent(in) :: op ! 'MAX', 'MIN' or 'SUM' logical, intent(in), optional :: all ! mimic mpi_allreduce instead of mpi_reduce logical, intent(in), optional :: row ! limit to PE on row_comm ! ! !OUTPUT PARAMETERS: ! real, intent(out) :: resultat integer, intent(out) :: status ! ! !REVISION HISTORY: ! 01 Nov 2011 - P. Le Sager - FIXME: NOT TESTED ! ! !REMARKS: ! ! (1) this is a convenient tool to get same results independently of the ! number of processors. If this is not what you want, you can do OP ! on each processor before calling par_reduce_r0 : it is less expensive. ! !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/PAR_REDUCE_R1' #ifdef MPI real, allocatable :: glb(:) logical :: l_row integer :: sz, l_root, l_comm, l_id sz = size(data) l_row=.false. if(present(row)) l_row=row if(l_row)then allocate( glb(npe_lon*sz) ) l_root = 0 l_comm = row_comm l_id = jcol else allocate( glb(npes*sz) ) l_root = root l_comm = localComm l_id = myid end if CALL MPI_GATHER( data, sz, my_real, glb, sz, my_real, l_root, l_comm, ierr) IF_NOTOK_MPI(status=1) if (l_root==l_id) then select case( op ) case('sum', 'Sum', 'SUM') resultat = sum(glb) case('max', 'Max', 'MAX') resultat = maxval(glb) case('min', 'Min', 'MIN') resultat = minval(glb) case default write(gol,*) 'UNSUPPORTED OPERATION'; status=1 IF_NOTOK_RETURN(status=1) end select end if deallocate(glb) if (present(all)) then if (all) call MPI_bcast(resultat, 1, my_real, l_root, l_comm, ierr) end if #else select case( op ) case('sum', 'Sum', 'SUM') resultat = sum(data) case('max', 'Max', 'MAX') resultat = maxval(data) case('min', 'Min', 'MIN') resultat = minval(data) case default write(gol,*) 'UNSUPPORTED OPERATION'; status=1 IF_NOTOK_RETURN(status=1) end select #endif status=0 END SUBROUTINE PAR_REDUCE_R1 !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: PAR_REDUCE_R2 ! ! !DESCRIPTION: !\\ !\\ ! !INTERFACE: ! ! subroutine Par_Reduce_r2( send, recv, oper, id_recv, status ) ! ! !REMARKS: ! in TM5 v3, this routine was a wrapper around MPI_REDUCE( ..., MPI_SUM, ...), and was called ! in project CLIMAQS/user_output_pdump.F90 (an updated retro output) and in ! various user_output_retro.F90, to sum budget... ! =====> in tm5-mp, this routine became PAR_REDUCE_ELEMENT_R2 () ! TO DEVELOP IF YOU NEED RESULT THAT ARE INDEPENDENT OF THE NUMBER OF ! PROCESSOR YOU ARE USING. ! end subroutine Par_Reduce_r2 !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: PAR_REDUCE_R3 ! ! !DESCRIPTION: Global reduction. Data are gathered on root, where OP is ! then done. Apply to a 3D REAL. !\\ !\\ ! !INTERFACE: ! SUBROUTINE PAR_REDUCE_R3( DATA, OP, RESULTAT, STATUS, ALL, ROW ) ! ! !INPUT PARAMETERS: ! real, intent(in) :: data(:,:,:) character(len=3), intent(in) :: op ! 'MAX', 'MIN' or 'SUM' logical, intent(in), optional :: all ! mimic mpi_allreduce instead of mpi_reduce logical, intent(in), optional :: row ! limit to PE on row_comm ! ! !OUTPUT PARAMETERS: ! real, intent(out) :: resultat integer, intent(out) :: status ! ! !REVISION HISTORY: ! ! !REMARKS: ! ! (1) this is a convenient tool to get same results independently of the ! number of processors. If this is not what you want, you can do OP ! on each processor before calling par_reduce_r0 : it is less expensive. ! !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/PAR_REDUCE_R3' #ifdef MPI real, allocatable :: glb(:) logical :: l_row integer :: sz, l_root, l_comm, l_id sz = size(data) l_row=.false. if(present(row)) l_row=row if(l_row)then allocate( glb(npe_lon*sz) ) l_root = 0 l_comm = row_comm l_id = jcol else allocate( glb(npes*sz) ) l_root = root l_comm = localComm l_id = myid end if CALL MPI_GATHER( data, sz, my_real, glb, sz, my_real, l_root, l_comm, ierr) IF_NOTOK_MPI(status=1) if (l_root==l_id) then select case( op ) case('sum', 'Sum', 'SUM') resultat = sum(glb) case('max', 'Max', 'MAX') resultat = maxval(glb) case('min', 'Min', 'MIN') resultat = minval(glb) case default write(gol,*) 'UNSUPPORTED OPERATION'; status=1 IF_NOTOK_RETURN(status=1) end select end if deallocate(glb) if (present(all)) then if (all) call MPI_bcast(resultat, 1, my_real, l_root, l_comm, ierr) end if #else select case( op ) case('sum', 'Sum', 'SUM') resultat = sum(data) case('max', 'Max', 'MAX') resultat = maxval(data) case('min', 'Min', 'MIN') resultat = minval(data) case default write(gol,*) 'UNSUPPORTED OPERATION'; status=1 IF_NOTOK_RETURN(status=1) end select #endif status=0 END SUBROUTINE PAR_REDUCE_R3 !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: PAR_REDUCE_ELEMENT_R1 ! ! !DESCRIPTION: reduce 1D arrays element-wise !\\ !\\ ! !INTERFACE: ! SUBROUTINE PAR_REDUCE_ELEMENT_R1 (DATA, OP, RESULTAT, STATUS, ALL, ROW) ! ! !USES: ! use dims, only : CheckShape ! ! !INPUT PARAMETERS: ! real, intent(in) :: data(:) character(len=3), intent(in) :: op ! 'MAX', 'MIN' or 'SUM' logical, intent(in), optional :: all ! mimic mpi_allreduce instead of mpi_reduce logical, intent(in), optional :: row ! limit to PE on row_comm ! ! !OUTPUT PARAMETERS: ! real, intent(out) :: resultat(:) integer, intent(out) :: status ! ! !REVISION HISTORY: ! 1 Mar 2012 - P. Le Sager - v0 ! ! !REMARKS: ! !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/PAR_REDUCE_ELEMENT_R2' logical :: l_row, l_all integer :: l_root, l_comm, l_id, l_op call CheckShape( shape(data), shape(resultat), status ) IF_NOTOK_RETURN(status=1) #ifdef MPI l_row=.false. ; if(present(row)) l_row=row l_all=.false. ; if(present(all)) l_all=all if(l_row)then l_root = 0 ! by our own design l_comm = row_comm l_id = jcol else l_root = root l_comm = localComm l_id = myid end if ! degenerate cases first if(l_row.and.(npe_lon==1)) then resultat = data else SELECT CASE( OP ) case('sum', 'Sum', 'SUM') l_op = MPI_SUM case('max', 'Max', 'MAX') l_op = MPI_MAX case('min', 'Min', 'MIN') l_op = MPI_MIN case default write(gol,*) 'UNSUPPORTED OPERATION :', OP; status=1 IF_NOTOK_RETURN(status=1) END SELECT if (l_all) then call MPI_ALLREDUCE(data, resultat, size(data), my_real, l_op, l_comm, ierr) IF_NOTOK_MPI(status=1) else call MPI_REDUCE(data, resultat, size(data), my_real, l_op, l_root, l_comm, ierr) IF_NOTOK_MPI(status=1) end if end if #else resultat = data #endif status=0 END SUBROUTINE PAR_REDUCE_ELEMENT_R1 !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: PAR_REDUCE_ELEMENT_R2 ! ! !DESCRIPTION: reduce 2D arrays element-wise !\\ !\\ ! !INTERFACE: ! SUBROUTINE PAR_REDUCE_ELEMENT_R2 (DATA, OP, RESULTAT, STATUS, ALL, ROW) ! ! !USES: ! use dims, only : CheckShape ! ! !INPUT PARAMETERS: ! real, intent(in) :: data(:,:) character(len=3), intent(in) :: op ! 'MAX', 'MIN' or 'SUM' logical, intent(in), optional :: all ! mimic mpi_allreduce instead of mpi_reduce logical, intent(in), optional :: row ! limit to PE on row_comm ! ! !OUTPUT PARAMETERS: ! real, intent(out) :: resultat(:,:) integer, intent(out) :: status ! ! !REVISION HISTORY: ! 1 Mar 2012 - P. Le Sager - v0 ! ! !REMARKS: ! !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/PAR_REDUCE_ELEMENT_R2' logical :: l_row, l_all integer :: l_root, l_comm, l_id, l_op call CheckShape( shape(data), shape(resultat), status ) IF_NOTOK_RETURN(status=1) #ifdef MPI l_row=.false. ; if(present(row)) l_row=row l_all=.false. ; if(present(all)) l_all=all if(l_row)then l_root = 0 ! by our own design l_comm = row_comm l_id = jcol else l_root = root l_comm = localComm l_id = myid end if ! degenerate cases first if(l_row.and.(npe_lon==1)) then resultat = data else SELECT CASE( OP ) case('sum', 'Sum', 'SUM') l_op = MPI_SUM case('max', 'Max', 'MAX') l_op = MPI_MAX case('min', 'Min', 'MIN') l_op = MPI_MIN case default write(gol,*) 'UNSUPPORTED OPERATION :', OP; status=1 IF_NOTOK_RETURN(status=1) END SELECT if (l_all) then call MPI_ALLREDUCE(data, resultat, size(data), my_real, l_op, l_comm, ierr) IF_NOTOK_MPI(status=1) else call MPI_REDUCE(data, resultat, size(data), my_real, l_op, l_root, l_comm, ierr) IF_NOTOK_MPI(status=1) end if end if #else resultat = data #endif status=0 END SUBROUTINE PAR_REDUCE_ELEMENT_R2 !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: PAR_REDUCE_ELEMENT_R3 ! ! !DESCRIPTION: reduce 3D arrays element-wise !\\ !\\ ! !INTERFACE: ! SUBROUTINE PAR_REDUCE_ELEMENT_R3 (DATA, OP, RESULTAT, STATUS, ALL, ROW) ! ! !USES: ! use dims, only : CheckShape ! ! !INPUT PARAMETERS: ! real, intent(in) :: data(:,:,:) character(len=3), intent(in) :: op ! 'MAX', 'MIN' or 'SUM' logical, intent(in), optional :: all ! mimic mpi_allreduce instead of mpi_reduce logical, intent(in), optional :: row ! limit to PE on row_comm ! ! !OUTPUT PARAMETERS: ! real, intent(out) :: resultat(:,:,:) integer, intent(out) :: status ! ! !REVISION HISTORY: ! 1 Mar 2012 - P. Le Sager - v0 ! ! !REMARKS: ! !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/PAR_REDUCE_ELEMENT_R3' logical :: l_row, l_all integer :: l_root, l_comm, l_id, l_op call CheckShape( shape(data), shape(resultat), status ) IF_NOTOK_RETURN(status=1) #ifdef MPI l_row=.false. ; if(present(row)) l_row=row l_all=.false. ; if(present(all)) l_all=all if(l_row)then l_root = 0 ! by our own design l_comm = row_comm l_id = jcol else l_root = root l_comm = localComm l_id = myid end if ! degenerate cases first if(l_row.and.(npe_lon==1)) then resultat = data else SELECT CASE( OP ) case('sum', 'Sum', 'SUM') l_op = MPI_SUM case('max', 'Max', 'MAX') l_op = MPI_MAX case('min', 'Min', 'MIN') l_op = MPI_MIN case default write(gol,*) 'UNSUPPORTED OPERATION :', OP; status=1 IF_NOTOK_RETURN(status=1) END SELECT if (l_all) then call MPI_ALLREDUCE(data, resultat, size(data), my_real, l_op, l_comm, ierr) IF_NOTOK_MPI(status=1) else call MPI_REDUCE(data, resultat, size(data), my_real, l_op, l_root, l_comm, ierr) IF_NOTOK_MPI(status=1) end if end if #else resultat = data #endif status=0 END SUBROUTINE PAR_REDUCE_ELEMENT_R3 !EOC END MODULE ParTools