TABLE OF CONTENTS
ABINIT/prep_wavef_sym_undo [ Functions ]
NAME
prep_wavef_sym_undo
FUNCTION
this routine dissociates each wave function in two waves functions as following C(G) = ( E*(-G) + E(G))/2 D(G) = i*( E*(-G) - E(G))/2 the values are redistributed on the processors in function of the value of mpi_enreg%distribfft%tab_fftwf2_distrib( (-kg_k_gather(2,i) )
COPYRIGHT
INPUTS
mpi_enreg = informations about mpi parallelization bandpp = number of groups of couple of waves functions nspinor = number of spin ndatarecv = number of values received by the processor and sended by the other processors band ndatarecv_tot = total number of received values (ndatarecv + number of received opposited planewave coordinates) ndatasend_sym = number of sended values to the processors fft to create opposited planewave coordinates idatarecv0 = position of the planewave coordinates (0,0,0) sendcounts_sym = number of sended values by the processor to each processor fft sdispls_sym = postions of the sended values by the processor to each processor fft recvcounts_sym = number of the received values by the processor to each processor fft gwavef_alltoall_sym = planewave coefficients of wavefunction initial of the processor + sended by other processors band + sended by other processors fft + and composited if bandpp >1 index_wavef_send = index to send the values by block to the other processor fft
OUTPUT
gwavef_alltoall = planewave coefficients of wavefunction ( for of the processor + to send to other processors band)
SIDE EFFECTS
PARENTS
prep_fourwf,prep_getghc
CHILDREN
xmpi_alltoallv
SOURCE
53 #if defined HAVE_CONFIG_H 54 #include "config.h" 55 #endif 56 57 #include "abi_common.h" 58 59 subroutine prep_wavef_sym_undo(mpi_enreg,bandpp,nspinor,& 60 & ndatarecv,& 61 & ndatarecv_tot,ndatasend_sym,idatarecv0,& 62 & gwavef_alltoall,& 63 & sendcounts_sym,sdispls_sym,& 64 & recvcounts_sym,rdispls_sym,& 65 & gwavef_alltoall_sym,& 66 & index_wavef_send) 67 68 use defs_basis 69 use defs_abitypes 70 use m_profiling_abi 71 use m_xmpi 72 73 !This section has been created automatically by the script Abilint (TD). 74 !Do not modify the following lines by hand. 75 #undef ABI_FUNC 76 #define ABI_FUNC 'prep_wavef_sym_undo' 77 !End of the abilint section 78 79 implicit none 80 81 !Arguments ------------------------------------ 82 !scalars 83 integer,intent(in) :: bandpp,idatarecv0,ndatarecv,ndatarecv_tot,ndatasend_sym 84 integer,intent(in) :: nspinor 85 type(mpi_type),intent(in) :: mpi_enreg 86 !arrays 87 integer,intent(in) :: index_wavef_send(:),rdispls_sym(:),recvcounts_sym(:) 88 integer,intent(in) :: sdispls_sym(:),sendcounts_sym(:) 89 real(dp),intent(inout) :: gwavef_alltoall(2,ndatarecv*nspinor*bandpp) 90 real(dp),intent(inout) :: gwavef_alltoall_sym(:,:) 91 92 !Local variables------------------------------- 93 !scalars 94 integer :: bandpp_sym,ibandpp,ideb_loc,idebc,idebd 95 integer :: idebe,ier,ifin_loc,ifinc,ifind,ifine 96 integer :: jbandpp,kbandpp,newspacecomm,nproc_fft 97 logical :: flag_compose 98 !arrays 99 integer,allocatable :: rdispls_sym_loc(:),recvcounts_sym_loc(:) 100 integer,allocatable :: sdispls_sym_loc(:),sendcounts_sym_loc(:) 101 real(dp),allocatable :: gwavef_alltoall_loc(:,:),gwavef_alltoall_rcv(:,:) 102 103 ! ********************************************************************* 104 105 !DEBUG 106 !write(std_out,*)' prep_wavef_sym_undo : enter ' 107 !ENDDEBUG 108 109 110 !--------------------------------------------- 111 !Initialisation 112 !--------------------------------------------- 113 nproc_fft = mpi_enreg%nproc_fft 114 115 newspacecomm = mpi_enreg%comm_fft 116 117 if (modulo(bandpp,2)==0) then 118 bandpp_sym = bandpp/2 119 flag_compose = .TRUE. 120 else 121 bandpp_sym = bandpp 122 flag_compose = .FALSE. 123 end if 124 125 !--------------------------------------------- 126 !Allocation 127 !--------------------------------------------- 128 ABI_ALLOCATE(gwavef_alltoall_loc ,(2,ndatarecv *bandpp_sym)) 129 ABI_ALLOCATE(gwavef_alltoall_rcv ,(2,ndatasend_sym *bandpp_sym)) 130 131 ABI_ALLOCATE(sendcounts_sym_loc ,(nproc_fft)) 132 ABI_ALLOCATE(sdispls_sym_loc ,(nproc_fft)) 133 ABI_ALLOCATE(recvcounts_sym_loc ,(nproc_fft)) 134 ABI_ALLOCATE(rdispls_sym_loc ,(nproc_fft)) 135 136 137 !--------------------------------------------- 138 !Initialisation 139 !--------------------------------------------- 140 gwavef_alltoall_loc(:,:) =0. 141 142 sendcounts_sym_loc(:) =0 143 sdispls_sym_loc(:) =0 144 recvcounts_sym_loc(:) =0 145 rdispls_sym_loc(:) =0 146 147 148 !------------------------------------------------- 149 !Calcul of number of the sended and received datas 150 !------------------------------------------------- 151 sendcounts_sym_loc = sendcounts_sym*2 152 recvcounts_sym_loc = recvcounts_sym*2 153 154 !---------------------------------------------------- 155 !Sending of the values 156 !---------------------------------------------------- 157 do ibandpp = 1,bandpp_sym 158 159 ! ------------------------------------------------- 160 ! Deplacment of the sended values because of bandpp 161 ! ------------------------------------------------- 162 sdispls_sym_loc(:) = sdispls_sym(:) + ndatasend_sym * (ibandpp-1) 163 sdispls_sym_loc = sdispls_sym_loc *2 164 165 ! --------------------------------------------------- 166 ! Deplacment of the received values because of bandpp 167 ! --------------------------------------------------- 168 rdispls_sym_loc(:) = rdispls_sym(:) + ndatarecv_tot * (ibandpp-1) 169 rdispls_sym_loc = rdispls_sym_loc *2 170 171 172 call xmpi_alltoallv(& 173 & gwavef_alltoall_sym(:,:) ,recvcounts_sym_loc,rdispls_sym_loc,& 174 & gwavef_alltoall_rcv(:,:) ,sendcounts_sym_loc,sdispls_sym_loc,& 175 & newspacecomm,ier) 176 177 end do 178 179 180 !---------------------- 181 !Dispatching the blocks 182 !---------------------- 183 gwavef_alltoall_loc(:,index_wavef_send(:)) = gwavef_alltoall_rcv(:,:) 184 185 !---------------------- 186 !Case -kg = [0 0 0] 187 !---------------------- 188 if (idatarecv0/=-1) then 189 do kbandpp=1,bandpp_sym 190 gwavef_alltoall_loc(:,(kbandpp-1)*ndatarecv + idatarecv0)= & 191 gwavef_alltoall_sym(:,(kbandpp-1)*ndatarecv_tot + idatarecv0) 192 end do 193 end if 194 195 !--------------------------------------------------- 196 !Build of hwavef_alltoall 197 ! 198 !We have got : 199 !bandpp_sym blocks to dissociate 200 !or bandpp_sym blokcs to not dissociate 201 !-------------------------------------------------- 202 do kbandpp=1,bandpp_sym 203 204 ! position of the 2 blocks 205 ! ---------------------------------- 206 ibandpp = (kbandpp-1) * 2 207 jbandpp = ibandpp + 1 208 209 idebe = (kbandpp-1) * ndatarecv_tot + 1 210 ifine = idebe + ndatarecv - 1 211 212 idebc = ibandpp * ndatarecv + 1 213 ifinc = idebc + ndatarecv - 1 214 215 idebd = jbandpp * ndatarecv + 1 216 ifind = idebd + ndatarecv - 1 217 218 ideb_loc = (kbandpp-1) * ndatarecv + 1 219 ifin_loc = ideb_loc + ndatarecv - 1 220 221 222 if (flag_compose) then 223 224 ! calcul cwavef(G) 225 ! ---------------- 226 gwavef_alltoall(1,idebc:ifinc) = gwavef_alltoall_sym(1,idebe:ifine) & 227 & + gwavef_alltoall_loc(1,ideb_loc:ifin_loc) 228 gwavef_alltoall(2,idebc:ifinc) = gwavef_alltoall_sym(2,idebe:ifine) & 229 & - gwavef_alltoall_loc(2,ideb_loc:ifin_loc) 230 231 ! calcul dwavef(G) 232 ! ------------------ 233 gwavef_alltoall(1,idebd:ifind) = gwavef_alltoall_sym(2,idebe:ifine) & 234 & + gwavef_alltoall_loc(2,ideb_loc:ifin_loc) 235 gwavef_alltoall(2,idebd:ifind) = - gwavef_alltoall_sym(1,idebe:ifine) & 236 & + gwavef_alltoall_loc(1,ideb_loc:ifin_loc) 237 else 238 239 ! calcul cwavef(G) 240 ! ---------------- 241 gwavef_alltoall(1,idebc:ifinc) = gwavef_alltoall_sym(1,idebe:ifine) & 242 & + gwavef_alltoall_loc(1,ideb_loc:ifin_loc) 243 gwavef_alltoall(2,idebc:ifinc) = gwavef_alltoall_sym(2,idebe:ifine) & 244 & - gwavef_alltoall_loc(2,ideb_loc:ifin_loc) 245 end if 246 247 end do 248 249 !We divise by two 250 gwavef_alltoall(:,:) = gwavef_alltoall(:,:)/2 251 252 !----------------------- 253 !Desallocation 254 !----------------------- 255 256 ABI_DEALLOCATE(sendcounts_sym_loc) 257 ABI_DEALLOCATE(recvcounts_sym_loc) 258 ABI_DEALLOCATE(sdispls_sym_loc) 259 ABI_DEALLOCATE(rdispls_sym_loc) 260 261 ABI_DEALLOCATE(gwavef_alltoall_loc) 262 ABI_DEALLOCATE(gwavef_alltoall_rcv) 263 264 end subroutine prep_wavef_sym_undo