TABLE OF CONTENTS
ABINIT/prep_wavef_sym_do [ Functions ]
NAME
prep_wavef_sym_do
FUNCTION
this routine associates waves functions by two as following E(G) = C(G) + D(G) E(-G) = C*(G) + iD*(G) the values are distributed 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 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 tab_proc = positions of opposited planewave coordinates in the list of the processors fft cwavef_alltoall = planewave coefficients of wavefunction ( initial of the processor + sended by other processors band) 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 from each processor fft rdispls_sym = postions of the received values by the processor from each processor fft
OUTPUT
ewavef_alltoall_sym = planewave coefficients of wavefunction initial of the processor + sended by other processors band + sended by other processors fft + and compisited if bandpp >1 index_wavef_send = index to send the values in blocks to the other processor fft
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_do(mpi_enreg,bandpp,nspinor,& 60 & ndatarecv,& 61 & ndatarecv_tot,ndatasend_sym,tab_proc,& 62 & cwavef_alltoall,& 63 & sendcounts_sym,sdispls_sym,& 64 & recvcounts_sym,rdispls_sym,& 65 & ewavef_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_do' 77 !End of the abilint section 78 79 implicit none 80 81 !Arguments ------------------------------------ 82 !scalars 83 integer,intent(in) :: bandpp,ndatarecv,ndatarecv_tot,ndatasend_sym 84 integer,intent(in) :: nspinor 85 type(mpi_type),intent(in) :: mpi_enreg 86 !arrays 87 integer,allocatable,intent(out) :: index_wavef_send(:) 88 integer,intent(in) :: rdispls_sym(:),recvcounts_sym(:) 89 integer,intent(in) :: sdispls_sym(:),sendcounts_sym(:) 90 integer,intent(in) :: tab_proc(:) 91 real(dp),intent(inout) :: cwavef_alltoall(2,ndatarecv*nspinor*bandpp) 92 real(dp),pointer :: ewavef_alltoall_sym(:,:) 93 94 !Local variables------------------------------- 95 !scalars 96 integer :: bandpp_sym,ibandpp,idatarecv,ideb_loc,idebc,idebd,idebe 97 integer :: ier,ifin_loc,ifinc,ifind,ifine,iproc,jbandpp,jsendloc 98 integer :: kbandpp,newspacecomm,nproc_fft 99 logical :: flag_compose 100 !arrays 101 integer,allocatable :: rdispls_sym_loc(:),recvcounts_sym_loc(:) 102 integer,allocatable :: sdispls_sym_loc(:),sendcounts_sym_loc(:) 103 real(dp),allocatable :: ewavef_alltoall_loc(:,:),ewavef_alltoall_send(:,:) 104 105 ! ********************************************************************* 106 107 !DEBUG 108 !write(std_out,*)' prep_wavef_sym_do : enter ' 109 !ENDDEBUG 110 111 !--------------------------------------------- 112 !Initialisation 113 !--------------------------------------------- 114 nproc_fft = mpi_enreg%nproc_fft 115 116 newspacecomm = mpi_enreg%comm_fft 117 118 if (modulo(bandpp,2)==0) then 119 bandpp_sym = bandpp/2 120 flag_compose = .TRUE. 121 else 122 bandpp_sym = bandpp 123 flag_compose = .FALSE. 124 end if 125 126 !--------------------------------------------- 127 !Allocation 128 !--------------------------------------------- 129 ABI_ALLOCATE(ewavef_alltoall_sym ,(2,ndatarecv_tot*bandpp_sym)) 130 ABI_ALLOCATE(ewavef_alltoall_loc ,(2,ndatarecv *bandpp_sym)) 131 ABI_ALLOCATE(ewavef_alltoall_send ,(2,ndatasend_sym*bandpp_sym)) 132 ABI_ALLOCATE(index_wavef_send ,( ndatasend_sym*bandpp_sym)) 133 134 ABI_ALLOCATE(sendcounts_sym_loc ,(nproc_fft)) 135 ABI_ALLOCATE(sdispls_sym_loc ,(nproc_fft)) 136 ABI_ALLOCATE(recvcounts_sym_loc ,(nproc_fft)) 137 ABI_ALLOCATE(rdispls_sym_loc ,(nproc_fft)) 138 139 140 !Initialisation 141 !-------------- 142 ewavef_alltoall_sym(:,:) =0. 143 ewavef_alltoall_loc(:,:) =0. 144 145 sendcounts_sym_loc(:) =0 146 sdispls_sym_loc(:) =0 147 recvcounts_sym_loc(:) =0 148 rdispls_sym_loc(:) =0 149 150 index_wavef_send(:) =0 151 152 153 !------------------------------------------------- 154 !We are bandpp blocks which we want to : 155 !associate by two (band_sym==bandpp/2) 156 !or not associate by two (band_sym==bandpp) 157 ! 158 !So We'll have got bandpp_sym blocks 159 !So we loop on the bandpp_sym blocks 160 !-------------------------------------------------- 161 162 do kbandpp=1,bandpp_sym 163 164 ! position of the two blocks 165 ! -------------------------- 166 ibandpp = (kbandpp-1) * 2 167 jbandpp = ibandpp + 1 168 169 idebe = (kbandpp-1) * ndatarecv_tot + 1 170 ifine = idebe + ndatarecv - 1 171 172 idebc = ibandpp * ndatarecv + 1 173 ifinc = idebc + ndatarecv - 1 174 175 idebd = jbandpp * ndatarecv + 1 176 ifind = idebd + ndatarecv - 1 177 178 ideb_loc = (kbandpp-1) * ndatarecv + 1 179 ifin_loc = ideb_loc + ndatarecv - 1 180 181 182 if (flag_compose) then 183 184 185 ! calcul ewavef(G) 186 ! ---------------- 187 ewavef_alltoall_sym(1,idebe:ifine) = & 188 & cwavef_alltoall(1,idebc:ifinc) & 189 & - cwavef_alltoall(2,idebd:ifind) 190 191 ewavef_alltoall_sym(2,idebe:ifine) = & 192 & cwavef_alltoall(2,idebc:ifinc) & 193 & + cwavef_alltoall(1,idebd:ifind) 194 195 ! calcul ewavef_loc(-G) 196 ! --------------------- 197 ewavef_alltoall_loc(1,ideb_loc:ifin_loc) = & 198 & cwavef_alltoall(1,idebc:ifinc) & 199 & + cwavef_alltoall(2,idebd:ifind) 200 201 ewavef_alltoall_loc(2,ideb_loc:ifin_loc) = & 202 & - cwavef_alltoall(2,idebc:ifinc) & 203 & + cwavef_alltoall(1,idebd:ifind) 204 else 205 206 ! calcul ewavef(G) 207 ! ---------------- 208 ewavef_alltoall_sym(1,idebe:ifine) = cwavef_alltoall(1,idebc:ifinc) 209 ewavef_alltoall_sym(2,idebe:ifine) = cwavef_alltoall(2,idebc:ifinc) 210 211 ! calcul ewavef_loc(-G) 212 ! --------------------- 213 ewavef_alltoall_loc(1,ideb_loc:ifin_loc) = cwavef_alltoall(1,idebc:ifinc) 214 ewavef_alltoall_loc(2,ideb_loc:ifin_loc) = - cwavef_alltoall(2,idebc:ifinc) 215 216 end if 217 218 end do 219 220 221 222 !------------------------------------------------------------------------ 223 !Creation of datas blocks for each processor fft from ewavef_alltoall_loc 224 !to send datas by blocks with a alltoall... 225 !------------------------------------------------------------------------ 226 227 !Position of the blocks 228 !---------------------- 229 jsendloc=0 230 do ibandpp=1,bandpp_sym 231 do iproc=1,nproc_fft 232 do idatarecv=1,ndatarecv 233 if (tab_proc(idatarecv)==(iproc-1)) then 234 jsendloc=jsendloc+1 235 index_wavef_send(jsendloc) = idatarecv + ndatarecv * (ibandpp-1) 236 end if 237 end do 238 end do 239 end do 240 241 !Classment 242 !---------- 243 ewavef_alltoall_send(:,:)=ewavef_alltoall_loc(:,index_wavef_send) 244 245 246 !------------------------------------------------- 247 !Calcul of the number of received and sended datas 248 !------------------------------------------------- 249 sendcounts_sym_loc = sendcounts_sym*2 250 recvcounts_sym_loc = recvcounts_sym*2 251 252 !------------------------------------------ 253 !Exchange of the datas ewavef_allto_all_loc 254 !------------------------------------------ 255 do ibandpp=1,bandpp_sym 256 257 ! ------------------------------------------------ 258 ! Deplacment of the sended datas because of bandpp 259 ! ------------------------------------------------ 260 sdispls_sym_loc(:) = sdispls_sym(:) + ndatasend_sym * (ibandpp-1) 261 sdispls_sym_loc = sdispls_sym_loc *2 262 263 ! -------------------------------------------------- 264 ! Deplacment of the received datas because of bandpp 265 ! -------------------------------------------------- 266 rdispls_sym_loc(:) = rdispls_sym(:) + ndatarecv_tot * (ibandpp-1) 267 rdispls_sym_loc = rdispls_sym_loc *2 268 269 270 call xmpi_alltoallv(& 271 & ewavef_alltoall_send(:,:) ,sendcounts_sym_loc,sdispls_sym_loc,& 272 & ewavef_alltoall_sym(:,:) ,recvcounts_sym_loc,rdispls_sym_loc,& 273 & newspacecomm,ier) 274 275 end do 276 277 !----------------------- 278 !Desallocation 279 !----------------------- 280 281 ABI_DEALLOCATE(sendcounts_sym_loc) 282 ABI_DEALLOCATE(recvcounts_sym_loc) 283 ABI_DEALLOCATE(sdispls_sym_loc) 284 ABI_DEALLOCATE(rdispls_sym_loc) 285 286 ABI_DEALLOCATE(ewavef_alltoall_loc) 287 ABI_DEALLOCATE(ewavef_alltoall_send) 288 289 end subroutine prep_wavef_sym_do