1!{\src2tex{textfont=tt}} 2!!****f* ABINIT/xmpi_alltoall 3!! NAME 4!! xmpi_alltoall 5!! 6!! FUNCTION 7!! This module contains functions that calls MPI routine, 8!! if we compile the code using the MPI CPP flags. 9!! xmpi_alltoall is the generic function. 10!! 11!! COPYRIGHT 12!! Copyright (C) 2001-2016 ABINIT group (AR,XG) 13!! This file is distributed under the terms of the 14!! GNU General Public License, see ~ABINIT/COPYING 15!! or http://www.gnu.org/copyleft/gpl.txt . 16!! 17!! SOURCE 18 19!!*** 20 21!!****f* ABINIT/xmpi_alltoall_int 22!! NAME 23!! xmpi_alltoall_int 24!! 25!! FUNCTION 26!! Sends data from all to all processes. 27!! Target: integer mono dimensional arrays. 28!! 29!! SOURCE 30 31subroutine xmpi_alltoall_int(xval, sendsize, recvbuf, recvsize, spaceComm, ier) 32 33 34 35!This section has been created automatically by the script Abilint (TD). 36!Do not modify the following lines by hand. 37#undef ABI_FUNC 38#define ABI_FUNC 'xmpi_alltoall_int' 39!End of the abilint section 40 41 implicit none 42 43!Arguments------------------------- 44 integer, DEV_CONTARRD intent(in) :: xval(:) 45 integer, DEV_CONTARRD intent(inout) :: recvbuf(:) 46 integer ,intent(in) :: sendsize, recvsize 47 integer ,intent(in) :: spaceComm 48 integer ,intent(out) :: ier 49 50!Local variables------------------- 51 52! ************************************************************************* 53 54 ier=0 55#if defined HAVE_MPI 56 if (spaceComm /= MPI_COMM_SELF .and. spaceComm /= MPI_COMM_NULL) then 57! allgather xval on all proc. in spaceComm 58 call MPI_ALLTOALL(xval, sendsize, MPI_INTEGER, recvbuf, & 59& recvsize, MPI_INTEGER, spaceComm, ier) 60 else if (spaceComm == MPI_COMM_SELF) then 61 recvbuf=xval 62 end if 63#else 64 recvbuf=xval 65#endif 66 67end subroutine xmpi_alltoall_int 68!!*** 69 70!!****f* ABINIT/xmpi_alltoall_dp2d 71!! NAME 72!! xmpi_alltoall_dp2d 73!! 74!! FUNCTION 75!! Sends data from all to all processes. 76!! Target: double precision two-dimensional arrays. 77!! 78!! SOURCE 79 80subroutine xmpi_alltoall_dp2d(xval, sendsize, recvbuf, recvsize, spaceComm, ier) 81 82 83 84!This section has been created automatically by the script Abilint (TD). 85!Do not modify the following lines by hand. 86#undef ABI_FUNC 87#define ABI_FUNC 'xmpi_alltoall_dp2d' 88!End of the abilint section 89 90 implicit none 91 92!Arguments------------------------- 93 real(dp), DEV_CONTARRD intent(in) :: xval(:,:) 94 real(dp), DEV_CONTARRD intent(inout) :: recvbuf(:,:) 95 integer ,intent(in) :: sendsize, recvsize 96 integer ,intent(in) :: spaceComm 97 integer ,intent(out) :: ier 98 99! ************************************************************************* 100 101 ier=0 102#if defined HAVE_MPI 103 if (spaceComm /= MPI_COMM_SELF .and. spaceComm /= MPI_COMM_NULL) then 104! allgather xval on all proc. in spaceComm 105 call MPI_ALLTOALL(xval, sendsize, MPI_DOUBLE_PRECISION, recvbuf, & 106& recvsize, MPI_DOUBLE_PRECISION, spaceComm, ier) 107 else if (spaceComm == MPI_COMM_SELF) then 108 recvbuf=xval 109 end if 110#else 111 recvbuf=xval 112#endif 113 114end subroutine xmpi_alltoall_dp2d 115!!*** 116 117!!****f* ABINIT/xmpi_alltoall_dp4d 118!! NAME 119!! xmpi_alltoall_dp4d 120!! 121!! FUNCTION 122!! Sends data from all to all processes. 123!! Target: double precision four-dimensional arrays. 124!! 125!! SOURCE 126 127subroutine xmpi_alltoall_dp4d(xval, sendsize, recvbuf, recvsize, spaceComm, ier) 128 129 130 131!This section has been created automatically by the script Abilint (TD). 132!Do not modify the following lines by hand. 133#undef ABI_FUNC 134#define ABI_FUNC 'xmpi_alltoall_dp4d' 135!End of the abilint section 136 137 implicit none 138 139!Arguments------------------------- 140 real(dp), DEV_CONTARRD intent(in) :: xval(:,:,:,:) 141 real(dp), DEV_CONTARRD intent(inout) :: recvbuf(:,:,:,:) 142 integer ,intent(in) :: sendsize, recvsize 143 integer ,intent(in) :: spaceComm 144 integer ,intent(out) :: ier 145 146! ************************************************************************* 147 148 ier=0 149#if defined HAVE_MPI 150 if (spaceComm /= MPI_COMM_SELF .and. spaceComm /= MPI_COMM_NULL) then 151! allgather xval on all proc. in spaceComm 152 call MPI_ALLTOALL(xval, sendsize, MPI_DOUBLE_PRECISION, recvbuf, & 153& recvsize, MPI_DOUBLE_PRECISION, spaceComm, ier) 154 else if (spaceComm == MPI_COMM_SELF) then 155 recvbuf=xval 156 end if 157#else 158 recvbuf=xval 159#endif 160 161end subroutine xmpi_alltoall_dp4d 162!!*** 163