TABLE OF CONTENTS
- ABINIT/xmpi_alltoallv_dp1d
- ABINIT/xmpi_alltoallv_dp1d2
- ABINIT/xmpi_alltoallv_dp2d
- ABINIT/xmpi_alltoallv_int2d
ABINIT/xmpi_alltoallv_dp1d [ Functions ]
NAME
xmpi_alltoallv_dp1d
FUNCTION
Sends data from all to all processes. Target: double precision one-dimensional arrays.
SOURCE
129 subroutine xmpi_alltoallv_dp1d(xval,sendcnts,sdispls,recvbuf,recvcnts,rdispls,comm,ier) 130 131 !Arguments------------------------- 132 real(dp), DEV_CONTARRD intent(in) :: xval(:) 133 real(dp), DEV_CONTARRD intent(inout) :: recvbuf(:) 134 integer, DEV_CONTARRD intent(in) :: sendcnts(:),sdispls(:),recvcnts(:) 135 integer,intent(in) :: comm, rdispls 136 integer,intent(out) :: ier 137 138 !Local variables------------------- 139 integer :: sc,sds,sdr 140 integer :: i 141 #if defined HAVE_MPI 142 integer, allocatable :: rdispls_on(:) 143 #endif 144 145 ! ********************************************************************* 146 147 ier=0 148 #if defined HAVE_MPI 149 if (comm /= MPI_COMM_SELF .and. comm /= MPI_COMM_NULL) then 150 ABI_STAT_MALLOC(rdispls_on,(size(sendcnts)), ier) 151 if (ier/= 0) call xmpi_abort(msg='error allocating rdispls_on in xmpi_alltoallv') 152 rdispls_on = 0 153 call MPI_ALLTOALLV(xval,sendcnts,sdispls,MPI_DOUBLE_PRECISION,recvbuf,& 154 & recvcnts,rdispls_on,MPI_DOUBLE_PRECISION,comm,ier) 155 ABI_FREE(rdispls_on) 156 else if (comm == MPI_COMM_SELF) then 157 #endif 158 sdr=rdispls;sds=0;if (size(sdispls)>0) sds=sdispls(1) 159 sc=size(xval);if (size(sendcnts)>0) sc=sendcnts(1) 160 !$OMP parallel do 161 do i = 1, sc 162 recvbuf(i)=xval(sds+i) 163 end do 164 #if defined HAVE_MPI 165 end if 166 #endif 167 168 end subroutine xmpi_alltoallv_dp1d
ABINIT/xmpi_alltoallv_dp1d2 [ Functions ]
NAME
xmpi_alltoallv_dp1d2
FUNCTION
Sends data from all to all processes. Target: double precision one-dimensional arrays.
SOURCE
181 subroutine xmpi_alltoallv_dp1d2(xval,sendcnts,sdispls,recvbuf,recvcnts,rdispls,comm,ier) 182 183 !Arguments------------------------- 184 real(dp), DEV_CONTARRD intent(in) :: xval(:) 185 real(dp), DEV_CONTARRD intent(inout) :: recvbuf(:) 186 integer, DEV_CONTARRD intent(in) :: sendcnts(:),sdispls(:),recvcnts(:),rdispls(:) 187 integer,intent(in) :: comm 188 integer,intent(out) :: ier 189 190 !Local variables------------------- 191 integer :: sc,sds,sdr 192 193 ! ********************************************************************* 194 195 ier=0 196 #if defined HAVE_MPI 197 if (comm /= MPI_COMM_SELF .and. comm /= MPI_COMM_NULL) then 198 call MPI_ALLTOALLV(xval,sendcnts,sdispls,MPI_DOUBLE_PRECISION,recvbuf,& 199 & recvcnts,rdispls,MPI_DOUBLE_PRECISION,comm,ier) 200 else if (comm == MPI_COMM_SELF) then 201 #endif 202 sds=0;if (size(sdispls)>0) sds=sdispls(1) 203 sdr=0;if (size(rdispls)>0) sdr=rdispls(1) 204 sc=size(xval);if (size(sendcnts)>0) sc=sendcnts(1) 205 recvbuf(sdr+1:sdr+sc)=xval(sds+1:sds+sc) 206 #if defined HAVE_MPI 207 end if 208 #endif 209 210 end subroutine xmpi_alltoallv_dp1d2
ABINIT/xmpi_alltoallv_dp2d [ Functions ]
NAME
xmpi_alltoallv_dp2d
FUNCTION
This module contains functions calling the MPI routine ALLTOALLV xmpi_alltoallv is the generic function.
COPYRIGHT
Copyright (C) 2001-2022 ABINIT group (AR,XG) This file is distributed under the terms of the GNU General Public License, see ~ABINIT/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
xval= buffer array sendcnts= number of sent elements sdispls= postions of values sent by the processor rdispls= positions of values received by the processor recvcnts= number of received elements comm= MPI communicator
OUTPUT
ier= exit status, a non-zero value meaning there is an error
SIDE EFFECTS
recvbuf= received buffer
PARENTS
CHILDREN
mpi_alltoallv
SOURCE
37 subroutine xmpi_alltoallv_dp2d(xval,sendcnts,sdispls,recvbuf,recvcnts,rdispls,comm,ier) 38 39 !Arguments------------------------- 40 real(dp), DEV_CONTARRD intent(in) :: xval(:,:) 41 real(dp), DEV_CONTARRD intent(inout) :: recvbuf(:,:) 42 integer , DEV_CONTARRD intent(in) :: sendcnts(:),sdispls(:),rdispls(:),recvcnts(:) 43 integer ,intent(in) :: comm 44 integer ,intent(out) :: ier 45 46 !Local variables------------------- 47 integer :: sc,sds,sdr,sz1 48 integer :: i 49 50 ! ********************************************************************* 51 52 ier=0 53 #if defined HAVE_MPI 54 if (comm /= MPI_COMM_SELF .and. comm /= MPI_COMM_NULL) then 55 ! allgather xval on all proc. in comm 56 call MPI_ALLTOALLV(xval,sendcnts,sdispls,MPI_DOUBLE_PRECISION,recvbuf,& 57 & recvcnts,rdispls,MPI_DOUBLE_PRECISION,comm,ier) 58 else if (comm == MPI_COMM_SELF) then 59 #endif 60 sz1=size(xval,1) 61 sds=0;if (size(sdispls)>0) sds=sdispls(1)/sz1 62 sdr=0;if (size(rdispls)>0) sdr=rdispls(1)/sz1 63 sc=size(xval,2);if (size(sendcnts)>0) sc=sendcnts(1)/sz1 64 !$OMP parallel do 65 do i=1,sc 66 recvbuf(:,sdr+i)=xval(:,sds+i) 67 end do 68 69 #if defined HAVE_MPI 70 end if 71 #endif 72 73 end subroutine xmpi_alltoallv_dp2d
ABINIT/xmpi_alltoallv_int2d [ Functions ]
NAME
xmpi_alltoallv_int2d
FUNCTION
Sends data from all to all processes. Target: two-dimensional integer arrays.
SOURCE
86 subroutine xmpi_alltoallv_int2d(xval,sendcnts,sdispls,recvbuf,recvcnts,rdispls,comm,ier) 87 88 !Arguments------------------------- 89 integer, DEV_CONTARRD intent(in) :: xval(:,:) 90 integer, DEV_CONTARRD intent(inout) :: recvbuf(:,:) 91 integer, DEV_CONTARRD intent(in) :: sendcnts(:),sdispls(:),rdispls(:),recvcnts(:) 92 integer,intent(in) :: comm 93 integer,intent(out) :: ier 94 95 !Local variables------------------- 96 integer :: sc,sds,sdr,sz1 97 98 ! ********************************************************************* 99 100 ier=0 101 #if defined HAVE_MPI 102 if (comm /= MPI_COMM_SELF .and. comm /= MPI_COMM_NULL) then 103 call MPI_ALLTOALLV(xval,sendcnts,sdispls,MPI_INTEGER,recvbuf,& 104 & recvcnts,rdispls,MPI_INTEGER,comm,ier) 105 else if (comm == MPI_COMM_SELF) then 106 #endif 107 sz1=size(xval,1) 108 sds=0;if (size(sdispls)>0) sds=sdispls(1)/sz1 109 sdr=0;if (size(rdispls)>0) sdr=rdispls(1)/sz1 110 sc=size(xval,2);if (size(sendcnts)>0) sc=sendcnts(1)/sz1 111 recvbuf(:,sdr+1:sdr+sc)=xval(:,sds+1:sds+sc) 112 #if defined HAVE_MPI 113 end if 114 #endif 115 116 end subroutine xmpi_alltoallv_int2d